Block-Structured AMR Software Framework
Loading...
Searching...
No Matches
AMReX_BaseFab.H
Go to the documentation of this file.
1#ifndef AMREX_BASEFAB_H_
2#define AMREX_BASEFAB_H_
3#include <AMReX_Config.H>
4
5#include <AMReX_Algorithm.H>
6#include <AMReX_Extension.H>
7#include <AMReX_BLassert.H>
8#include <AMReX_Array.H>
9#include <AMReX_Box.H>
10#include <AMReX_Loop.H>
11#include <AMReX_BoxList.H>
12#include <AMReX_BArena.H>
13#include <AMReX_CArena.H>
14#include <AMReX_DataAllocator.H>
15#include <AMReX_REAL.H>
16#include <AMReX_BLProfiler.H>
17#include <AMReX_BoxIterator.H>
18#include <AMReX_MakeType.H>
19#include <AMReX_Utility.H>
20#include <AMReX_Reduce.H>
21#include <AMReX_Gpu.H>
22#include <AMReX_Scan.H>
23#include <AMReX_Math.H>
24#include <AMReX_OpenMP.H>
25#include <AMReX_MemPool.H>
26#include <AMReX_TypeTraits.H>
27
28#include <cmath>
29#include <cstdlib>
30#include <algorithm>
31#include <limits>
32#include <climits>
33#include <array>
34#include <type_traits>
35#include <memory>
36#include <atomic>
37#include <utility>
38
39namespace amrex
40{
41
43
44extern std::atomic<Long> atomic_total_bytes_allocated_in_fabs;
45extern std::atomic<Long> atomic_total_bytes_allocated_in_fabs_hwm;
46extern std::atomic<Long> atomic_total_cells_allocated_in_fabs;
47extern std::atomic<Long> atomic_total_cells_allocated_in_fabs_hwm;
48extern Long private_total_bytes_allocated_in_fabs;
49extern Long private_total_bytes_allocated_in_fabs_hwm;
50extern Long private_total_cells_allocated_in_fabs;
51extern Long private_total_cells_allocated_in_fabs_hwm;
52#ifdef AMREX_USE_OMP
53#pragma omp threadprivate(private_total_bytes_allocated_in_fabs)
54#pragma omp threadprivate(private_total_bytes_allocated_in_fabs_hwm)
55#pragma omp threadprivate(private_total_cells_allocated_in_fabs)
56#pragma omp threadprivate(private_total_cells_allocated_in_fabs_hwm)
57#endif
58
60
66void update_fab_stats (Long n, Long s, std::size_t szt) noexcept;
67
68void BaseFab_Initialize ();
69void BaseFab_Finalize ();
70
71struct SrcComp {
73 explicit SrcComp (int ai) noexcept : i(ai) {}
74 int i;
75};
76
77struct DestComp {
79 explicit DestComp (int ai) noexcept : i(ai) {}
80 int i;
81};
82
83struct NumComps {
85 explicit NumComps (int an) noexcept : n(an) {}
86 int n;
87};
88
89template <typename T>
92makeArray4 (T* p, Box const& bx, int ncomp) noexcept
93{
94 return Array4<T>{p, amrex::begin(bx), amrex::end(bx), ncomp};
95}
96
97template <typename T>
98std::enable_if_t<std::is_arithmetic_v<T>>
99placementNew (T* const /*ptr*/, Long /*n*/)
100{}
101
102template <typename T>
103std::enable_if_t<std::is_trivially_default_constructible_v<T>
104 && !std::is_arithmetic_v<T>>
105placementNew (T* const ptr, Long n)
106{
107 for (Long i = 0; i < n; ++i) {
108 new (ptr+i) T;
109 }
110}
111
112template <typename T>
113std::enable_if_t<!std::is_trivially_default_constructible_v<T>>
114placementNew (T* const ptr, Long n)
115{
117 {
118 new (ptr+i) T;
119 });
120}
121
122template <typename T>
123std::enable_if_t<std::is_trivially_destructible_v<T>>
124placementDelete (T* const /*ptr*/, Long /*n*/)
125{}
126
127template <typename T>
128std::enable_if_t<!std::is_trivially_destructible_v<T>>
129placementDelete (T* const ptr, Long n)
130{
132 {
133 (ptr+i)->~T();
134 });
135}
136
187template <class T>
189 : public DataAllocator
190{
191public:
192
193 template <class U> friend class BaseFab;
194
195 using value_type = T;
196
198 BaseFab () noexcept = default;
199
200 explicit BaseFab (Arena* ar) noexcept;
201
202 BaseFab (const Box& bx, int n, Arena* ar);
203
205 explicit BaseFab (const Box& bx, int n = 1, bool alloc = true,
206 bool shared = false, Arena* ar = nullptr);
207
208 BaseFab (const BaseFab<T>& rhs, MakeType make_type, int scomp, int ncomp);
209
215 BaseFab (const Box& bx, int ncomp, T* p);
216 BaseFab (const Box& bx, int ncomp, T const* p);
217
218 explicit BaseFab (Array4<T> const& a) noexcept;
219
220 explicit BaseFab (Array4<T> const& a, IndexType t) noexcept;
221
222 explicit BaseFab (Array4<T const> const& a) noexcept;
223
224 explicit BaseFab (Array4<T const> const& a, IndexType t) noexcept;
225
227 virtual ~BaseFab () noexcept;
228
229 BaseFab (const BaseFab<T>& rhs) = delete;
230 BaseFab<T>& operator= (const BaseFab<T>& rhs) = delete;
231
232 BaseFab (BaseFab<T>&& rhs) noexcept;
233 BaseFab<T>& operator= (BaseFab<T>&& rhs) noexcept;
234
235 template <RunOn run_on AMREX_DEFAULT_RUNON>
236 BaseFab& operator= (T const&) noexcept;
237
238 static void Initialize();
239 static void Finalize();
240
254 void resize (const Box& b, int N = 1, Arena* ar = nullptr);
255
256 template <class U=T, std::enable_if_t<std::is_trivially_destructible_v<U>,int> = 0>
257 [[nodiscard]] Elixir elixir () noexcept;
258
263 void clear ();
264
266 [[nodiscard]] std::unique_ptr<T,DataDeleter> release () noexcept;
267
269 [[nodiscard]] std::size_t nBytes () const noexcept { return this->truesize*sizeof(T); }
270
271 [[nodiscard]] std::size_t nBytesOwned () const noexcept {
272 return (this->ptr_owner) ? nBytes() : 0;
273 }
274
276 [[nodiscard]] std::size_t nBytes (const Box& bx, int ncomps) const noexcept
277 { return bx.numPts() * sizeof(T) * ncomps; }
278
280 [[nodiscard]] int nComp () const noexcept { return this->nvar; }
281
283 [[nodiscard]] const int* nCompPtr() const noexcept {
284 return &(this->nvar);
285 }
286
288 [[nodiscard]] Long numPts () const noexcept { return this->domain.numPts(); }
289
291 [[nodiscard]] Long size () const noexcept { return this->nvar*this->domain.numPts(); }
292
294 [[nodiscard]] const Box& box () const noexcept { return this->domain; }
295
300 [[nodiscard]] IntVect length () const noexcept { return this->domain.length(); }
301
306 [[nodiscard]] const IntVect& smallEnd () const noexcept { return this->domain.smallEnd(); }
307
309 [[nodiscard]] const IntVect& bigEnd () const noexcept { return this->domain.bigEnd(); }
310
319 [[nodiscard]] const int* loVect () const noexcept { return this->domain.loVect(); }
320
329 [[nodiscard]] const int* hiVect () const noexcept { return this->domain.hiVect(); }
330
335 [[nodiscard]] bool contains (const BaseFab<T>& fab) const noexcept
336 {
337 return box().contains(fab.box()) && this->nvar >= fab.nvar;
338 }
339
344 [[nodiscard]] bool contains (const Box& bx) const noexcept { return box().contains(bx); }
345
355 [[nodiscard]] T* dataPtr (int n = 0) noexcept {
356 if (this->dptr) {
357 return &(this->dptr[n*this->domain.numPts()]);
358 } else {
359 return nullptr;
360 }
361 }
362
364 [[nodiscard]] const T* dataPtr (int n = 0) const noexcept {
365 if (this->dptr) {
366 return &(this->dptr[n*this->domain.numPts()]);
367 } else {
368 return nullptr;
369 }
370 }
371
372 [[nodiscard]] T* dataPtr (const IntVect& p, int n = 0) noexcept;
373
374 [[nodiscard]] const T* dataPtr (const IntVect& p, int n = 0) const noexcept;
375
376 void setPtr (T* p, Long sz) noexcept { AMREX_ASSERT(this->dptr == nullptr && this->truesize == 0); this->dptr = p; this->truesize = sz; }
377
378 void prefetchToHost () const noexcept;
379 void prefetchToDevice () const noexcept;
380
381 [[nodiscard]] AMREX_FORCE_INLINE
382 Array4<T const> array () const noexcept
383 {
384 return makeArray4<T const>(this->dptr, this->domain, this->nvar);
385 }
386
387 [[nodiscard]] AMREX_FORCE_INLINE
388 Array4<T const> array (int start_comp) const noexcept
389 {
390 return Array4<T const>(makeArray4<T const>(this->dptr, this->domain, this->nvar),start_comp);
391 }
392
393 [[nodiscard]] AMREX_FORCE_INLINE
394 Array4<T const> array (int start_comp, int num_comps) const noexcept
395 {
396 return Array4<T const>(makeArray4<T const>(this->dptr, this->domain, this->nvar), start_comp, num_comps);
397 }
398
399 [[nodiscard]] AMREX_FORCE_INLINE
400 Array4<T> array () noexcept
401 {
402 return makeArray4<T>(this->dptr, this->domain, this->nvar);
403 }
404
405 [[nodiscard]] AMREX_FORCE_INLINE
406 Array4<T> array (int start_comp) noexcept
407 {
408 return Array4<T>(makeArray4<T>(this->dptr, this->domain, this->nvar),start_comp);
409 }
410
411 [[nodiscard]] AMREX_FORCE_INLINE
412 Array4<T> array (int start_comp, int num_comps) noexcept
413 {
414 return Array4<T>(makeArray4<T>(this->dptr, this->domain, this->nvar), start_comp, num_comps);
415 }
416
417 [[nodiscard]] AMREX_FORCE_INLINE
418 Array4<T const> const_array () const noexcept
419 {
420 return makeArray4<T const>(this->dptr, this->domain, this->nvar);
421 }
422
423 [[nodiscard]] AMREX_FORCE_INLINE
424 Array4<T const> const_array (int start_comp) const noexcept
425 {
426 return Array4<T const>(makeArray4<T const>(this->dptr, this->domain, this->nvar),start_comp);
427 }
428
429 [[nodiscard]] AMREX_FORCE_INLINE
430 Array4<T const> const_array (int start_comp, int num_comps) const noexcept
431 {
432 return Array4<T const>(makeArray4<T const>(this->dptr, this->domain, this->nvar), start_comp, num_comps);
433 }
434
436 [[nodiscard]] bool isAllocated () const noexcept { return this->dptr != nullptr; }
437
444 [[nodiscard]] T& operator() (const IntVect& p, int N) noexcept;
445
447 [[nodiscard]] T& operator() (const IntVect& p) noexcept;
448
450 [[nodiscard]] const T& operator() (const IntVect& p, int N) const noexcept;
451
453 [[nodiscard]] const T& operator() (const IntVect& p) const noexcept;
454
460 void getVal (T* data, const IntVect& pos, int N, int numcomp) const noexcept;
462 void getVal (T* data, const IntVect& pos) const noexcept;
463
464 template <RunOn run_on AMREX_DEFAULT_RUNON,
465 class U=T, std::enable_if_t<std::is_same_v<U,float> || std::is_same_v<U,double>,int> FOO = 0>
466 void fill_snan () noexcept;
467
474 template <RunOn run_on AMREX_DEFAULT_RUNON>
475 void setVal (T const& x, const Box& bx, int dcomp, int ncomp) noexcept;
477 template <RunOn run_on AMREX_DEFAULT_RUNON>
478 void setVal (T const& x, const Box& bx, int N = 0) noexcept;
480 template <RunOn run_on AMREX_DEFAULT_RUNON>
481 void setVal (T const& x, int N) noexcept;
482
483 template <RunOn run_on AMREX_DEFAULT_RUNON>
484 void setValIfNot (T const& val, const Box& bx, const BaseFab<int>& mask, int nstart, int num) noexcept;
485
491 template <RunOn run_on AMREX_DEFAULT_RUNON>
492 void setComplement (T const& x, const Box& b, int ns, int num) noexcept;
493
510 template <RunOn run_on AMREX_DEFAULT_RUNON>
511 BaseFab<T>& copy (const BaseFab<T>& src, const Box& srcbox, int srccomp,
512 const Box& destbox, int destcomp, int numcomp) noexcept;
513
520 template <RunOn run_on AMREX_DEFAULT_RUNON>
521 BaseFab<T>& copy (const BaseFab<T>& src, int srccomp, int destcomp,
522 int numcomp = 1) noexcept;
529 template <RunOn run_on AMREX_DEFAULT_RUNON>
530 BaseFab<T>& copy (const BaseFab<T>& src, const Box& destbox) noexcept;
531
533 template <RunOn run_on AMREX_DEFAULT_RUNON>
534 std::size_t copyToMem (const Box& srcbox, int srccomp,
535 int numcomp, void* dst) const noexcept;
536
538 template <RunOn run_on AMREX_DEFAULT_RUNON, typename BUF = T>
539 std::size_t copyFromMem (const Box& dstbox, int dstcomp,
540 int numcomp, const void* src) noexcept;
541
543 template <RunOn run_on AMREX_DEFAULT_RUNON, typename BUF = T>
544 std::size_t addFromMem (const Box& dstbox, int dstcomp,
545 int numcomp, const void* src) noexcept;
546
552 BaseFab<T>& shift (const IntVect& v) noexcept;
558 BaseFab<T>& shift (int idir, int n_cell) noexcept;
564 BaseFab<T>& shiftHalf (int dir, int n_cell) noexcept;
570 BaseFab<T>& shiftHalf (const IntVect& v) noexcept;
571
572 template <RunOn run_on AMREX_DEFAULT_RUNON>
573 [[nodiscard]] Real norminfmask (const Box& subbox, const BaseFab<int>& mask, int scomp=0, int ncomp=1) const noexcept;
574
581 template <RunOn run_on AMREX_DEFAULT_RUNON>
582 [[nodiscard]] Real norm (int p, int scomp = 0, int numcomp = 1) const;
583
585 template <RunOn run_on AMREX_DEFAULT_RUNON>
586 [[nodiscard]] Real norm (const Box& subbox, int p, int scomp = 0, int numcomp = 1) const;
588 template <RunOn run_on AMREX_DEFAULT_RUNON>
589 void abs () noexcept;
591 template <RunOn run_on AMREX_DEFAULT_RUNON>
592 void abs (int comp, int numcomp=1) noexcept;
596 template <RunOn run_on AMREX_DEFAULT_RUNON>
597 void abs (const Box& subbox, int comp = 0, int numcomp=1) noexcept;
601 template <RunOn run_on AMREX_DEFAULT_RUNON>
602 [[nodiscard]] T min (int comp = 0) const noexcept;
606 template <RunOn run_on AMREX_DEFAULT_RUNON>
607 [[nodiscard]] T min (const Box& subbox, int comp = 0) const noexcept;
611 template <RunOn run_on AMREX_DEFAULT_RUNON>
612 [[nodiscard]] T max (int comp = 0) const noexcept;
616 template <RunOn run_on AMREX_DEFAULT_RUNON>
617 [[nodiscard]] T max (const Box& subbox, int comp = 0) const noexcept;
621 template <RunOn run_on AMREX_DEFAULT_RUNON>
622 [[nodiscard]] std::pair<T,T> minmax (int comp = 0) const noexcept;
626 template <RunOn run_on AMREX_DEFAULT_RUNON>
627 [[nodiscard]] std::pair<T,T> minmax (const Box& subbox, int comp = 0) const noexcept;
631 template <RunOn run_on AMREX_DEFAULT_RUNON>
632 [[nodiscard]] T maxabs (int comp = 0) const noexcept;
636 template <RunOn run_on AMREX_DEFAULT_RUNON>
637 [[nodiscard]] T maxabs (const Box& subbox, int comp = 0) const noexcept;
638
639 /*(
640 * \return location of given value
641 */
642 template <RunOn run_on AMREX_DEFAULT_RUNON>
643 [[nodiscard]] IntVect indexFromValue (const Box& subbox, int comp, T const& value) const noexcept;
644
648 template <RunOn run_on AMREX_DEFAULT_RUNON>
649 [[nodiscard]] IntVect minIndex (int comp = 0) const noexcept;
654 template <RunOn run_on AMREX_DEFAULT_RUNON>
655 [[nodiscard]] IntVect minIndex (const Box& subbox, int comp = 0) const noexcept;
660 template <RunOn run_on AMREX_DEFAULT_RUNON>
661 void minIndex (const Box& subbox, Real& min_val, IntVect& min_idx, int comp = 0) const noexcept;
662
666 template <RunOn run_on AMREX_DEFAULT_RUNON>
667 [[nodiscard]] IntVect maxIndex (int comp = 0) const noexcept;
672 template <RunOn run_on AMREX_DEFAULT_RUNON>
673 [[nodiscard]] IntVect maxIndex (const Box& subbox, int comp = 0) const noexcept;
678 template <RunOn run_on AMREX_DEFAULT_RUNON>
679 void maxIndex (const Box& subbox, Real& max_value, IntVect& max_idx, int comp = 0) const noexcept;
680
687 template <RunOn run_on AMREX_DEFAULT_RUNON>
688 int maskLT (BaseFab<int>& mask, T const& val, int comp = 0) const noexcept;
690 template <RunOn run_on AMREX_DEFAULT_RUNON>
691 int maskLE (BaseFab<int>& mask, T const& val, int comp = 0) const noexcept;
692
694 template <RunOn run_on AMREX_DEFAULT_RUNON>
695 int maskEQ (BaseFab<int>& mask, T const& val, int comp = 0) const noexcept;
697 template <RunOn run_on AMREX_DEFAULT_RUNON>
698 int maskGT (BaseFab<int>& mask, T const& val, int comp = 0) const noexcept;
700 template <RunOn run_on AMREX_DEFAULT_RUNON>
701 int maskGE (BaseFab<int>& mask, T const& val, int comp = 0) const noexcept;
702
704 template <RunOn run_on AMREX_DEFAULT_RUNON>
705 [[nodiscard]] T sum (int comp, int numcomp = 1) const noexcept;
707 template <RunOn run_on AMREX_DEFAULT_RUNON>
708 [[nodiscard]] T sum (const Box& subbox, int comp, int numcomp = 1) const noexcept;
709
711 template <RunOn run_on AMREX_DEFAULT_RUNON>
712 BaseFab<T>& invert (T const& r, const Box& b, int comp=0, int numcomp=1) noexcept;
714 template <RunOn run_on AMREX_DEFAULT_RUNON>
715 BaseFab<T>& invert (T const& r, int comp, int numcomp=1) noexcept;
716
718 template <RunOn run_on AMREX_DEFAULT_RUNON>
719 BaseFab<T>& negate (const Box& b, int comp=0, int numcomp=1) noexcept;
721 template <RunOn run_on AMREX_DEFAULT_RUNON>
722 BaseFab<T>& negate (int comp, int numcomp=1) noexcept;
723
725 template <RunOn run_on AMREX_DEFAULT_RUNON>
726 BaseFab<T>& plus (T const& r, const Box& b, int comp=0, int numcomp=1) noexcept;
727
729 template <RunOn run_on AMREX_DEFAULT_RUNON>
730 BaseFab<T>& plus (T const& r, int comp, int numcomp=1) noexcept;
731
737 template <RunOn run_on AMREX_DEFAULT_RUNON>
738 BaseFab<T>& plus (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp=1) noexcept;
744 template <RunOn run_on AMREX_DEFAULT_RUNON>
745 BaseFab<T>& plus (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp, int numcomp=1) noexcept;
750 template <RunOn run_on AMREX_DEFAULT_RUNON>
751 BaseFab<T>& plus (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
752 int srccomp, int destcomp, int numcomp=1) noexcept;
753
755 template <RunOn run_on AMREX_DEFAULT_RUNON>
756 BaseFab<T>& atomicAdd (const BaseFab<T>& x) noexcept;
757
763 template <RunOn run_on AMREX_DEFAULT_RUNON>
764 BaseFab<T>& atomicAdd (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp=1) noexcept;
770 template <RunOn run_on AMREX_DEFAULT_RUNON>
771 BaseFab<T>& atomicAdd (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
772 int numcomp=1) noexcept;
777 template <RunOn run_on AMREX_DEFAULT_RUNON>
778 BaseFab<T>& atomicAdd (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
779 int srccomp, int destcomp, int numcomp=1) noexcept;
780
786 template <RunOn run_on AMREX_DEFAULT_RUNON>
787 BaseFab<T>& lockAdd (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
788 int srccomp, int destcomp, int numcomp) noexcept;
789
791 template <RunOn run_on AMREX_DEFAULT_RUNON>
792 BaseFab<T>& saxpy (T a, const BaseFab<T>& x, const Box& srcbox, const Box& destbox,
793 int srccomp, int destcomp, int numcomp=1) noexcept;
795 template <RunOn run_on AMREX_DEFAULT_RUNON>
796 BaseFab<T>& saxpy (T a, const BaseFab<T>& x) noexcept;
797
799 template <RunOn run_on AMREX_DEFAULT_RUNON>
800 BaseFab<T>& xpay (T a, const BaseFab<T>& x, const Box& srcbox, const Box& destbox,
801 int srccomp, int destcomp, int numcomp=1) noexcept;
802
804 template <RunOn run_on AMREX_DEFAULT_RUNON>
805 BaseFab<T>& addproduct (const Box& destbox, int destcomp, int numcomp,
806 const BaseFab<T>& src1, int comp1,
807 const BaseFab<T>& src2, int comp2) noexcept;
808
814 template <RunOn run_on AMREX_DEFAULT_RUNON>
815 BaseFab<T>& minus (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp=1) noexcept;
821 template <RunOn run_on AMREX_DEFAULT_RUNON>
822 BaseFab<T>& minus (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
823 int numcomp=1) noexcept;
828 template <RunOn run_on AMREX_DEFAULT_RUNON>
829 BaseFab<T>& minus (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
830 int srccomp, int destcomp, int numcomp=1) noexcept;
831
833 template <RunOn run_on AMREX_DEFAULT_RUNON>
834 BaseFab<T>& mult (T const& r, int comp, int numcomp=1) noexcept;
838 template <RunOn run_on AMREX_DEFAULT_RUNON>
839 BaseFab<T>& mult (T const& r, const Box& b, int comp=0, int numcomp=1) noexcept;
840
846 template <RunOn run_on AMREX_DEFAULT_RUNON>
847 BaseFab<T>& mult (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp=1) noexcept;
848
854 template <RunOn run_on AMREX_DEFAULT_RUNON>
855 BaseFab<T>& mult (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
856 int numcomp=1) noexcept;
857
862 template <RunOn run_on AMREX_DEFAULT_RUNON>
863 BaseFab<T>& mult (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
864 int srccomp, int destcomp, int numcomp=1) noexcept;
865
867 template <RunOn run_on AMREX_DEFAULT_RUNON>
868 BaseFab<T>& divide (T const& r, int comp, int numcomp=1) noexcept;
869
871 template <RunOn run_on AMREX_DEFAULT_RUNON>
872 BaseFab<T>& divide (T const& r, const Box& b, int comp=0, int numcomp=1) noexcept;
873
880 template <RunOn run_on AMREX_DEFAULT_RUNON>
881 BaseFab<T>& divide (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp=1) noexcept;
887 template <RunOn run_on AMREX_DEFAULT_RUNON>
888 BaseFab<T>& divide (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
889 int numcomp=1) noexcept;
894 template <RunOn run_on AMREX_DEFAULT_RUNON>
895 BaseFab<T>& divide (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
896 int srccomp, int destcomp, int numcomp=1) noexcept;
900 template <RunOn run_on AMREX_DEFAULT_RUNON>
901 BaseFab<T>& protected_divide (const BaseFab<T>& src) noexcept;
902
910 template <RunOn run_on AMREX_DEFAULT_RUNON>
911 BaseFab<T>& protected_divide (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp=1) noexcept;
912
919 template <RunOn run_on AMREX_DEFAULT_RUNON>
920 BaseFab<T>& protected_divide (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
921 int numcomp=1) noexcept;
922
928 template <RunOn run_on AMREX_DEFAULT_RUNON>
929 BaseFab<T>& protected_divide (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
930 int srccomp, int destcomp, int numcomp=1) noexcept;
931
942 template <RunOn run_on AMREX_DEFAULT_RUNON>
943 BaseFab<T>& linInterp (const BaseFab<T>& f1, const Box& b1, int comp1,
944 const BaseFab<T>& f2, const Box& b2, int comp2,
945 Real t1, Real t2, Real t,
946 const Box& b, int comp, int numcomp = 1) noexcept;
947
949 template <RunOn run_on AMREX_DEFAULT_RUNON>
950 BaseFab<T>& linInterp (const BaseFab<T>& f1, int comp1,
951 const BaseFab<T>& f2, int comp2,
952 Real t1, Real t2, Real t,
953 const Box& b, int comp, int numcomp = 1) noexcept;
954
964 template <RunOn run_on AMREX_DEFAULT_RUNON>
965 BaseFab<T>& linComb (const BaseFab<T>& f1, const Box& b1, int comp1,
966 const BaseFab<T>& f2, const Box& b2, int comp2,
967 Real alpha, Real beta, const Box& b,
968 int comp, int numcomp = 1) noexcept;
969
971 template <RunOn run_on AMREX_DEFAULT_RUNON>
972 [[nodiscard]] T dot (const Box& xbx, int xcomp, const BaseFab<T>& y, const Box& ybx, int ycomp,
973 int numcomp = 1) const noexcept;
974
975 template <RunOn run_on AMREX_DEFAULT_RUNON>
976 [[nodiscard]] T dotmask (const BaseFab<int>& mask, const Box& xbx, int xcomp,
977 const BaseFab<T>& y, const Box& ybx, int ycomp,
978 int numcomp) const noexcept;
979
981 void SetBoxType (const IndexType& typ) noexcept { this->domain.setType(typ); }
982
983 //
984 // New interfaces
985 //
986
988 template <RunOn run_on AMREX_DEFAULT_RUNON>
989 void setVal (T const& val) noexcept;
990 //
992 template <RunOn run_on AMREX_DEFAULT_RUNON>
993 void setVal (T const& x, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept;
994
995 template <RunOn run_on AMREX_DEFAULT_RUNON>
996 void setValIf (T const& val, const BaseFab<int>& mask) noexcept;
997 //
999 template <RunOn run_on AMREX_DEFAULT_RUNON>
1000 void setValIf (T const& val, Box const& bx, const BaseFab<int>& mask, DestComp dcomp, NumComps ncomp) noexcept;
1001
1002 template <RunOn run_on AMREX_DEFAULT_RUNON>
1003 void setValIfNot (T const& val, const BaseFab<int>& mask) noexcept;
1004 //
1006 template <RunOn run_on AMREX_DEFAULT_RUNON>
1007 void setValIfNot (T const& val, Box const& bx, const BaseFab<int>& mask, DestComp dcomp, NumComps ncomp) noexcept;
1008
1010 template <RunOn run_on AMREX_DEFAULT_RUNON>
1011 void setComplement (T const& x, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept;
1012
1018 template <RunOn run_on AMREX_DEFAULT_RUNON>
1019 BaseFab<T>& copy (const BaseFab<T>& src) noexcept;
1020 //
1022 template <RunOn run_on AMREX_DEFAULT_RUNON>
1023 BaseFab<T>& copy (const BaseFab<T>& src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept;
1024
1026 template <RunOn run_on AMREX_DEFAULT_RUNON>
1027 BaseFab<T>& plus (T const& val) noexcept;
1028 //
1029 template <RunOn run_on AMREX_DEFAULT_RUNON>
1030 BaseFab<T>& operator+= (T const& val) noexcept;
1031 //
1033 template <RunOn run_on AMREX_DEFAULT_RUNON>
1034 BaseFab<T>& plus (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept;
1040 template <RunOn run_on AMREX_DEFAULT_RUNON>
1041 BaseFab<T>& plus (const BaseFab<T>& src) noexcept;
1042 //
1043 template <RunOn run_on AMREX_DEFAULT_RUNON>
1044 BaseFab<T>& operator+= (const BaseFab<T>& src) noexcept;
1045 //
1047 template <RunOn run_on AMREX_DEFAULT_RUNON>
1048 BaseFab<T>& plus (const BaseFab<T>& src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept;
1049
1051 template <RunOn run_on AMREX_DEFAULT_RUNON>
1052 BaseFab<T>& minus (T const& val) noexcept;
1053 //
1054 template <RunOn run_on AMREX_DEFAULT_RUNON>
1055 BaseFab<T>& operator-= (T const& val) noexcept;
1056 //
1058 template <RunOn run_on AMREX_DEFAULT_RUNON>
1059 BaseFab<T>& minus (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept;
1065 template <RunOn run_on AMREX_DEFAULT_RUNON>
1066 BaseFab<T>& minus (const BaseFab<T>& src) noexcept;
1067 //
1068 template <RunOn run_on AMREX_DEFAULT_RUNON>
1069 BaseFab<T>& operator-= (const BaseFab<T>& src) noexcept;
1070 //
1072 template <RunOn run_on AMREX_DEFAULT_RUNON>
1073 BaseFab<T>& minus (const BaseFab<T>& src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept;
1074
1076 template <RunOn run_on AMREX_DEFAULT_RUNON>
1077 BaseFab<T>& mult (T const& val) noexcept;
1078 //
1079 template <RunOn run_on AMREX_DEFAULT_RUNON>
1080 BaseFab<T>& operator*= (T const& val) noexcept;
1081 //
1083 template <RunOn run_on AMREX_DEFAULT_RUNON>
1084 BaseFab<T>& mult (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept;
1090 template <RunOn run_on AMREX_DEFAULT_RUNON>
1091 BaseFab<T>& mult (const BaseFab<T>& src) noexcept;
1092 //
1093 template <RunOn run_on AMREX_DEFAULT_RUNON>
1094 BaseFab<T>& operator*= (const BaseFab<T>& src) noexcept;
1095 //
1097 template <RunOn run_on AMREX_DEFAULT_RUNON>
1098 BaseFab<T>& mult (const BaseFab<T>& src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept;
1099
1101 template <RunOn run_on AMREX_DEFAULT_RUNON>
1102 BaseFab<T>& divide (T const& val) noexcept;
1103 //
1104 template <RunOn run_on AMREX_DEFAULT_RUNON>
1105 BaseFab<T>& operator/= (T const& val) noexcept;
1106 //
1108 template <RunOn run_on AMREX_DEFAULT_RUNON>
1109 BaseFab<T>& divide (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept;
1115 template <RunOn run_on AMREX_DEFAULT_RUNON>
1116 BaseFab<T>& divide (const BaseFab<T>& src) noexcept;
1117 //
1118 template <RunOn run_on AMREX_DEFAULT_RUNON>
1119 BaseFab<T>& operator/= (const BaseFab<T>& src) noexcept;
1120 //
1122 template <RunOn run_on AMREX_DEFAULT_RUNON>
1123 BaseFab<T>& divide (const BaseFab<T>& src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept;
1124
1126 template <RunOn run_on AMREX_DEFAULT_RUNON>
1127 BaseFab<T>& negate () noexcept;
1128 //
1129 template <RunOn run_on AMREX_DEFAULT_RUNON>
1130 BaseFab<T>& negate (const Box& bx, DestComp dcomp, NumComps ncomp) noexcept;
1131
1133 template <RunOn run_on AMREX_DEFAULT_RUNON>
1134 BaseFab<T>& invert (T const& r) noexcept;
1135 //
1136 template <RunOn run_on AMREX_DEFAULT_RUNON>
1137 BaseFab<T>& invert (T const& r, const Box& bx, DestComp dcomp, NumComps ncomp) noexcept;
1138
1140 template <RunOn run_on AMREX_DEFAULT_RUNON>
1141 [[nodiscard]] T sum (const Box& bx, DestComp dcomp, NumComps ncomp) const noexcept;
1142
1144 template <RunOn run_on AMREX_DEFAULT_RUNON>
1145 [[nodiscard]] T dot (const BaseFab<T>& src, const Box& bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) const noexcept;
1146
1148 template <RunOn run_on AMREX_DEFAULT_RUNON>
1149 [[nodiscard]] T dot (const Box& bx, int destcomp, int numcomp) const noexcept;
1150
1152 template <RunOn run_on AMREX_DEFAULT_RUNON>
1153 [[nodiscard]] T dot (const Box& bx, DestComp dcomp, NumComps ncomp) const noexcept;
1154
1156 template <RunOn run_on AMREX_DEFAULT_RUNON>
1157 [[nodiscard]] T dotmask (const BaseFab<T>& src, const Box& bx, const BaseFab<int>& mask,
1158 SrcComp scomp, DestComp dcomp, NumComps ncomp) const noexcept;
1159
1160protected:
1162 void define ();
1163
1164 T* dptr = nullptr;
1166 int nvar = 0;
1168 bool ptr_owner = false;
1169 bool shared_memory = false;
1170#ifdef AMREX_USE_GPU
1172#endif
1173};
1174
1175template <class T>
1177T*
1178BaseFab<T>::dataPtr (const IntVect& p, int n) noexcept
1179{
1180 AMREX_ASSERT(n >= 0);
1181 AMREX_ASSERT(n < this->nvar);
1182 AMREX_ASSERT(!(this->dptr == nullptr));
1183 AMREX_ASSERT(this->domain.contains(p));
1184
1185 return this->dptr + (this->domain.index(p)+n*this->domain.numPts());
1186}
1187
1188template <class T>
1190const T*
1191BaseFab<T>::dataPtr (const IntVect& p, int n) const noexcept
1192{
1193 AMREX_ASSERT(n >= 0);
1194 AMREX_ASSERT(n < this->nvar);
1195 AMREX_ASSERT(!(this->dptr == nullptr));
1196 AMREX_ASSERT(this->domain.contains(p));
1197
1198 return this->dptr + (this->domain.index(p)+n*this->domain.numPts());
1199}
1200
1201template <class T>
1202void
1204{
1205#ifdef AMREX_USE_GPU
1206 if (this->arena()->isManaged()) {
1207#if defined(AMREX_USE_SYCL)
1208 // xxxxx SYCL todo: prefetchToHost
1209 // std::size_t s = sizeof(T)*this->nvar*this->domain.numPts();
1210 // auto& q = Gpu::Device::streamQueue();
1211 // q.submit([&] (sycl::handler& h) { h.prefetch(this->dptr, s); });
1212#elif defined(AMREX_USE_CUDA) && !defined(_WIN32)
1213 if (Gpu::Device::devicePropMajor() >= 6) {
1214 std::size_t s = sizeof(T)*this->nvar*this->domain.numPts();
1215#if defined(CUDART_VERSION) && (CUDART_VERSION >= 13000)
1216 cudaMemLocation location = {};
1217 location.type = cudaMemLocationTypeHost;
1218 AMREX_CUDA_SAFE_CALL(cudaMemPrefetchAsync(this->dptr, s, location, 0,
1219 Gpu::gpuStream()));
1220#else
1221 AMREX_CUDA_SAFE_CALL(cudaMemPrefetchAsync(this->dptr, s,
1222 cudaCpuDeviceId,
1223 Gpu::gpuStream()));
1224#endif
1225 }
1226#elif defined(AMREX_USE_HIP)
1227 // xxxxx HIP FIX HERE after managed memory is supported
1228#endif
1229 }
1230#endif
1231}
1232
1233template <class T>
1234void
1236{
1237#ifdef AMREX_USE_GPU
1238 if (this->arena()->isManaged()) {
1239#if defined(AMREX_USE_SYCL)
1240 std::size_t s = sizeof(T)*this->nvar*this->domain.numPts();
1241 auto& q = Gpu::Device::streamQueue();
1242 q.submit([&] (sycl::handler& h) { h.prefetch(this->dptr, s); });
1243#elif defined(AMREX_USE_CUDA) && !defined(_WIN32)
1244 if (Gpu::Device::devicePropMajor() >= 6) {
1245 std::size_t s = sizeof(T)*this->nvar*this->domain.numPts();
1246#if defined(CUDART_VERSION) && (CUDART_VERSION >= 13000)
1247 cudaMemLocation location = {};
1248 location.type = cudaMemLocationTypeDevice;
1249 location.id = Gpu::Device::deviceId();
1250 AMREX_CUDA_SAFE_CALL(cudaMemPrefetchAsync(this->dptr, s, location, 0,
1251 Gpu::gpuStream()));
1252#else
1253 AMREX_CUDA_SAFE_CALL(cudaMemPrefetchAsync(this->dptr, s,
1255 Gpu::gpuStream()));
1256#endif
1257 }
1258#elif defined(AMREX_USE_HIP)
1259 // xxxxx HIP FIX HERE after managed memory is supported
1260#endif
1261 }
1262#endif
1263}
1264
1265template <class T>
1267T&
1268BaseFab<T>::operator() (const IntVect& p, int n) noexcept
1269{
1270 AMREX_ASSERT(n >= 0);
1271 AMREX_ASSERT(n < this->nvar);
1272 AMREX_ASSERT(!(this->dptr == nullptr));
1273 AMREX_ASSERT(this->domain.contains(p));
1274
1275 return this->dptr[this->domain.index(p)+n*this->domain.numPts()];
1276}
1277
1278template <class T>
1280T&
1282{
1283 AMREX_ASSERT(!(this->dptr == nullptr));
1284 AMREX_ASSERT(this->domain.contains(p));
1285
1286 return this->dptr[this->domain.index(p)];
1287}
1288
1289template <class T>
1291const T&
1292BaseFab<T>::operator() (const IntVect& p, int n) const noexcept
1293{
1294 AMREX_ASSERT(n >= 0);
1295 AMREX_ASSERT(n < this->nvar);
1296 AMREX_ASSERT(!(this->dptr == nullptr));
1297 AMREX_ASSERT(this->domain.contains(p));
1298
1299 return this->dptr[this->domain.index(p)+n*this->domain.numPts()];
1300}
1301
1302template <class T>
1304const T&
1305BaseFab<T>::operator() (const IntVect& p) const noexcept
1306{
1307 AMREX_ASSERT(!(this->dptr == nullptr));
1308 AMREX_ASSERT(this->domain.contains(p));
1309
1310 return this->dptr[this->domain.index(p)];
1311}
1312
1313template <class T>
1314void
1316 const IntVect& pos,
1317 int n,
1318 int numcomp) const noexcept
1319{
1320 const int loc = this->domain.index(pos);
1321 const Long sz = this->domain.numPts();
1322
1323 AMREX_ASSERT(!(this->dptr == nullptr));
1324 AMREX_ASSERT(n >= 0 && n + numcomp <= this->nvar);
1325
1326 for (int k = 0; k < numcomp; k++) {
1327 data[k] = this->dptr[loc+(n+k)*sz];
1328 }
1329}
1330
1331template <class T>
1332void
1334 const IntVect& pos) const noexcept
1335{
1336 getVal(data,pos,0,this->nvar);
1337}
1338
1339template <class T>
1341BaseFab<T>::shift (const IntVect& v) noexcept
1342{
1343 this->domain += v;
1344 return *this;
1345}
1346
1347template <class T>
1349BaseFab<T>::shift (int idir, int n_cell) noexcept
1350{
1351 this->domain.shift(idir,n_cell);
1352 return *this;
1353}
1354
1355template <class T>
1356BaseFab<T> &
1358{
1359 this->domain.shiftHalf(v);
1360 return *this;
1361}
1362
1363template <class T>
1364BaseFab<T> &
1365BaseFab<T>::shiftHalf (int idir, int n_cell) noexcept
1366{
1367 this->domain.shiftHalf(idir,n_cell);
1368 return *this;
1369}
1370
1371template <class T>
1372template <RunOn run_on, class U,
1373 std::enable_if_t<std::is_same_v<U,float> || std::is_same_v<U,double>, int> FOO>
1374void
1376{
1377 amrex::fill_snan<run_on>(this->dptr, this->truesize);
1378}
1379
1380template <class T>
1381template <RunOn run_on>
1382void
1383BaseFab<T>::setVal (T const& x, const Box& bx, int n) noexcept
1384{
1385 this->setVal<run_on>(x, bx, DestComp{n}, NumComps{1});
1386}
1387
1388template <class T>
1389template <RunOn run_on>
1390void
1391BaseFab<T>::setVal (T const& x, int n) noexcept
1392{
1393 this->setVal<run_on>(x, this->domain, DestComp{n}, NumComps{1});
1394}
1395
1396template <class T>
1397template <RunOn run_on>
1398void
1399BaseFab<T>::setVal (T const& x, const Box& bx, int dcomp, int ncomp) noexcept
1400{
1401 this->setVal<run_on>(x, bx, DestComp{dcomp}, NumComps{ncomp});
1402}
1403
1404template <class T>
1405template <RunOn run_on>
1406void
1407BaseFab<T>::setValIfNot (T const& val, const Box& bx, const BaseFab<int>& mask, int ns, int num) noexcept
1408{
1409 this->setValIfNot<run_on>(val, bx, mask, DestComp{ns}, NumComps{num});
1410}
1411
1412template <class T>
1413template <RunOn run_on>
1415BaseFab<T>::copy (const BaseFab<T>& src, const Box& srcbox, int srccomp,
1416 const Box& destbox, int destcomp, int numcomp) noexcept
1417{
1418 AMREX_ASSERT(destbox.ok());
1419 AMREX_ASSERT(srcbox.sameSize(destbox));
1420 AMREX_ASSERT(src.box().contains(srcbox));
1421 AMREX_ASSERT(this->domain.contains(destbox));
1422 AMREX_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
1423 AMREX_ASSERT(destcomp >= 0 && destcomp+numcomp <= this->nvar);
1424
1425 Array4<T> const& d = this->array();
1426 Array4<T const> const& s = src.const_array();
1427 const auto dlo = amrex::lbound(destbox);
1428 const auto slo = amrex::lbound(srcbox);
1429 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
1430
1431 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
1432 {
1433 d(i,j,k,n+destcomp) = s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
1434 });
1435
1436 return *this;
1437}
1438
1439template <class T>
1440template <RunOn run_on>
1442BaseFab<T>::copy (const BaseFab<T>& src, const Box& destbox) noexcept
1443{
1444 return this->copy<run_on>(src, destbox, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
1445}
1446
1447template <class T>
1448template <RunOn run_on>
1450BaseFab<T>::copy (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
1451{
1452 return copy<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
1453}
1454
1455template <class T>
1456void
1458{
1459 AMREX_ASSERT(this->dptr == nullptr);
1460 AMREX_ASSERT(this->domain.numPts() > 0);
1461 AMREX_ASSERT(this->nvar >= 0);
1462 if (this->nvar == 0) { return; }
1463 AMREX_ASSERT(std::numeric_limits<Long>::max()/this->nvar > this->domain.numPts());
1464
1465 this->truesize = this->nvar*this->domain.numPts();
1466 this->ptr_owner = true;
1467 this->dptr = static_cast<T*>(this->alloc(this->truesize*sizeof(T)));
1468#ifdef AMREX_USE_GPU
1469 this->alloc_stream = Gpu::gpuStream();
1470#endif
1471
1472 placementNew(this->dptr, this->truesize);
1473
1474 amrex::update_fab_stats(this->domain.numPts(), this->truesize, sizeof(T));
1475
1476 if constexpr (std::is_same_v<T,float> || std::is_same_v<T,double>) {
1477 if (amrex::InitSNaN() && this->truesize > 0) {
1478#ifdef AMREX_USE_GPU
1479 if (Gpu::inLaunchRegion() && arena()->isDeviceAccessible()) {
1480 this->template fill_snan<RunOn::Device>();
1482 } else
1483#endif
1484 {
1485 this->template fill_snan<RunOn::Host>();
1486 }
1487 }
1488 }
1489}
1490
1491template <class T>
1493 : DataAllocator{ar}
1494{}
1495
1496template <class T>
1497BaseFab<T>::BaseFab (const Box& bx, int n, Arena* ar)
1498 : DataAllocator{ar}, domain(bx), nvar(n)
1499{
1500 define();
1501}
1502
1503template <class T>
1504BaseFab<T>::BaseFab (const Box& bx, int n, bool alloc, bool shared, Arena* ar)
1505 : DataAllocator{ar}, domain(bx), nvar(n), shared_memory(shared)
1506{
1507 if (!this->shared_memory && alloc) { define(); }
1508}
1509
1510template <class T>
1511BaseFab<T>::BaseFab (const BaseFab<T>& rhs, MakeType make_type, int scomp, int ncomp)
1512 : DataAllocator{rhs.arena()},
1513 dptr(const_cast<T*>(rhs.dataPtr(scomp))),
1514 domain(rhs.domain), nvar(ncomp),
1515 truesize(ncomp*rhs.domain.numPts())
1516{
1517 AMREX_ASSERT(scomp+ncomp <= rhs.nComp());
1518 if (make_type == amrex::make_deep_copy)
1519 {
1520 this->dptr = nullptr;
1521 define();
1522 this->copy<RunOn::Device>(rhs, this->domain, scomp, this->domain, 0, ncomp);
1523 } else if (make_type == amrex::make_alias) {
1524 ; // nothing to do
1525 } else {
1526 amrex::Abort("BaseFab: unknown MakeType");
1527 }
1528}
1529
1530template<class T>
1531BaseFab<T>::BaseFab (const Box& bx, int ncomp, T* p)
1532 : dptr(p), domain(bx), nvar(ncomp), truesize(bx.numPts()*ncomp)
1533{
1534}
1535
1536template<class T>
1537BaseFab<T>::BaseFab (const Box& bx, int ncomp, T const* p)
1538 : dptr(const_cast<T*>(p)), domain(bx), nvar(ncomp), truesize(bx.numPts()*ncomp)
1539{
1540}
1541
1542template<class T>
1544 : dptr(a.p),
1545 domain(IntVect(AMREX_D_DECL(a.begin[0],a.begin[1],a.begin[2])),
1546 IntVect(AMREX_D_DECL(a.end[0]-1,a.end[1]-1,a.end[2]-1))),
1547 nvar(a.nComp()), truesize(a.size())
1548{}
1549
1550template<class T>
1552 : dptr(a.p),
1553 domain(IntVect(AMREX_D_DECL(a.begin[0],a.begin[1],a.begin[2])),
1554 IntVect(AMREX_D_DECL(a.end[0]-1,a.end[1]-1,a.end[2]-1)), t),
1555 nvar(a.nComp()), truesize(a.size())
1556{}
1557
1558template<class T>
1560 : dptr(const_cast<T*>(a.p)),
1561 domain(IntVect(AMREX_D_DECL(a.begin[0],a.begin[1],a.begin[2])),
1562 IntVect(AMREX_D_DECL(a.end[0]-1,a.end[1]-1,a.end[2]-1))),
1563 nvar(a.nComp()), truesize(a.size())
1564{}
1565
1566template<class T>
1568 : dptr(const_cast<T*>(a.p)),
1569 domain(IntVect(AMREX_D_DECL(a.begin[0],a.begin[1],a.begin[2])),
1570 IntVect(AMREX_D_DECL(a.end[0]-1,a.end[1]-1,a.end[2]-1)), t),
1571 nvar(a.nComp()), truesize(a.size())
1572{}
1573
1574template <class T>
1576{
1577 clear();
1578}
1579
1580template <class T>
1582 : DataAllocator{rhs.arena()},
1583 dptr(rhs.dptr), domain(rhs.domain),
1584 nvar(rhs.nvar), truesize(rhs.truesize),
1585 ptr_owner(rhs.ptr_owner), shared_memory(rhs.shared_memory)
1586#ifdef AMREX_USE_GPU
1587 , alloc_stream(rhs.alloc_stream)
1588#endif
1589{
1590 rhs.dptr = nullptr;
1591 rhs.ptr_owner = false;
1592}
1593
1594template <class T>
1595BaseFab<T>&
1597{
1598 if (this != &rhs) {
1599 clear();
1600 DataAllocator::operator=(rhs);
1601 dptr = rhs.dptr;
1602 domain = rhs.domain;
1603 nvar = rhs.nvar;
1604 truesize = rhs.truesize;
1605 ptr_owner = rhs.ptr_owner;
1606 shared_memory = rhs.shared_memory;
1607#ifdef AMREX_USE_GPU
1608 alloc_stream = rhs.alloc_stream;
1609#endif
1610
1611 rhs.dptr = nullptr;
1612 rhs.ptr_owner = false;
1613 }
1614 return *this;
1615}
1616
1617template <class T>
1618template <RunOn run_on>
1620BaseFab<T>::operator= (T const& t) noexcept
1621{
1622 setVal<run_on>(t);
1623 return *this;
1624}
1625
1626template <class T>
1627void
1628BaseFab<T>::resize (const Box& b, int n, Arena* ar)
1629{
1630 this->nvar = n;
1631 this->domain = b;
1632
1633 if (ar == nullptr) {
1634 ar = m_arena;
1635 }
1636
1637 if (arena() != DataAllocator(ar).arena()) {
1638 clear();
1639 m_arena = ar;
1640 define();
1641 }
1642 else if (this->dptr == nullptr || !this->ptr_owner)
1643 {
1644 if (this->shared_memory) {
1645 amrex::Abort("BaseFab::resize: BaseFab in shared memory cannot increase size");
1646 }
1647
1648 this->dptr = nullptr;
1649 define();
1650 }
1651 else if (this->nvar*this->domain.numPts() > this->truesize
1652#ifdef AMREX_USE_GPU
1653 || (arena()->isStreamOrderedArena() && alloc_stream != Gpu::gpuStream())
1654#endif
1655 )
1656 {
1657 if (this->shared_memory) {
1658 amrex::Abort("BaseFab::resize: BaseFab in shared memory cannot increase size");
1659 }
1660
1661 clear();
1662
1663 define();
1664 }
1665}
1666
1667template <class T>
1668template <class U, std::enable_if_t<std::is_trivially_destructible_v<U>,int>>
1669Elixir
1671{
1672 bool o;
1673 if (Gpu::inLaunchRegion()) {
1674 o = this->ptr_owner;
1675 this->ptr_owner = false;
1676 if (o && this->dptr) {
1677 if (this->nvar > 1) {
1678 amrex::update_fab_stats(-this->truesize/this->nvar, -this->truesize, sizeof(T));
1679 } else {
1680 amrex::update_fab_stats(0, -this->truesize, sizeof(T));
1681 }
1682 }
1683 } else {
1684 o = false;
1685 }
1686 return Elixir((o ? this->dptr : nullptr), this->arena());
1687}
1688
1689template <class T>
1690void
1692{
1693 if (this->dptr)
1694 {
1695 //
1696 // Call T::~T() on the to-be-destroyed memory.
1697 //
1698 if (this->ptr_owner)
1699 {
1700 if (this->shared_memory)
1701 {
1702 amrex::Abort("BaseFab::clear: BaseFab cannot be owner of shared memory");
1703 }
1704
1705 placementDelete(this->dptr, this->truesize);
1706
1707#ifdef AMREX_USE_GPU
1708 this->arena()->streamOrderedFree(this->dptr, alloc_stream);
1709#else
1710 this->free(this->dptr);
1711#endif
1712
1713 if (this->nvar > 1) {
1714 amrex::update_fab_stats(-this->truesize/this->nvar, -this->truesize, sizeof(T));
1715 } else {
1716 amrex::update_fab_stats(0, -this->truesize, sizeof(T));
1717 }
1718 }
1719
1720 this->dptr = nullptr;
1721 this->truesize = 0;
1722 }
1723}
1724
1725template <class T>
1726std::unique_ptr<T,DataDeleter>
1728{
1729 std::unique_ptr<T,DataDeleter> r(nullptr, DataDeleter{this->arena()});
1730 if (this->dptr && this->ptr_owner) {
1731 r.reset(this->dptr);
1732 this->ptr_owner = false;
1733 if (this->nvar > 1) {
1734 amrex::update_fab_stats(-this->truesize/this->nvar, -this->truesize, sizeof(T));
1735 } else {
1736 amrex::update_fab_stats(0, -this->truesize, sizeof(T));
1737 }
1738 }
1739 return r;
1740}
1741
1742template <class T>
1743template <RunOn run_on>
1744std::size_t
1746 int srccomp,
1747 int numcomp,
1748 void* dst) const noexcept
1749{
1750 BL_ASSERT(box().contains(srcbox));
1751 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= nComp());
1752
1753 if (srcbox.ok())
1754 {
1755 Array4<T> d(static_cast<T*>(dst),amrex::begin(srcbox),amrex::end(srcbox),numcomp);
1756 Array4<T const> const& s = this->const_array();
1757 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, srcbox, numcomp, i, j, k, n,
1758 {
1759 d(i,j,k,n) = s(i,j,k,n+srccomp);
1760 });
1761 return sizeof(T)*d.size();
1762 }
1763 else
1764 {
1765 return 0;
1766 }
1767}
1768
1769template <class T>
1770template <RunOn run_on, typename BUF>
1771std::size_t
1773 int dstcomp,
1774 int numcomp,
1775 const void* src) noexcept
1776{
1777 BL_ASSERT(box().contains(dstbox));
1778 BL_ASSERT(dstcomp >= 0 && dstcomp+numcomp <= nComp());
1779
1780 if (dstbox.ok())
1781 {
1782 Array4<BUF const> s(static_cast<BUF const*>(src), amrex::begin(dstbox),
1783 amrex::end(dstbox), numcomp);
1784 Array4<T> const& d = this->array();
1785 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, dstbox, numcomp, i, j, k, n,
1786 {
1787 d(i,j,k,n+dstcomp) = static_cast<T>(s(i,j,k,n));
1788 });
1789 return sizeof(BUF)*s.size();
1790 }
1791 else
1792 {
1793 return 0;
1794 }
1795}
1796
1797template <class T>
1798template <RunOn run_on, typename BUF>
1799std::size_t
1801 int dstcomp,
1802 int numcomp,
1803 const void* src) noexcept
1804{
1805 BL_ASSERT(box().contains(dstbox));
1806 BL_ASSERT(dstcomp >= 0 && dstcomp+numcomp <= nComp());
1807
1808 if (dstbox.ok())
1809 {
1810 Array4<BUF const> s(static_cast<BUF const*>(src), amrex::begin(dstbox),
1811 amrex::end(dstbox), numcomp);
1812 Array4<T> const& d = this->array();
1813 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, dstbox, numcomp, i, j, k, n,
1814 {
1815 d(i,j,k,n+dstcomp) += static_cast<T>(s(i,j,k,n));
1816 });
1817 return sizeof(BUF)*s.size();
1818 }
1819 else
1820 {
1821 return 0;
1822 }
1823}
1824
1825template <class T>
1826template <RunOn run_on>
1827void
1828BaseFab<T>::setComplement (T const& x, const Box& b, int ns, int num) noexcept
1829{
1830 this->setComplement<run_on>(x, b, DestComp{ns}, NumComps{num});
1831}
1832
1833template <class T>
1834template <RunOn run_on>
1835void
1837{
1838 this->abs<run_on>(this->domain,0,this->nvar);
1839}
1840
1841template <class T>
1842template <RunOn run_on>
1843void
1844BaseFab<T>::abs (int comp, int numcomp) noexcept
1845{
1846 this->abs<run_on>(this->domain,comp,numcomp);
1847}
1848
1849template <class T>
1850template <RunOn run_on>
1851void
1852BaseFab<T>::abs (const Box& subbox, int comp, int numcomp) noexcept
1853{
1854 Array4<T> const& a = this->array();
1855 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG (run_on, subbox, numcomp, i, j, k, n,
1856 {
1857 a(i,j,k,n+comp) = std::abs(a(i,j,k,n+comp));
1858 });
1859}
1860
1861template <class T>
1862template <RunOn run_on>
1863Real
1865 int scomp, int ncomp) const noexcept
1866{
1867 BL_ASSERT(this->domain.contains(subbox));
1868 BL_ASSERT(scomp >= 0 && scomp + ncomp <= this->nvar);
1869
1870 Array4<T const> const& a = this->const_array();
1871 Array4<int const> const& m = mask.const_array();
1872 Real r = 0.0;
1873#ifdef AMREX_USE_GPU
1874 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
1875 ReduceOps<ReduceOpMax> reduce_op;
1876 ReduceData<Real> reduce_data(reduce_op);
1877 using ReduceTuple = ReduceData<Real>::Type;
1878 reduce_op.eval(subbox, reduce_data,
1879 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
1880 {
1881 Real t = 0.0;
1882 if (m(i,j,k)) {
1883 for (int n = 0; n < ncomp; ++n) {
1884 t = amrex::max(t,static_cast<Real>(std::abs(a(i,j,k,n+scomp))));
1885 }
1886 }
1887 return {t};
1888 });
1889 ReduceTuple hv = reduce_data.value(reduce_op);
1890 r = amrex::get<0>(hv);
1891 } else
1892#endif
1893 {
1894 amrex::LoopOnCpu(subbox, ncomp, [=,&r] (int i, int j, int k, int n) noexcept
1895 {
1896 if (m(i,j,k)) {
1897 Real t = static_cast<Real>(std::abs(a(i,j,k,n+scomp)));
1898 r = amrex::max(r,t);
1899 }
1900 });
1901 }
1902 return r;
1903}
1904
1905template <class T>
1906template <RunOn run_on>
1907Real
1908BaseFab<T>::norm (int p, int comp, int numcomp) const
1909{
1910 return norm<run_on>(this->domain,p,comp,numcomp);
1911}
1912
1913template <class T>
1914template <RunOn run_on>
1915Real
1916BaseFab<T>::norm (const Box& subbox, int p, int comp, int numcomp) const
1917{
1918 BL_ASSERT(this->domain.contains(subbox));
1919 BL_ASSERT(comp >= 0 && comp + numcomp <= this->nvar);
1920
1921 Array4<T const> const& a = this->const_array();
1922 Real nrm = 0.;
1923#ifdef AMREX_USE_GPU
1924 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
1925 if (p == 0) {
1926 ReduceOps<ReduceOpMax> reduce_op;
1927 ReduceData<Real> reduce_data(reduce_op);
1928 using ReduceTuple = ReduceData<Real>::Type;
1929 reduce_op.eval(subbox, reduce_data,
1930 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
1931 {
1932 Real t = 0.0;
1933 for (int n = 0; n < numcomp; ++n) {
1934 t = amrex::max(t, static_cast<Real>(std::abs(a(i,j,k,n+comp))));
1935 }
1936 return {t};
1937 });
1938 ReduceTuple hv = reduce_data.value(reduce_op);
1939 nrm = amrex::get<0>(hv);
1940 } else if (p == 1) {
1941 ReduceOps<ReduceOpSum> reduce_op;
1942 ReduceData<Real> reduce_data(reduce_op);
1943 using ReduceTuple = ReduceData<Real>::Type;
1944 reduce_op.eval(subbox, reduce_data,
1945 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
1946 {
1947 Real t = 0.0;
1948 for (int n = 0; n < numcomp; ++n) {
1949 t += static_cast<Real>(std::abs(a(i,j,k,n+comp)));
1950 }
1951 return {t};
1952 });
1953 ReduceTuple hv = reduce_data.value(reduce_op);
1954 nrm = amrex::get<0>(hv);
1955 } else if (p == 2) {
1956 ReduceOps<ReduceOpSum> reduce_op;
1957 ReduceData<Real> reduce_data(reduce_op);
1958 using ReduceTuple = ReduceData<Real>::Type;
1959 reduce_op.eval(subbox, reduce_data,
1960 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
1961 {
1962 Real t = 0.0;
1963 for (int n = 0; n < numcomp; ++n) {
1964 t += static_cast<Real>(a(i,j,k,n+comp)*a(i,j,k,n+comp));
1965 }
1966 return {t};
1967 });
1968 ReduceTuple hv = reduce_data.value(reduce_op);
1969 nrm = amrex::get<0>(hv);
1970 } else {
1971 amrex::Error("BaseFab<T>::norm: wrong p");
1972 }
1973 } else
1974#endif
1975 {
1976 if (p == 0) {
1977 amrex::LoopOnCpu(subbox, numcomp, [=,&nrm] (int i, int j, int k, int n) noexcept
1978 {
1979 Real t = static_cast<Real>(std::abs(a(i,j,k,n+comp)));
1980 nrm = amrex::max(nrm,t);
1981 });
1982 } else if (p == 1) {
1983 amrex::LoopOnCpu(subbox, numcomp, [=,&nrm] (int i, int j, int k, int n) noexcept
1984 {
1985 nrm += std::abs(a(i,j,k,n+comp));
1986 });
1987 } else if (p == 2) {
1988 amrex::LoopOnCpu(subbox, numcomp, [=,&nrm] (int i, int j, int k, int n) noexcept
1989 {
1990 nrm += a(i,j,k,n+comp)*a(i,j,k,n+comp);
1991 });
1992 } else {
1993 amrex::Error("BaseFab<T>::norm: wrong p");
1994 }
1995 }
1996
1997 return nrm;
1998}
1999
2000template <class T>
2001template <RunOn run_on>
2002T
2003BaseFab<T>::min (int comp) const noexcept
2004{
2005 return this->min<run_on>(this->domain,comp);
2006}
2007
2008template <class T>
2009template <RunOn run_on>
2010T
2011BaseFab<T>::min (const Box& subbox, int comp) const noexcept
2012{
2013 Array4<T const> const& a = this->const_array(comp);
2014#ifdef AMREX_USE_GPU
2015 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2016 ReduceOps<ReduceOpMin> reduce_op;
2017 ReduceData<T> reduce_data(reduce_op);
2018 using ReduceTuple = typename decltype(reduce_data)::Type;
2019 reduce_op.eval(subbox, reduce_data,
2020 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2021 {
2022 return { a(i,j,k) };
2023 });
2024 ReduceTuple hv = reduce_data.value(reduce_op);
2025 return amrex::get<0>(hv);
2026 } else
2027#endif
2028 {
2029 T r = std::numeric_limits<T>::max();
2030 amrex::LoopOnCpu(subbox, [=,&r] (int i, int j, int k) noexcept
2031 {
2032 r = amrex::min(r, a(i,j,k));
2033 });
2034 return r;
2035 }
2036}
2037
2038template <class T>
2039template <RunOn run_on>
2040T
2041BaseFab<T>::max (int comp) const noexcept
2042{
2043 return this->max<run_on>(this->domain,comp);
2044}
2045
2046template <class T>
2047template <RunOn run_on>
2048T
2049BaseFab<T>::max (const Box& subbox, int comp) const noexcept
2050{
2051 Array4<T const> const& a = this->const_array(comp);
2052#ifdef AMREX_USE_GPU
2053 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2054 ReduceOps<ReduceOpMax> reduce_op;
2055 ReduceData<T> reduce_data(reduce_op);
2056 using ReduceTuple = typename decltype(reduce_data)::Type;
2057 reduce_op.eval(subbox, reduce_data,
2058 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2059 {
2060 return { a(i,j,k) };
2061 });
2062 ReduceTuple hv = reduce_data.value(reduce_op);
2063 return amrex::get<0>(hv);
2064 } else
2065#endif
2066 {
2067 T r = std::numeric_limits<T>::lowest();
2068 amrex::LoopOnCpu(subbox, [=,&r] (int i, int j, int k) noexcept
2069 {
2070 r = amrex::max(r, a(i,j,k));
2071 });
2072 return r;
2073 }
2074}
2075
2076template <class T>
2077template <RunOn run_on>
2078std::pair<T,T>
2079BaseFab<T>::minmax (int comp) const noexcept
2080{
2081 return this->minmax<run_on>(this->domain,comp);
2082}
2083
2084template <class T>
2085template <RunOn run_on>
2086std::pair<T,T>
2087BaseFab<T>::minmax (const Box& subbox, int comp) const noexcept
2088{
2089 Array4<T const> const& a = this->const_array(comp);
2090#ifdef AMREX_USE_GPU
2091 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2093 ReduceData<T,T> reduce_data(reduce_op);
2094 using ReduceTuple = typename decltype(reduce_data)::Type;
2095 reduce_op.eval(subbox, reduce_data,
2096 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2097 {
2098 auto const x = a(i,j,k);
2099 return { x, x };
2100 });
2101 ReduceTuple hv = reduce_data.value(reduce_op);
2102 return std::make_pair(amrex::get<0>(hv), amrex::get<1>(hv));
2103 } else
2104#endif
2105 {
2106 T rmax = std::numeric_limits<T>::lowest();
2107 T rmin = std::numeric_limits<T>::max();
2108 amrex::LoopOnCpu(subbox, [=,&rmin,&rmax] (int i, int j, int k) noexcept
2109 {
2110 auto const x = a(i,j,k);
2111 rmin = amrex::min(rmin, x);
2112 rmax = amrex::max(rmax, x);
2113 });
2114 return std::make_pair(rmin,rmax);
2115 }
2116}
2117
2118template <class T>
2119template <RunOn run_on>
2120T
2121BaseFab<T>::maxabs (int comp) const noexcept
2122{
2123 return this->maxabs<run_on>(this->domain,comp);
2124}
2125
2126template <class T>
2127template <RunOn run_on>
2128T
2129BaseFab<T>::maxabs (const Box& subbox, int comp) const noexcept
2130{
2131 Array4<T const> const& a = this->const_array(comp);
2132#ifdef AMREX_USE_GPU
2133 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2134 ReduceOps<ReduceOpMax> reduce_op;
2135 ReduceData<T> reduce_data(reduce_op);
2136 using ReduceTuple = typename decltype(reduce_data)::Type;
2137 reduce_op.eval(subbox, reduce_data,
2138 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2139 {
2140 return { std::abs(a(i,j,k)) };
2141 });
2142 ReduceTuple hv = reduce_data.value(reduce_op);
2143 return amrex::get<0>(hv);
2144 } else
2145#endif
2146 {
2147 T r = 0;
2148 amrex::LoopOnCpu(subbox, [=,&r] (int i, int j, int k) noexcept
2149 {
2150 r = amrex::max(r, std::abs(a(i,j,k)));
2151 });
2152 return r;
2153 }
2154}
2155
2156
2157template <class T>
2158template <RunOn run_on>
2159IntVect
2160BaseFab<T>::indexFromValue (Box const& subbox, int comp, T const& value) const noexcept
2161{
2162 Array4<T const> const& a = this->const_array(comp);
2163#ifdef AMREX_USE_GPU
2164 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2165 Array<int,1+AMREX_SPACEDIM> ha{0,AMREX_D_DECL(std::numeric_limits<int>::lowest(),
2166 std::numeric_limits<int>::lowest(),
2167 std::numeric_limits<int>::lowest())};
2168 Gpu::DeviceVector<int> dv(1+AMREX_SPACEDIM);
2169 int* p = dv.data();
2170 Gpu::htod_memcpy_async(p, ha.data(), sizeof(int)*ha.size());
2171 amrex::ParallelFor(subbox, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
2172 {
2173 int* flag = p;
2174 if ((*flag == 0) && (a(i,j,k) == value)) {
2175 if (Gpu::Atomic::Exch(flag,1) == 0) {
2176 AMREX_D_TERM(p[1] = i;,
2177 p[2] = j;,
2178 p[3] = k;);
2179 }
2180 }
2181 });
2182 Gpu::dtoh_memcpy_async(ha.data(), p, sizeof(int)*ha.size());
2184 return IntVect(AMREX_D_DECL(ha[1],ha[2],ha[3]));
2185 } else
2186#endif
2187 {
2188 AMREX_LOOP_3D(subbox, i, j, k,
2189 {
2190 if (a(i,j,k) == value) { return IntVect(AMREX_D_DECL(i,j,k)); }
2191 });
2192 return IntVect::TheMinVector();
2193 }
2194}
2195
2196template <class T>
2197template <RunOn run_on>
2198IntVect
2199BaseFab<T>::minIndex (int comp) const noexcept
2200{
2201 return this->minIndex<run_on>(this->domain,comp);
2202}
2203
2204template <class T>
2205template <RunOn run_on>
2206IntVect
2207BaseFab<T>::minIndex (const Box& subbox, int comp) const noexcept
2208{
2209 T min_val = this->min<run_on>(subbox, comp);
2210 return this->indexFromValue<run_on>(subbox, comp, min_val);
2211}
2212
2213template <class T>
2214template <RunOn run_on>
2215void
2216BaseFab<T>::minIndex (const Box& subbox, Real& min_val, IntVect& min_idx, int comp) const noexcept
2217{
2218 min_val = this->min<run_on>(subbox, comp);
2219 min_idx = this->indexFromValue<run_on>(subbox, comp, min_val);
2220}
2221
2222template <class T>
2223template <RunOn run_on>
2224IntVect
2225BaseFab<T>::maxIndex (int comp) const noexcept
2226{
2227 return this->maxIndex<run_on>(this->domain,comp);
2228}
2229
2230template <class T>
2231template <RunOn run_on>
2232IntVect
2233BaseFab<T>::maxIndex (const Box& subbox, int comp) const noexcept
2234{
2235 T max_val = this->max<run_on>(subbox, comp);
2236 return this->indexFromValue<run_on>(subbox, comp, max_val);
2237}
2238
2239template <class T>
2240template <RunOn run_on>
2241void
2242BaseFab<T>::maxIndex (const Box& subbox, Real& max_val, IntVect& max_idx, int comp) const noexcept
2243{
2244 max_val = this->max<run_on>(subbox, comp);
2245 max_idx = this->indexFromValue<run_on>(subbox, comp, max_val);
2246}
2247
2248template <class T>
2249template <RunOn run_on>
2250int
2251BaseFab<T>::maskLT (BaseFab<int>& mask, T const& val, int comp) const noexcept
2252{
2253 mask.resize(this->domain,1);
2254 int cnt = 0;
2255 Array4<int> const& m = mask.array();
2256 Array4<T const> const& a = this->const_array(comp);
2257#ifdef AMREX_USE_GPU
2258 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2259 ReduceOps<ReduceOpSum> reduce_op;
2260 ReduceData<int> reduce_data(reduce_op);
2261 using ReduceTuple = typename decltype(reduce_data)::Type;
2262 reduce_op.eval(this->domain, reduce_data,
2263 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2264 {
2265 int t;
2266 if (a(i,j,k) < val) {
2267 m(i,j,k) = 1;
2268 t = 1;
2269 } else {
2270 m(i,j,k) = 0;
2271 t = 0;
2272 }
2273 return {t};
2274 });
2275 ReduceTuple hv = reduce_data.value(reduce_op);
2276 cnt = amrex::get<0>(hv);
2277 } else
2278#endif
2279 {
2280 AMREX_LOOP_3D(this->domain, i, j, k,
2281 {
2282 if (a(i,j,k) < val) {
2283 m(i,j,k) = 1;
2284 ++cnt;
2285 } else {
2286 m(i,j,k) = 0;
2287 }
2288 });
2289 }
2290
2291 return cnt;
2292}
2293
2294template <class T>
2295template <RunOn run_on>
2296int
2297BaseFab<T>::maskLE (BaseFab<int>& mask, T const& val, int comp) const noexcept
2298{
2299 mask.resize(this->domain,1);
2300 int cnt = 0;
2301 Array4<int> const& m = mask.array();
2302 Array4<T const> const& a = this->const_array(comp);
2303#ifdef AMREX_USE_GPU
2304 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2305 ReduceOps<ReduceOpSum> reduce_op;
2306 ReduceData<int> reduce_data(reduce_op);
2307 using ReduceTuple = typename decltype(reduce_data)::Type;
2308 reduce_op.eval(this->domain, reduce_data,
2309 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2310 {
2311 int t;
2312 if (a(i,j,k) <= val) {
2313 m(i,j,k) = 1;
2314 t = 1;
2315 } else {
2316 m(i,j,k) = 0;
2317 t = 0;
2318 }
2319 return {t};
2320 });
2321 ReduceTuple hv = reduce_data.value(reduce_op);
2322 cnt = amrex::get<0>(hv);
2323 } else
2324#endif
2325 {
2326 AMREX_LOOP_3D(this->domain, i, j, k,
2327 {
2328 if (a(i,j,k) <= val) {
2329 m(i,j,k) = 1;
2330 ++cnt;
2331 } else {
2332 m(i,j,k) = 0;
2333 }
2334 });
2335 }
2336
2337 return cnt;
2338}
2339
2340template <class T>
2341template <RunOn run_on>
2342int
2343BaseFab<T>::maskEQ (BaseFab<int>& mask, T const& val, int comp) const noexcept
2344{
2345 mask.resize(this->domain,1);
2346 int cnt = 0;
2347 Array4<int> const& m = mask.array();
2348 Array4<T const> const& a = this->const_array(comp);
2349#ifdef AMREX_USE_GPU
2350 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2351 ReduceOps<ReduceOpSum> reduce_op;
2352 ReduceData<int> reduce_data(reduce_op);
2353 using ReduceTuple = typename decltype(reduce_data)::Type;
2354 reduce_op.eval(this->domain, reduce_data,
2355 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2356 {
2357 int t;
2358 if (a(i,j,k) == val) {
2359 m(i,j,k) = 1;
2360 t = 1;
2361 } else {
2362 m(i,j,k) = 0;
2363 t = 0;
2364 }
2365 return {t};
2366 });
2367 ReduceTuple hv = reduce_data.value(reduce_op);
2368 cnt = amrex::get<0>(hv);
2369 } else
2370#endif
2371 {
2372 AMREX_LOOP_3D(this->domain, i, j, k,
2373 {
2374 if (a(i,j,k) == val) {
2375 m(i,j,k) = 1;
2376 ++cnt;
2377 } else {
2378 m(i,j,k) = 0;
2379 }
2380 });
2381 }
2382
2383 return cnt;
2384}
2385
2386template <class T>
2387template <RunOn run_on>
2388int
2389BaseFab<T>::maskGT (BaseFab<int>& mask, T const& val, int comp) const noexcept
2390{
2391 mask.resize(this->domain,1);
2392 int cnt = 0;
2393 Array4<int> const& m = mask.array();
2394 Array4<T const> const& a = this->const_array(comp);
2395#ifdef AMREX_USE_GPU
2396 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2397 ReduceOps<ReduceOpSum> reduce_op;
2398 ReduceData<int> reduce_data(reduce_op);
2399 using ReduceTuple = typename decltype(reduce_data)::Type;
2400 reduce_op.eval(this->domain, reduce_data,
2401 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2402 {
2403 int t;
2404 if (a(i,j,k) > val) {
2405 m(i,j,k) = 1;
2406 t = 1;
2407 } else {
2408 m(i,j,k) = 0;
2409 t = 0;
2410 }
2411 return {t};
2412 });
2413 ReduceTuple hv = reduce_data.value(reduce_op);
2414 cnt = amrex::get<0>(hv);
2415 } else
2416#endif
2417 {
2418 AMREX_LOOP_3D(this->domain, i, j, k,
2419 {
2420 if (a(i,j,k) > val) {
2421 m(i,j,k) = 1;
2422 ++cnt;
2423 } else {
2424 m(i,j,k) = 0;
2425 }
2426 });
2427 }
2428
2429 return cnt;
2430}
2431
2432template <class T>
2433template <RunOn run_on>
2434int
2435BaseFab<T>::maskGE (BaseFab<int>& mask, T const& val, int comp) const noexcept
2436{
2437 mask.resize(this->domain,1);
2438 int cnt = 0;
2439 Array4<int> const& m = mask.array();
2440 Array4<T const> const& a = this->const_array(comp);
2441#ifdef AMREX_USE_GPU
2442 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2443 ReduceOps<ReduceOpSum> reduce_op;
2444 ReduceData<int> reduce_data(reduce_op);
2445 using ReduceTuple = typename decltype(reduce_data)::Type;
2446 reduce_op.eval(this->domain, reduce_data,
2447 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2448 {
2449 int t;
2450 if (a(i,j,k) >= val) {
2451 m(i,j,k) = 1;
2452 t = 1;
2453 } else {
2454 m(i,j,k) = 0;
2455 t = 0;
2456 }
2457 return {t};
2458 });
2459 ReduceTuple hv = reduce_data.value(reduce_op);
2460 cnt = amrex::get<0>(hv);
2461 } else
2462#endif
2463 {
2464 AMREX_LOOP_3D(this->domain, i, j, k,
2465 {
2466 if (a(i,j,k) >= val) {
2467 m(i,j,k) = 1;
2468 ++cnt;
2469 } else {
2470 m(i,j,k) = 0;
2471 }
2472 });
2473 }
2474
2475 return cnt;
2476}
2477
2478template <class T>
2479template <RunOn run_on>
2482{
2483 Box ovlp(this->domain);
2484 ovlp &= x.domain;
2485 return ovlp.ok() ? this->atomicAdd<run_on>(x,ovlp,ovlp,0,0,this->nvar) : *this;
2486}
2487
2488template <class T>
2489template <RunOn run_on>
2491BaseFab<T>::saxpy (T a, const BaseFab<T>& x, const Box& srcbox, const Box& destbox,
2492 int srccomp, int destcomp, int numcomp) noexcept
2493{
2494 BL_ASSERT(srcbox.ok());
2495 BL_ASSERT(x.box().contains(srcbox));
2496 BL_ASSERT(destbox.ok());
2497 BL_ASSERT(box().contains(destbox));
2498 BL_ASSERT(destbox.sameSize(srcbox));
2499 BL_ASSERT( srccomp >= 0 && srccomp+numcomp <= x.nComp());
2500 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2501
2502 Array4<T> const& d = this->array();
2503 Array4<T const> const& s = x.const_array();
2504 const auto dlo = amrex::lbound(destbox);
2505 const auto slo = amrex::lbound(srcbox);
2506 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
2507 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
2508 {
2509 d(i,j,k,n+destcomp) += a * s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
2510 });
2511
2512 return *this;
2513}
2514
2515template <class T>
2516template <RunOn run_on>
2518BaseFab<T>::saxpy (T a, const BaseFab<T>& x) noexcept
2519{
2520 Box ovlp(this->domain);
2521 ovlp &= x.domain;
2522 return ovlp.ok() ? saxpy<run_on>(a,x,ovlp,ovlp,0,0,this->nvar) : *this;
2523}
2524
2525template <class T>
2526template <RunOn run_on>
2529 const Box& srcbox, const Box& destbox,
2530 int srccomp, int destcomp, int numcomp) noexcept
2531{
2532 BL_ASSERT(srcbox.ok());
2533 BL_ASSERT(x.box().contains(srcbox));
2534 BL_ASSERT(destbox.ok());
2535 BL_ASSERT(box().contains(destbox));
2536 BL_ASSERT(destbox.sameSize(srcbox));
2537 BL_ASSERT( srccomp >= 0 && srccomp+numcomp <= x.nComp());
2538 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2539
2540 Array4<T> const& d = this->array();
2541 Array4<T const> const& s = x.const_array();
2542 const auto dlo = amrex::lbound(destbox);
2543 const auto slo = amrex::lbound(srcbox);
2544 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
2545 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
2546 {
2547 d(i,j,k,n+destcomp) = s(i+offset.x,j+offset.y,k+offset.z,n+srccomp) + a*d(i,j,k,n+destcomp);
2548 });
2549
2550 return *this;
2551}
2552
2553template <class T>
2554template <RunOn run_on>
2556BaseFab<T>::addproduct (const Box& destbox, int destcomp, int numcomp,
2557 const BaseFab<T>& src1, int comp1,
2558 const BaseFab<T>& src2, int comp2) noexcept
2559{
2560 BL_ASSERT(destbox.ok());
2561 BL_ASSERT(box().contains(destbox));
2562 BL_ASSERT( comp1 >= 0 && comp1+numcomp <= src1.nComp());
2563 BL_ASSERT( comp2 >= 0 && comp2+numcomp <= src2.nComp());
2564 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2565
2566 Array4<T> const& d = this->array();
2567 Array4<T const> const& s1 = src1.const_array();
2568 Array4<T const> const& s2 = src2.const_array();
2569 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
2570 {
2571 d(i,j,k,n+destcomp) += s1(i,j,k,n+comp1) * s2(i,j,k,n+comp2);
2572 });
2573
2574 return *this;
2575}
2576
2577template <class T>
2578template <RunOn run_on>
2580BaseFab<T>::linComb (const BaseFab<T>& f1, const Box& b1, int comp1,
2581 const BaseFab<T>& f2, const Box& b2, int comp2,
2582 Real alpha, Real beta, const Box& b,
2583 int comp, int numcomp) noexcept
2584{
2585 BL_ASSERT(b1.ok());
2586 BL_ASSERT(f1.box().contains(b1));
2587 BL_ASSERT(b2.ok());
2588 BL_ASSERT(f2.box().contains(b2));
2589 BL_ASSERT(b.ok());
2590 BL_ASSERT(box().contains(b));
2591 BL_ASSERT(b.sameSize(b1));
2592 BL_ASSERT(b.sameSize(b2));
2593 BL_ASSERT(comp1 >= 0 && comp1+numcomp <= f1.nComp());
2594 BL_ASSERT(comp2 >= 0 && comp2+numcomp <= f2.nComp());
2595 BL_ASSERT(comp >= 0 && comp +numcomp <= nComp());
2596
2597 Array4<T> const& d = this->array();
2598 Array4<T const> const& s1 = f1.const_array();
2599 Array4<T const> const& s2 = f2.const_array();
2600 const auto dlo = amrex::lbound(b);
2601 const auto slo1 = amrex::lbound(b1);
2602 const auto slo2 = amrex::lbound(b2);
2603 const Dim3 off1{.x = slo1.x-dlo.x, .y = slo1.y-dlo.y, .z = slo1.z-dlo.z};
2604 const Dim3 off2{.x = slo2.x-dlo.x, .y = slo2.y-dlo.y, .z = slo2.z-dlo.z};
2605
2606 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, b, numcomp, i, j, k, n,
2607 {
2608 d(i,j,k,n+comp) = alpha*s1(i+off1.x,j+off1.y,k+off1.z,n+comp1)
2609 + beta*s2(i+off2.x,j+off2.y,k+off2.z,n+comp2);
2610 });
2611 return *this;
2612}
2613
2614template <class T>
2615template <RunOn run_on>
2616T
2617BaseFab<T>::dot (const Box& xbx, int xcomp,
2618 const BaseFab<T>& y, const Box& ybx, int ycomp,
2619 int numcomp) const noexcept
2620{
2621 BL_ASSERT(xbx.ok());
2622 BL_ASSERT(box().contains(xbx));
2623 BL_ASSERT(y.box().contains(ybx));
2624 BL_ASSERT(xbx.sameSize(ybx));
2625 BL_ASSERT(xcomp >= 0 && xcomp+numcomp <= nComp());
2626 BL_ASSERT(ycomp >= 0 && ycomp+numcomp <= y.nComp());
2627
2628 T r = 0;
2629
2630 const auto xlo = amrex::lbound(xbx);
2631 const auto ylo = amrex::lbound(ybx);
2632 const Dim3 offset{.x = ylo.x-xlo.x, .y = ylo.y-xlo.y, .z = ylo.z-xlo.z};
2633 Array4<T const> const& xa = this->const_array();
2634 Array4<T const> const& ya = y.const_array();
2635
2636#ifdef AMREX_USE_GPU
2637 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2638 ReduceOps<ReduceOpSum> reduce_op;
2639 ReduceData<T> reduce_data(reduce_op);
2640 using ReduceTuple = typename decltype(reduce_data)::Type;
2641 reduce_op.eval(xbx, reduce_data,
2642 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2643 {
2644 T t = 0;
2645 for (int n = 0; n < numcomp; ++n) {
2646 t += xa(i,j,k,n+xcomp) * ya(i+offset.x,j+offset.y,k+offset.z,n+ycomp);
2647 }
2648 return {t};
2649 });
2650 ReduceTuple hv = reduce_data.value(reduce_op);
2651 r = amrex::get<0>(hv);
2652 } else
2653#endif
2654 {
2655 AMREX_LOOP_4D(xbx, numcomp, i, j, k, n,
2656 {
2657 r += xa(i,j,k,n+xcomp) * ya(i+offset.x,j+offset.y,k+offset.z,n+ycomp);
2658 });
2659 }
2660
2661 return r;
2662}
2663
2664template <class T>
2665template <RunOn run_on>
2666T
2667BaseFab<T>::dotmask (const BaseFab<int>& mask, const Box& xbx, int xcomp,
2668 const BaseFab<T>& y, const Box& ybx, int ycomp,
2669 int numcomp) const noexcept
2670{
2671 BL_ASSERT(xbx.ok());
2672 BL_ASSERT(box().contains(xbx));
2673 BL_ASSERT(y.box().contains(ybx));
2674 BL_ASSERT(xbx.sameSize(ybx));
2675 BL_ASSERT(xcomp >= 0 && xcomp+numcomp <= nComp());
2676 BL_ASSERT(ycomp >= 0 && ycomp+numcomp <= y.nComp());
2677
2678 T r = 0;
2679
2680 const auto xlo = amrex::lbound(xbx);
2681 const auto ylo = amrex::lbound(ybx);
2682 const Dim3 offset{.x = ylo.x-xlo.x, .y = ylo.y-xlo.y, .z = ylo.z-xlo.z};
2683
2684 Array4<T const> const& xa = this->const_array();
2685 Array4<T const> const& ya = y.const_array();
2686 Array4<int const> const& ma = mask.const_array();
2687
2688#ifdef AMREX_USE_GPU
2689 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2690 ReduceOps<ReduceOpSum> reduce_op;
2691 ReduceData<T> reduce_data(reduce_op);
2692 using ReduceTuple = typename decltype(reduce_data)::Type;
2693 reduce_op.eval(xbx, reduce_data,
2694 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2695 {
2696 int m = static_cast<int>(static_cast<bool>(ma(i,j,k)));
2697 T t = 0;
2698 for (int n = 0; n < numcomp; ++n) {
2699 t += xa(i,j,k,n+xcomp) * ya(i+offset.x,j+offset.y,k+offset.z,n+ycomp) * m;
2700 }
2701 return {t};
2702 });
2703 ReduceTuple hv = reduce_data.value(reduce_op);
2704 r = amrex::get<0>(hv);
2705 } else
2706#endif
2707 {
2708 AMREX_LOOP_4D(xbx, numcomp, i, j, k, n,
2709 {
2710 int m = static_cast<int>(static_cast<bool>(ma(i,j,k)));
2711 r += xa(i,j,k,n+xcomp) * ya(i+offset.x,j+offset.y,k+offset.z,n+ycomp) * m;
2712 });
2713 }
2714
2715 return r;
2716}
2717
2718template <class T>
2719template <RunOn run_on>
2720T
2721BaseFab<T>::sum (int comp, int numcomp) const noexcept
2722{
2723 return this->sum<run_on>(this->domain, DestComp{comp}, NumComps{numcomp});
2724}
2725
2726template <class T>
2727template <RunOn run_on>
2728T
2729BaseFab<T>::sum (const Box& subbox, int comp, int numcomp) const noexcept
2730{
2731 return this->sum<run_on>(subbox, DestComp{comp}, NumComps{numcomp});
2732}
2733
2734template <class T>
2735template <RunOn run_on>
2737BaseFab<T>::negate (int comp, int numcomp) noexcept
2738{
2739 return this->negate<run_on>(this->domain, DestComp{comp}, NumComps{numcomp});
2740}
2741
2742template <class T>
2743template <RunOn run_on>
2745BaseFab<T>::negate (const Box& b, int comp, int numcomp) noexcept
2746{
2747 return this->negate<run_on>(b, DestComp{comp}, NumComps{numcomp});
2748}
2749
2750template <class T>
2751template <RunOn run_on>
2753BaseFab<T>::invert (T const& r, int comp, int numcomp) noexcept
2754{
2755 return this->invert<run_on>(r, this->domain, DestComp{comp}, NumComps{numcomp});
2756}
2757
2758template <class T>
2759template <RunOn run_on>
2761BaseFab<T>::invert (T const& r, const Box& b, int comp, int numcomp) noexcept
2762{
2763 return this->invert<run_on>(r, b, DestComp{comp}, NumComps{numcomp});
2764}
2765
2766template <class T>
2767template <RunOn run_on>
2769BaseFab<T>::plus (T const& r, int comp, int numcomp) noexcept
2770{
2771 return this->plus<run_on>(r, this->domain, DestComp{comp}, NumComps{numcomp});
2772}
2773
2774template <class T>
2775template <RunOn run_on>
2777BaseFab<T>::plus (T const& r, const Box& b, int comp, int numcomp) noexcept
2778{
2779 return this->plus<run_on>(r, b, DestComp{comp}, NumComps{numcomp});
2780}
2781
2782template <class T>
2783template <RunOn run_on>
2785BaseFab<T>::plus (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
2786{
2787 return this->plus<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
2788}
2789
2790template <class T>
2791template <RunOn run_on>
2793BaseFab<T>::atomicAdd (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
2794{
2795 Box ovlp(this->domain);
2796 ovlp &= src.domain;
2797 return ovlp.ok() ? this->atomicAdd<run_on>(src,ovlp,ovlp,srccomp,destcomp,numcomp) : *this;
2798}
2799
2800template <class T>
2801template <RunOn run_on>
2803BaseFab<T>::plus (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
2804 int numcomp) noexcept
2805{
2806 return this->plus<run_on>(src, subbox, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
2807}
2808
2809template <class T>
2810template <RunOn run_on>
2812BaseFab<T>::atomicAdd (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
2813 int numcomp) noexcept
2814{
2815 Box ovlp(this->domain);
2816 ovlp &= src.domain;
2817 ovlp &= subbox;
2818 return ovlp.ok() ? this->atomicAdd<run_on>(src,ovlp,ovlp,srccomp,destcomp,numcomp) : *this;
2819}
2820
2821template <class T>
2822template <RunOn run_on>
2824BaseFab<T>::plus (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
2825 int srccomp, int destcomp, int numcomp) noexcept
2826{
2827 BL_ASSERT(destbox.ok());
2828 BL_ASSERT(src.box().contains(srcbox));
2829 BL_ASSERT(box().contains(destbox));
2830 BL_ASSERT(destbox.sameSize(srcbox));
2831 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
2832 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2833
2834 Array4<T> const& d = this->array();
2835 Array4<T const> const& s = src.const_array();
2836 const auto dlo = amrex::lbound(destbox);
2837 const auto slo = amrex::lbound(srcbox);
2838 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
2839 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
2840 {
2841 d(i,j,k,n+destcomp) += s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
2842 });
2843
2844 return *this;
2845}
2846
2847namespace detail
2848{
2849
2850template <RunOn run_on, typename T,
2851 std::enable_if_t<HasAtomicAdd<T>::value,int> FOO = 0>
2852void basefab_atomic_add (BaseFab<T>& dfab, const BaseFab<T>& sfab,
2853 const Box& srcbox, const Box& destbox,
2854 int srccomp, int destcomp, int numcomp) noexcept
2855{
2856 Array4<T> const& d = dfab.array();
2857 Array4<T const> const& s = sfab.const_array();
2858 const auto dlo = amrex::lbound(destbox);
2859 const auto slo = amrex::lbound(srcbox);
2860 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
2861 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
2862 {
2863 T* p = d.ptr(i,j,k,n+destcomp);
2864 HostDevice::Atomic::Add(p, s(i+offset.x,j+offset.y,k+offset.z,n+srccomp));
2865 });
2866}
2867
2868template <RunOn run_on, typename T,
2869 std::enable_if_t<! HasAtomicAdd<T>::value,int> FOO = 0>
2870void basefab_atomic_add (BaseFab<T>& dfab, const BaseFab<T>& sfab,
2871 const Box& srcbox, const Box& destbox,
2872 int srccomp, int destcomp, int numcomp)
2873{
2874 amrex::ignore_unused(dfab, sfab, srcbox, destbox, srccomp, destcomp, numcomp);
2875 amrex::Abort("BaseFab: atomicAdd not supported");
2876}
2877
2878}
2879
2880template <class T>
2881template <RunOn run_on>
2882BaseFab<T>&
2883BaseFab<T>::atomicAdd (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
2884 int srccomp, int destcomp, int numcomp) noexcept
2885{
2886 BL_ASSERT(destbox.ok());
2887 BL_ASSERT(src.box().contains(srcbox));
2888 BL_ASSERT(box().contains(destbox));
2889 BL_ASSERT(destbox.sameSize(srcbox));
2890 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
2891 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2892
2893 detail::basefab_atomic_add<run_on>(*this, src, srcbox, destbox,
2894 srccomp, destcomp, numcomp);
2895
2896 return *this;
2897}
2898
2899template <class T>
2900template <RunOn run_on>
2902BaseFab<T>::lockAdd (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
2903 int srccomp, int destcomp, int numcomp) noexcept
2904{
2905#if defined(AMREX_USE_OMP) && (AMREX_SPACEDIM > 1)
2906#if defined(AMREX_USE_GPU)
2907 if (run_on == RunOn::Host || Gpu::notInLaunchRegion()) {
2908#endif
2909 BL_ASSERT(destbox.ok());
2910 BL_ASSERT(src.box().contains(srcbox));
2911 BL_ASSERT(box().contains(destbox));
2912 BL_ASSERT(destbox.sameSize(srcbox));
2913 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
2914 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2915
2916 Array4<T> const& d = this->array();
2917 Array4<T const> const& s = src.const_array();
2918 auto const& dlo = amrex::lbound(destbox);
2919 auto const& dhi = amrex::ubound(destbox);
2920 auto const& len = amrex::length(destbox);
2921 auto const& slo = amrex::lbound(srcbox);
2922 Dim3 const offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
2923
2924 int planedim;
2925 int nplanes;
2926 int plo;
2927 if (len.z == 1) {
2928 planedim = 1;
2929 nplanes = len.y;
2930 plo = dlo.y;
2931 } else {
2932 planedim = 2;
2933 nplanes = len.z;
2934 plo = dlo.z;
2935 }
2936
2937 auto* mask = (bool*) amrex_mempool_alloc(sizeof(bool)*nplanes);
2938 for (int ip = 0; ip < nplanes; ++ip) {
2939 mask[ip] = false;
2940 }
2941
2942 int mm = 0;
2943 int planes_left = nplanes;
2944 while (planes_left > 0) {
2945 AMREX_ASSERT(mm < nplanes);
2946 auto const m = mm + plo;
2947 auto* lock = OpenMP::get_lock(m);
2948 if (omp_test_lock(lock))
2949 {
2950 auto lo = dlo;
2951 auto hi = dhi;
2952 if (planedim == 1) {
2953 lo.y = m;
2954 hi.y = m;
2955 } else {
2956 lo.z = m;
2957 hi.z = m;
2958 }
2959
2960 for (int n = 0; n < numcomp; ++n) {
2961 for (int k = lo.z; k <= hi.z; ++k) {
2962 for (int j = lo.y; j <= hi.y; ++j) {
2963 auto * pdst = d.ptr(dlo.x,j ,k ,n+destcomp);
2964 auto const* psrc = s.ptr(slo.x,j+offset.y,k+offset.z,n+ srccomp);
2965#pragma omp simd
2966 for (int ii = 0; ii < len.x; ++ii) {
2967 pdst[ii] += psrc[ii];
2968 }
2969 }
2970 }
2971 }
2972
2973 mask[mm] = true;
2974 --planes_left;
2975 omp_unset_lock(lock);
2976 if (planes_left == 0) { break; }
2977 }
2978
2979 ++mm;
2980 for (int ip = 0; ip < nplanes; ++ip) {
2981 int new_mm = (mm+ip) % nplanes;
2982 if ( ! mask[new_mm] ) {
2983 mm = new_mm;
2984 break;
2985 }
2986 }
2987 }
2988
2990
2991 return *this;
2992
2993#if defined(AMREX_USE_GPU)
2994 } else {
2995 return this->template atomicAdd<run_on>(src, srcbox, destbox, srccomp, destcomp, numcomp);
2996 }
2997#endif
2998#else
2999 return this->template atomicAdd<run_on>(src, srcbox, destbox, srccomp, destcomp, numcomp);
3000#endif
3001}
3002
3003template <class T>
3004template <RunOn run_on>
3006BaseFab<T>::minus (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
3007{
3008 return this->minus<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3009}
3010
3011template <class T>
3012template <RunOn run_on>
3014BaseFab<T>::minus (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp, int numcomp) noexcept
3015{
3016 return this->minus<run_on>(src, subbox, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3017}
3018
3019template <class T>
3020template <RunOn run_on>
3022BaseFab<T>::minus (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3023 int srccomp, int destcomp, int numcomp) noexcept
3024{
3025 BL_ASSERT(destbox.ok());
3026 BL_ASSERT(src.box().contains(srcbox));
3027 BL_ASSERT(box().contains(destbox));
3028 BL_ASSERT(destbox.sameSize(srcbox));
3029 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3030 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3031
3032 Array4<T> const& d = this->array();
3033 Array4<T const> const& s = src.const_array();
3034 const auto dlo = amrex::lbound(destbox);
3035 const auto slo = amrex::lbound(srcbox);
3036 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3037 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3038 {
3039 d(i,j,k,n+destcomp) -= s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
3040 });
3041
3042 return *this;
3043}
3044
3045template <class T>
3046template <RunOn run_on>
3048BaseFab<T>::mult (T const& r, int comp, int numcomp) noexcept
3049{
3050 return this->mult<run_on>(r, this->domain, DestComp{comp}, NumComps{numcomp});
3051}
3052
3053template <class T>
3054template <RunOn run_on>
3056BaseFab<T>::mult (T const& r, const Box& b, int comp, int numcomp) noexcept
3057{
3058 return this->mult<run_on>(r, b, DestComp{comp}, NumComps{numcomp});
3059}
3060
3061template <class T>
3062template <RunOn run_on>
3064BaseFab<T>::mult (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
3065{
3066 return this->mult<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3067}
3068
3069template <class T>
3070template <RunOn run_on>
3072BaseFab<T>::mult (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp, int numcomp) noexcept
3073{
3074 return this->mult<run_on>(src, subbox, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3075}
3076
3077template <class T>
3078template <RunOn run_on>
3080BaseFab<T>::mult (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3081 int srccomp, int destcomp, int numcomp) noexcept
3082{
3083 BL_ASSERT(destbox.ok());
3084 BL_ASSERT(src.box().contains(srcbox));
3085 BL_ASSERT(box().contains(destbox));
3086 BL_ASSERT(destbox.sameSize(srcbox));
3087 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3088 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3089
3090 Array4<T> const& d = this->array();
3091 Array4<T const> const& s = src.const_array();
3092 const auto dlo = amrex::lbound(destbox);
3093 const auto slo = amrex::lbound(srcbox);
3094 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3095 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3096 {
3097 d(i,j,k,n+destcomp) *= s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
3098 });
3099
3100 return *this;
3101}
3102
3103template <class T>
3104template <RunOn run_on>
3106BaseFab<T>::divide (T const& r, int comp, int numcomp) noexcept
3107{
3108 return this->divide<run_on>(r, this->domain, DestComp{comp}, NumComps{numcomp});
3109}
3110
3111template <class T>
3112template <RunOn run_on>
3114BaseFab<T>::divide (T const& r, const Box& b, int comp, int numcomp) noexcept
3115{
3116 return this->divide<run_on>(r, b, DestComp{comp}, NumComps{numcomp});
3117}
3118
3119template <class T>
3120template <RunOn run_on>
3122BaseFab<T>::divide (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
3123{
3124 return this->divide<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3125}
3126
3127template <class T>
3128template <RunOn run_on>
3130BaseFab<T>::divide (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp, int numcomp) noexcept
3131{
3132 return this->divide<run_on>(src, subbox, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3133}
3134
3135template <class T>
3136template <RunOn run_on>
3138BaseFab<T>::divide (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3139 int srccomp, int destcomp, int numcomp) noexcept
3140{
3141 BL_ASSERT(destbox.ok());
3142 BL_ASSERT(src.box().contains(srcbox));
3143 BL_ASSERT(box().contains(destbox));
3144 BL_ASSERT(destbox.sameSize(srcbox));
3145 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3146 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3147
3148 Array4<T> const& d = this->array();
3149 Array4<T const> const& s = src.const_array();
3150 const auto dlo = amrex::lbound(destbox);
3151 const auto slo = amrex::lbound(srcbox);
3152 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3153 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3154 {
3155 d(i,j,k,n+destcomp) /= s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
3156 });
3157
3158 return *this;
3159}
3160
3161template <class T>
3162template <RunOn run_on>
3165{
3166 Box ovlp(this->domain);
3167 ovlp &= src.domain;
3168 return ovlp.ok() ? this->protected_divide<run_on>(src,ovlp,ovlp,0,0,this->nvar) : *this;
3169}
3170
3171template <class T>
3172template <RunOn run_on>
3174BaseFab<T>::protected_divide (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
3175{
3176 Box ovlp(this->domain);
3177 ovlp &= src.domain;
3178 return ovlp.ok() ? this->protected_divide<run_on>(src,ovlp,ovlp,srccomp,destcomp,numcomp) : *this;
3179}
3180
3181template <class T>
3182template <RunOn run_on>
3184BaseFab<T>::protected_divide (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
3185 int numcomp) noexcept
3186{
3187 Box ovlp(this->domain);
3188 ovlp &= src.domain;
3189 ovlp &= subbox;
3190 return ovlp.ok() ? this->protected_divide<run_on>(src,ovlp,ovlp,srccomp,destcomp,numcomp) : *this;
3191}
3192
3193template <class T>
3194template <RunOn run_on>
3196BaseFab<T>::protected_divide (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3197 int srccomp, int destcomp, int numcomp) noexcept
3198{
3199 BL_ASSERT(destbox.ok());
3200 BL_ASSERT(src.box().contains(srcbox));
3201 BL_ASSERT(box().contains(destbox));
3202 BL_ASSERT(destbox.sameSize(srcbox));
3203 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3204 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3205
3206 Array4<T> const& d = this->array();
3207 Array4<T const> const& s = src.const_array();
3208 const auto dlo = amrex::lbound(destbox);
3209 const auto slo = amrex::lbound(srcbox);
3210 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3211 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3212 {
3213 if (s(i+offset.x,j+offset.y,k+offset.z,n+srccomp) != 0) {
3214 d(i,j,k,n+destcomp) /= s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
3215 }
3216 });
3217
3218 return *this;
3219}
3220
3231template <class T>
3232template <RunOn run_on>
3234BaseFab<T>::linInterp (const BaseFab<T>& f1, const Box& b1, int comp1,
3235 const BaseFab<T>& f2, const Box& b2, int comp2,
3236 Real t1, Real t2, Real t,
3237 const Box& b, int comp, int numcomp) noexcept
3238{
3239 if (amrex::almostEqual(t1,t2) || amrex::almostEqual(t1,t)) {
3240 return copy<run_on>(f1,b1,comp1,b,comp,numcomp);
3241 } else if (amrex::almostEqual(t2,t)) {
3242 return copy<run_on>(f2,b2,comp2,b,comp,numcomp);
3243 } else {
3244 Real alpha = (t2-t)/(t2-t1);
3245 Real beta = (t-t1)/(t2-t1);
3246 return linComb<run_on>(f1,b1,comp1,f2,b2,comp2,alpha,beta,b,comp,numcomp);
3247 }
3248}
3249
3250template <class T>
3251template <RunOn run_on>
3253BaseFab<T>::linInterp (const BaseFab<T>& f1, int comp1,
3254 const BaseFab<T>& f2, int comp2,
3255 Real t1, Real t2, Real t,
3256 const Box& b, int comp, int numcomp) noexcept
3257{
3258 if (amrex::almostEqual(t1,t2) || amrex::almostEqual(t1,t)) {
3259 return copy<run_on>(f1,b,comp1,b,comp,numcomp);
3260 } else if (amrex::almostEqual(t2,t)) {
3261 return copy<run_on>(f2,b,comp2,b,comp,numcomp);
3262 } else {
3263 Real alpha = (t2-t)/(t2-t1);
3264 Real beta = (t-t1)/(t2-t1);
3265 return linComb<run_on>(f1,b,comp1,f2,b,comp2,alpha,beta,b,comp,numcomp);
3266 }
3267}
3268
3269//
3270// New interfaces
3271//
3272
3273template <class T>
3274template <RunOn run_on>
3275void
3276BaseFab<T>::setVal (T const& val) noexcept
3277{
3278 this->setVal<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3279}
3280
3281template <class T>
3282template <RunOn run_on>
3283void
3284BaseFab<T>::setVal (T const& x, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept
3285{
3286 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3287 Array4<T> const& a = this->array();
3288 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG (run_on, bx, ncomp.n, i, j, k, n,
3289 {
3290 a(i,j,k,n+dcomp.i) = x;
3291 });
3292}
3293
3294template <class T>
3295template <RunOn run_on>
3296void
3297BaseFab<T>::setValIf (T const& val, const BaseFab<int>& mask) noexcept
3298{
3299 this->setValIf<run_on>(val, this->domain, mask, DestComp{0}, NumComps{this->nvar});
3300}
3301
3302template <class T>
3303template <RunOn run_on>
3304void
3305BaseFab<T>::setValIf (T const& val, Box const& bx, const BaseFab<int>& mask, DestComp dcomp, NumComps ncomp) noexcept
3306{
3307 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3308 Array4<T> const& a = this->array();
3309 Array4<int const> const& m = mask.const_array();
3310 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG (run_on, bx, ncomp.n, i, j, k, n,
3311 {
3312 if (m(i,j,k)) { a(i,j,k,n+dcomp.i) = val; }
3313 });
3314}
3315
3316template <class T>
3317template <RunOn run_on>
3318void
3319BaseFab<T>::setValIfNot (T const& val, const BaseFab<int>& mask) noexcept
3320{
3321 this->setValIfNot<run_on>(val, this->domain, mask, DestComp{0}, NumComps{this->nvar});
3322}
3323
3324template <class T>
3325template <RunOn run_on>
3326void
3327BaseFab<T>::setValIfNot (T const& val, Box const& bx, const BaseFab<int>& mask, DestComp dcomp, NumComps ncomp) noexcept
3328{
3329 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3330 Array4<T> const& a = this->array();
3331 Array4<int const> const& m = mask.const_array();
3332 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG (run_on, bx, ncomp.n, i, j, k, n,
3333 {
3334 if (!m(i,j,k)) { a(i,j,k,n+dcomp.i) = val; }
3335 });
3336}
3337
3338template <class T>
3339template <RunOn run_on>
3340void
3341BaseFab<T>::setComplement (T const& x, const Box& bx, DestComp dcomp, NumComps ncomp) noexcept
3342{
3343#ifdef AMREX_USE_GPU
3344 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
3345 Array4<T> const& a = this->array();
3346 amrex::ParallelFor(this->domain, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
3347 {
3348 if (! bx.contains(IntVect(AMREX_D_DECL(i,j,k)))) {
3349 for (int n = dcomp.i; n < ncomp.n+dcomp.i; ++n) {
3350 a(i,j,k,n) = x;
3351 }
3352 }
3353 });
3354 } else
3355#endif
3356 {
3357 const BoxList b_lst = amrex::boxDiff(this->domain,bx);
3358 for (auto const& b : b_lst) {
3359 this->setVal<RunOn::Host>(x, b, dcomp, ncomp);
3360 }
3361 }
3362}
3363
3364template <class T>
3365template <RunOn run_on>
3367BaseFab<T>::copy (const BaseFab<T>& src) noexcept
3368{
3369 this->copy<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3370 return *this;
3371}
3372
3373template <class T>
3374template <RunOn run_on>
3377 SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
3378{
3379 AMREX_ASSERT(this->domain.sameType(src.domain));
3380 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3381 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3382
3383 bx &= src.domain;
3384
3385 Array4<T> const& d = this->array();
3386 Array4<T const> const& s = src.const_array();
3387 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3388 {
3389 d(i,j,k,n+dcomp.i) = s(i,j,k,n+scomp.i);
3390 });
3391
3392 return *this;
3393}
3394
3395template <class T>
3396template <RunOn run_on>
3398BaseFab<T>::plus (T const& val) noexcept
3399{
3400 return this->plus<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3401}
3402
3403template <class T>
3404template <RunOn run_on>
3406BaseFab<T>::operator+= (T const& val) noexcept
3407{
3408 return this->plus<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3409}
3410
3411template <class T>
3412template <RunOn run_on>
3414BaseFab<T>::plus (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept
3415{
3416 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3417
3418 Array4<T> const& a = this->array();
3419 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3420 {
3421 a(i,j,k,n+dcomp.i) += val;
3422 });
3423
3424 return *this;
3425}
3426
3427template <class T>
3428template <RunOn run_on>
3430BaseFab<T>::plus (const BaseFab<T>& src) noexcept
3431{
3432 return this->plus<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3433}
3434
3435template <class T>
3436template <RunOn run_on>
3439{
3440 return this->plus<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3441}
3442
3443template <class T>
3444template <RunOn run_on>
3447 SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
3448{
3449 AMREX_ASSERT(this->domain.sameType(src.domain));
3450 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3451 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3452
3453 bx &= src.domain;
3454
3455 Array4<T> const& d = this->array();
3456 Array4<T const> const& s = src.const_array();
3457 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3458 {
3459 d(i,j,k,n+dcomp.i) += s(i,j,k,n+scomp.i);
3460 });
3461
3462 return *this;
3463}
3464
3465template <class T>
3466template <RunOn run_on>
3468BaseFab<T>::minus (T const& val) noexcept
3469{
3470 return this->minus<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3471}
3472
3473template <class T>
3474template <RunOn run_on>
3476BaseFab<T>::operator-= (T const& val) noexcept
3477{
3478 return this->minus<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3479}
3480
3481template <class T>
3482template <RunOn run_on>
3484BaseFab<T>::minus (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept
3485{
3486 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3487
3488 Array4<T> const& a = this->array();
3489 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3490 {
3491 a(i,j,k,n+dcomp.i) -= val;
3492 });
3493
3494 return *this;
3495}
3496
3497template <class T>
3498template <RunOn run_on>
3500BaseFab<T>::minus (const BaseFab<T>& src) noexcept
3501{
3502 return this->minus<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3503}
3504
3505template <class T>
3506template <RunOn run_on>
3509{
3510 return this->minus<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3511}
3512
3513template <class T>
3514template <RunOn run_on>
3517 SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
3518{
3519 AMREX_ASSERT(this->domain.sameType(src.domain));
3520 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3521 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3522
3523 bx &= src.domain;
3524
3525 Array4<T> const& d = this->array();
3526 Array4<T const> const& s = src.const_array();
3527 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3528 {
3529 d(i,j,k,n+dcomp.i) -= s(i,j,k,n+scomp.i);
3530 });
3531
3532 return *this;
3533}
3534
3535template <class T>
3536template <RunOn run_on>
3538BaseFab<T>::mult (T const& val) noexcept
3539{
3540 return this->mult<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3541}
3542
3543template <class T>
3544template <RunOn run_on>
3546BaseFab<T>::operator*= (T const& val) noexcept
3547{
3548 return this->mult<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3549}
3550
3551template <class T>
3552template <RunOn run_on>
3554BaseFab<T>::mult (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept
3555{
3556 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3557
3558 Array4<T> const& a = this->array();
3559 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3560 {
3561 a(i,j,k,n+dcomp.i) *= val;
3562 });
3563
3564 return *this;
3565}
3566
3567template <class T>
3568template <RunOn run_on>
3570BaseFab<T>::mult (const BaseFab<T>& src) noexcept
3571{
3572 return this->mult<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3573}
3574
3575template <class T>
3576template <RunOn run_on>
3579{
3580 return this->mult<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3581}
3582
3583template <class T>
3584template <RunOn run_on>
3587 SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
3588{
3589 AMREX_ASSERT(this->domain.sameType(src.domain));
3590 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3591 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3592
3593 bx &= src.domain;
3594
3595 Array4<T> const& d = this->array();
3596 Array4<T const> const& s = src.const_array();
3597 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3598 {
3599 d(i,j,k,n+dcomp.i) *= s(i,j,k,n+scomp.i);
3600 });
3601
3602 return *this;
3603}
3604
3605template <class T>
3606template <RunOn run_on>
3608BaseFab<T>::divide (T const& val) noexcept
3609{
3610 return this->divide<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3611}
3612
3613template <class T>
3614template <RunOn run_on>
3616BaseFab<T>::operator/= (T const& val) noexcept
3617{
3618 return this->divide<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3619}
3620
3621template <class T>
3622template <RunOn run_on>
3624BaseFab<T>::divide (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept
3625{
3626 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3627
3628 Array4<T> const& a = this->array();
3629 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3630 {
3631 a(i,j,k,n+dcomp.i) /= val;
3632 });
3633
3634 return *this;
3635}
3636
3637template <class T>
3638template <RunOn run_on>
3640BaseFab<T>::divide (const BaseFab<T>& src) noexcept
3641{
3642 return this->divide<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3643}
3644
3645template <class T>
3646template <RunOn run_on>
3649{
3650 return this->divide<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3651}
3652
3653template <class T>
3654template <RunOn run_on>
3657 SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
3658{
3659 AMREX_ASSERT(this->domain.sameType(src.domain));
3660 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3661 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3662
3663 bx &= src.domain;
3664
3665 Array4<T> const& d = this->array();
3666 Array4<T const> const& s = src.const_array();
3667 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3668 {
3669 d(i,j,k,n+dcomp.i) /= s(i,j,k,n+scomp.i);
3670 });
3671
3672 return *this;
3673}
3674
3675template <class T>
3676template <RunOn run_on>
3679{
3680 return this->negate<run_on>(this->domain, DestComp{0}, NumComps{this->nvar});
3681}
3682
3683template <class T>
3684template <RunOn run_on>
3686BaseFab<T>::negate (const Box& bx, DestComp dcomp, NumComps ncomp) noexcept
3687{
3688 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3689
3690 Array4<T> const& a = this->array();
3691 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3692 {
3693 a(i,j,k,n+dcomp.i) = -a(i,j,k,n+dcomp.i);
3694 });
3695
3696 return *this;
3697}
3698
3699template <class T>
3700template <RunOn run_on>
3702BaseFab<T>::invert (T const& r) noexcept
3703{
3704 return this->invert<run_on>(r, this->domain, DestComp{0}, NumComps{this->nvar});
3705}
3706
3707template <class T>
3708template <RunOn run_on>
3710BaseFab<T>::invert (T const& r, const Box& bx, DestComp dcomp, NumComps ncomp) noexcept
3711{
3712 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3713
3714 Array4<T> const& a = this->array();
3715 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3716 {
3717 a(i,j,k,n+dcomp.i) = r / a(i,j,k,n+dcomp.i);
3718 });
3719
3720 return *this;
3721}
3722
3723template <class T>
3724template <RunOn run_on>
3725T
3726BaseFab<T>::sum (const Box& bx, DestComp dcomp, NumComps ncomp) const noexcept
3727{
3728 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3729
3730 T r = 0;
3731 Array4<T const> const& a = this->const_array();
3732#ifdef AMREX_USE_GPU
3733 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
3734 ReduceOps<ReduceOpSum> reduce_op;
3735 ReduceData<T> reduce_data(reduce_op);
3736 using ReduceTuple = typename decltype(reduce_data)::Type;
3737 reduce_op.eval(bx, reduce_data,
3738 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
3739 {
3740 T t = 0;
3741 for (int n = 0; n < ncomp.n; ++n) {
3742 t += a(i,j,k,n+dcomp.i);
3743 }
3744 return { t };
3745 });
3746 ReduceTuple hv = reduce_data.value(reduce_op);
3747 r = amrex::get<0>(hv);
3748 } else
3749#endif
3750 {
3751 amrex::LoopOnCpu(bx, ncomp.n, [=,&r] (int i, int j, int k, int n) noexcept
3752 {
3753 r += a(i,j,k,n+dcomp.i);
3754 });
3755 }
3756
3757 return r;
3758}
3759
3760template <class T>
3761template <RunOn run_on>
3762T
3763BaseFab<T>::dot (const BaseFab<T>& src, const Box& bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) const noexcept
3764{
3765 AMREX_ASSERT(this->domain.sameType(src.domain));
3766 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3767 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3768
3769 T r = 0;
3770 Array4<T const> const& d = this->const_array();
3771 Array4<T const> const& s = src.const_array();
3772#ifdef AMREX_USE_GPU
3773 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
3774 ReduceOps<ReduceOpSum> reduce_op;
3775 ReduceData<T> reduce_data(reduce_op);
3776 using ReduceTuple = typename decltype(reduce_data)::Type;
3777 reduce_op.eval(bx, reduce_data,
3778 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
3779 {
3780 T t = 0;
3781 for (int n = 0; n < ncomp.n; ++n) {
3782 t += d(i,j,k,n+dcomp.i) * s(i,j,k,n+scomp.i);
3783 }
3784 return { t };
3785 });
3786 ReduceTuple hv = reduce_data.value(reduce_op);
3787 r = amrex::get<0>(hv);
3788 } else
3789#endif
3790 {
3791 amrex::LoopOnCpu(bx, ncomp.n, [=,&r] (int i, int j, int k, int n) noexcept
3792 {
3793 r += d(i,j,k,n+dcomp.i) * s(i,j,k,n+scomp.i);
3794 });
3795 }
3796
3797 return r;
3798}
3799
3800template <class T>
3801template <RunOn run_on>
3802T
3803BaseFab<T>::dot (const Box& bx, int destcomp, int numcomp) const noexcept
3804{
3805 return dot<run_on>(bx, DestComp{destcomp}, NumComps{numcomp});
3806}
3807
3808
3809template <class T>
3810template <RunOn run_on>
3811T
3812BaseFab<T>::dot (const Box& bx, DestComp dcomp, NumComps ncomp) const noexcept
3813{
3814 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3815
3816 T r = 0;
3817 Array4<T const> const& a = this->const_array();
3818#ifdef AMREX_USE_GPU
3819 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
3820 ReduceOps<ReduceOpSum> reduce_op;
3821 ReduceData<T> reduce_data(reduce_op);
3822 using ReduceTuple = typename decltype(reduce_data)::Type;
3823 reduce_op.eval(bx, reduce_data,
3824 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
3825 {
3826 T t = 0;
3827 for (int n = 0; n < ncomp.n; ++n) {
3828 t += a(i,j,k,n+dcomp.i)*a(i,j,k,n+dcomp.i);
3829 }
3830 return { t };
3831 });
3832 ReduceTuple hv = reduce_data.value(reduce_op);
3833 r = amrex::get<0>(hv);
3834 } else
3835#endif
3836 {
3837 amrex::LoopOnCpu(bx, ncomp.n, [=,&r] (int i, int j, int k, int n) noexcept
3838 {
3839 r += a(i,j,k,n+dcomp.i)*a(i,j,k,n+dcomp.i);
3840 });
3841 }
3842
3843 return r;
3844}
3845
3846template <class T>
3847template <RunOn run_on>
3848T
3849BaseFab<T>::dotmask (const BaseFab<T>& src, const Box& bx, const BaseFab<int>& mask,
3850 SrcComp scomp, DestComp dcomp, NumComps ncomp) const noexcept
3851{
3852 AMREX_ASSERT(this->domain.sameType(src.domain));
3853 AMREX_ASSERT(this->domain.sameType(mask.domain));
3854 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3855 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3856
3857 T r = 0;
3858 Array4<T const> const& d = this->const_array();
3859 Array4<T const> const& s = src.const_array();
3860 Array4<int const> const& m = mask.const_array();
3861#ifdef AMREX_USE_GPU
3862 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
3863 ReduceOps<ReduceOpSum> reduce_op;
3864 ReduceData<T> reduce_data(reduce_op);
3865 using ReduceTuple = typename decltype(reduce_data)::Type;
3866 reduce_op.eval(bx, reduce_data,
3867 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
3868 {
3869 T t = 0;
3870 T mi = static_cast<T>(static_cast<int>(static_cast<bool>(m(i,j,k))));
3871 for (int n = 0; n < ncomp.n; ++n) {
3872 t += d(i,j,k,n+dcomp.i)*s(i,j,k,n+scomp.i)*mi;
3873 }
3874 return { t };
3875 });
3876 ReduceTuple hv = reduce_data.value(reduce_op);
3877 r = amrex::get<0>(hv);
3878 } else
3879#endif
3880 {
3881 amrex::LoopOnCpu(bx, ncomp.n, [=,&r] (int i, int j, int k, int n) noexcept
3882 {
3883 int mi = static_cast<int>(static_cast<bool>(m(i,j,k)));
3884 r += d(i,j,k,n+dcomp.i)*s(i,j,k,n+scomp.i)*mi;
3885 });
3886 }
3887
3888 return r;
3889}
3890
3891}
3892
3893#endif /*BL_BASEFAB_H*/
#define BL_ASSERT(EX)
Definition AMReX_BLassert.H:39
#define AMREX_ASSERT(EX)
Definition AMReX_BLassert.H:38
#define AMREX_FORCE_INLINE
Definition AMReX_Extension.H:119
#define AMREX_DEFAULT_RUNON
Definition AMReX_GpuControl.H:69
#define AMREX_CUDA_SAFE_CALL(call)
Definition AMReX_GpuError.H:73
#define AMREX_HOST_DEVICE_FOR_1D(...)
Definition AMReX_GpuLaunchMacrosC.nolint.H:105
#define AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(where_to_run, box, nc, i, j, k, n, block)
Definition AMReX_GpuLaunch.nolint.H:75
#define AMREX_GPU_DEVICE
Definition AMReX_GpuQualifiers.H:18
#define AMREX_GPU_HOST_DEVICE
Definition AMReX_GpuQualifiers.H:20
int idir
Definition AMReX_HypreMLABecLap.cpp:1143
Array4< int const > offset
Definition AMReX_HypreMLABecLap.cpp:1139
Real * pdst
Definition AMReX_HypreMLABecLap.cpp:1140
Array4< int const > mask
Definition AMReX_InterpFaceRegister.cpp:93
#define AMREX_LOOP_3D(bx, i, j, k, block)
Definition AMReX_Loop.nolint.H:4
#define AMREX_LOOP_4D(bx, ncomp, i, j, k, n, block)
Definition AMReX_Loop.nolint.H:16
void amrex_mempool_free(void *p)
Definition AMReX_MemPool.cpp:80
void * amrex_mempool_alloc(size_t nbytes)
Definition AMReX_MemPool.cpp:74
#define AMREX_D_TERM(a, b, c)
Definition AMReX_SPACE.H:172
#define AMREX_D_DECL(a, b, c)
Definition AMReX_SPACE.H:171
A virtual base class for objects that manage their own dynamic memory allocation.
Definition AMReX_Arena.H:132
A FortranArrayBox(FAB)-like object.
Definition AMReX_BaseFab.H:190
int maskGE(BaseFab< int > &mask, T const &val, int comp=0) const noexcept
Same as above except mark cells with value greater than or equal to val.
Definition AMReX_BaseFab.H:2435
BaseFab< T > & saxpy(T a, const BaseFab< T > &x, const Box &srcbox, const Box &destbox, int srccomp, int destcomp, int numcomp=1) noexcept
FAB SAXPY (y[i] <- y[i] + a * x[i]), in place.
Definition AMReX_BaseFab.H:2491
Array4< T const > const_array() const noexcept
Definition AMReX_BaseFab.H:418
BaseFab< T > & copy(const BaseFab< T > &src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx does not intersect with src fab.
Definition AMReX_BaseFab.H:3376
T sum(int comp, int numcomp=1) const noexcept
Returns sum of given component of FAB state vector.
Definition AMReX_BaseFab.H:2721
gpuStream_t alloc_stream
Definition AMReX_BaseFab.H:1171
Real norminfmask(const Box &subbox, const BaseFab< int > &mask, int scomp=0, int ncomp=1) const noexcept
Definition AMReX_BaseFab.H:1864
BaseFab< T > & divide(T const &val) noexcept
Scalar division on the whole domain and all components.
Definition AMReX_BaseFab.H:3608
const int * hiVect() const noexcept
Returns the upper corner of the domain.
Definition AMReX_BaseFab.H:329
BaseFab< T > & minus(const BaseFab< T > &src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx does not intersect with src fab.
Definition AMReX_BaseFab.H:3516
BaseFab< T > & lockAdd(const BaseFab< T > &src, const Box &srcbox, const Box &destbox, int srccomp, int destcomp, int numcomp) noexcept
Atomically add srcbox region of src FAB to destbox region of this FAB. The srcbox and destbox must be...
Definition AMReX_BaseFab.H:2902
static void Initialize()
std::size_t copyToMem(const Box &srcbox, int srccomp, int numcomp, void *dst) const noexcept
Copy from the srcbox of this Fab to raw memory and return the number of bytes copied.
Definition AMReX_BaseFab.H:1745
std::size_t addFromMem(const Box &dstbox, int dstcomp, int numcomp, const void *src) noexcept
Add from raw memory to the dstbox of this Fab and return the number of bytes copied.
Definition AMReX_BaseFab.H:1800
std::size_t nBytesOwned() const noexcept
Definition AMReX_BaseFab.H:271
BaseFab< T > & copy(const BaseFab< T > &src) noexcept
Definition AMReX_BaseFab.H:3367
BaseFab< T > & addproduct(const Box &destbox, int destcomp, int numcomp, const BaseFab< T > &src1, int comp1, const BaseFab< T > &src2, int comp2) noexcept
y[i] <- y[i] + x1[i] * x2[i])
Definition AMReX_BaseFab.H:2556
BaseFab< T > & minus(T const &val) noexcept
Scalar subtraction on the whole domain and all components.
Definition AMReX_BaseFab.H:3468
int maskLT(BaseFab< int > &mask, T const &val, int comp=0) const noexcept
Compute mask array with value of 1 in cells where BaseFab has value less than val,...
Definition AMReX_BaseFab.H:2251
BaseFab< T > & plus(T const &val) noexcept
Scalar addition on the whole domain and all components.
Definition AMReX_BaseFab.H:3398
BaseFab< T > & mult(T const &val, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3554
BaseFab< T > & mult(const BaseFab< T > &src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx does not intersect with src fab.
Definition AMReX_BaseFab.H:3586
std::size_t nBytes(const Box &bx, int ncomps) const noexcept
Returns bytes used in the Box for those components.
Definition AMReX_BaseFab.H:276
void setPtr(T *p, Long sz) noexcept
Definition AMReX_BaseFab.H:376
BaseFab< T > & linComb(const BaseFab< T > &f1, const Box &b1, int comp1, const BaseFab< T > &f2, const Box &b2, int comp2, Real alpha, Real beta, const Box &b, int comp, int numcomp=1) noexcept
Linear combination. Result is alpha*f1 + beta*f2. Data is taken from b1 region of f1,...
Definition AMReX_BaseFab.H:2580
void define()
Allocates memory for the BaseFab<T>.
Definition AMReX_BaseFab.H:1457
BaseFab< T > & operator*=(T const &val) noexcept
Definition AMReX_BaseFab.H:3546
void resize(const Box &b, int N=1, Arena *ar=nullptr)
This function resizes a BaseFab so it covers the Box b with N components.
Definition AMReX_BaseFab.H:1628
void clear()
The function returns the BaseFab to the invalid state. The memory is freed.
Definition AMReX_BaseFab.H:1691
const IntVect & smallEnd() const noexcept
Returns the lower corner of the domain See class Box for analogue.
Definition AMReX_BaseFab.H:306
BaseFab< T > & mult(T const &r, int comp, int numcomp=1) noexcept
Scalar multiplication, except control which components are multiplied.
Definition AMReX_BaseFab.H:3048
BaseFab< T > & atomicAdd(const BaseFab< T > &x) noexcept
Atomic FAB addition (a[i] <- a[i] + b[i]).
Definition AMReX_BaseFab.H:2481
int maskEQ(BaseFab< int > &mask, T const &val, int comp=0) const noexcept
Same as above except mark cells with value equal to val.
Definition AMReX_BaseFab.H:2343
const int * loVect() const noexcept
Returns the lower corner of the domain.
Definition AMReX_BaseFab.H:319
bool contains(const BaseFab< T > &fab) const noexcept
Returns true if the domain of fab is totally contained within the domain of this BaseFab.
Definition AMReX_BaseFab.H:335
bool isAllocated() const noexcept
Returns true if the data for the FAB has been allocated.
Definition AMReX_BaseFab.H:436
std::unique_ptr< T, DataDeleter > release() noexcept
Release ownership of memory.
Definition AMReX_BaseFab.H:1727
void setVal(T const &x, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3284
BaseFab< T > & operator-=(T const &val) noexcept
Definition AMReX_BaseFab.H:3476
const Box & box() const noexcept
Returns the domain (box) where the array is defined.
Definition AMReX_BaseFab.H:294
void setValIf(T const &val, Box const &bx, const BaseFab< int > &mask, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3305
static void Finalize()
Array4< T > array() noexcept
Definition AMReX_BaseFab.H:400
IntVect indexFromValue(const Box &subbox, int comp, T const &value) const noexcept
Definition AMReX_BaseFab.H:2160
BaseFab< T > & mult(const BaseFab< T > &src) noexcept
Definition AMReX_BaseFab.H:3570
bool shared_memory
Is the memory allocated in shared memory?
Definition AMReX_BaseFab.H:1169
int maskLE(BaseFab< int > &mask, T const &val, int comp=0) const noexcept
Same as above except mark cells with value less than or equal to val.
Definition AMReX_BaseFab.H:2297
void setValIf(T const &val, const BaseFab< int > &mask) noexcept
Definition AMReX_BaseFab.H:3297
BaseFab< T > & plus(const BaseFab< T > &src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx does not intersect with src fab.
Definition AMReX_BaseFab.H:3446
void setValIfNot(T const &val, const BaseFab< int > &mask) noexcept
Definition AMReX_BaseFab.H:3319
BaseFab< T > & xpay(T a, const BaseFab< T > &x, const Box &srcbox, const Box &destbox, int srccomp, int destcomp, int numcomp=1) noexcept
FAB XPAY (y[i] <- x[i] + a * y[i])
Definition AMReX_BaseFab.H:2528
T & operator()(const IntVect &p, int N) noexcept
Returns a reference to the Nth component value defined at position p in the domain....
Definition AMReX_BaseFab.H:1268
std::size_t nBytes() const noexcept
Returns how many bytes used.
Definition AMReX_BaseFab.H:269
std::size_t copyFromMem(const Box &dstbox, int dstcomp, int numcomp, const void *src) noexcept
Copy from raw memory to the dstbox of this Fab and return the number of bytes copied.
Definition AMReX_BaseFab.H:1772
BaseFab< T > & negate() noexcept
on the whole domain and all components
Definition AMReX_BaseFab.H:3678
BaseFab< T > & minus(const BaseFab< T > &src) noexcept
Definition AMReX_BaseFab.H:3500
BaseFab< T > & minus(const BaseFab< T > &src, int srccomp, int destcomp, int numcomp=1) noexcept
Subtract src components (srccomp:srccomp+numcomp-1) to this FABs components (destcomp:destcomp+numcom...
Definition AMReX_BaseFab.H:3006
T value_type
Definition AMReX_BaseFab.H:195
void SetBoxType(const IndexType &typ) noexcept
Change the Box type without change the length.
Definition AMReX_BaseFab.H:981
Array4< T const > array() const noexcept
Definition AMReX_BaseFab.H:382
T maxabs(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2121
BaseFab< T > & operator+=(T const &val) noexcept
Definition AMReX_BaseFab.H:3406
void setComplement(T const &x, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
setVal on the complement of bx in the fab's domain
Definition AMReX_BaseFab.H:3341
BaseFab< T > & minus(T const &val, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3484
Long truesize
nvar*numpts that was allocated on heap.
Definition AMReX_BaseFab.H:1167
void setVal(T const &val) noexcept
Set value on the whole domain and all components.
Definition AMReX_BaseFab.H:3276
const int * nCompPtr() const noexcept
for calls to fortran.
Definition AMReX_BaseFab.H:283
Array4< T const > const_array(int start_comp, int num_comps) const noexcept
Definition AMReX_BaseFab.H:430
Box domain
My index space.
Definition AMReX_BaseFab.H:1165
bool contains(const Box &bx) const noexcept
Returns true if bx is totally contained within the domain of this BaseFab.
Definition AMReX_BaseFab.H:344
T * dptr
The data pointer.
Definition AMReX_BaseFab.H:1164
BaseFab< T > & shift(const IntVect &v) noexcept
Perform shifts upon the domain of the BaseFab. They are completely analogous to the corresponding Box...
Definition AMReX_BaseFab.H:1341
BaseFab< T > & divide(T const &val, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3624
T max(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2041
BaseFab< T > & copy(const BaseFab< T > &src, const Box &srcbox, int srccomp, const Box &destbox, int destcomp, int numcomp) noexcept
The copy functions copy the contents of one BaseFab into another. The destination BaseFab is always t...
Definition AMReX_BaseFab.H:1415
Array4< T > array(int start_comp, int num_comps) noexcept
Definition AMReX_BaseFab.H:412
int nvar
Number components.
Definition AMReX_BaseFab.H:1166
T dot(const Box &xbx, int xcomp, const BaseFab< T > &y, const Box &ybx, int ycomp, int numcomp=1) const noexcept
Dot product of x (i.e.,this) and y.
Definition AMReX_BaseFab.H:2617
BaseFab< T > & divide(T const &r, int comp, int numcomp=1) noexcept
As above except specify which components.
Definition AMReX_BaseFab.H:3106
Array4< T const > array(int start_comp) const noexcept
Definition AMReX_BaseFab.H:388
BaseFab< T > & operator/=(T const &val) noexcept
Definition AMReX_BaseFab.H:3616
std::pair< T, T > minmax(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2079
Array4< T const > const_array(int start_comp) const noexcept
Definition AMReX_BaseFab.H:424
void fill_snan() noexcept
Definition AMReX_BaseFab.H:1375
void setVal(T const &x, const Box &bx, int dcomp, int ncomp) noexcept
The setVal functions set sub-regions in the BaseFab to a constant value. This most general form speci...
Definition AMReX_BaseFab.H:1399
Long size() const noexcept
Returns the total number of points of all components.
Definition AMReX_BaseFab.H:291
BaseFab< T > & plus(T const &r, const Box &b, int comp=0, int numcomp=1) noexcept
Scalar addition (a[i] <- a[i] + r), most general.
Definition AMReX_BaseFab.H:2777
BaseFab< T > & operator=(const BaseFab< T > &rhs)=delete
void getVal(T *data, const IntVect &pos, int N, int numcomp) const noexcept
This function puts numcomp component values, starting at component N, from position pos in the domain...
Definition AMReX_BaseFab.H:1315
const IntVect & bigEnd() const noexcept
Returns the upper corner of the domain. See class Box for analogue.
Definition AMReX_BaseFab.H:309
Array4< T > array(int start_comp) noexcept
Definition AMReX_BaseFab.H:406
Elixir elixir() noexcept
Definition AMReX_BaseFab.H:1670
Long numPts() const noexcept
Returns the number of points.
Definition AMReX_BaseFab.H:288
const T * dataPtr(int n=0) const noexcept
Same as above except works on const FABs.
Definition AMReX_BaseFab.H:364
void setValIfNot(T const &val, Box const &bx, const BaseFab< int > &mask, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3327
BaseFab< T > & mult(T const &val) noexcept
Scalar multiplication on the whole domain and all components.
Definition AMReX_BaseFab.H:3538
BaseFab< T > & divide(const BaseFab< T > &src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx does not intersect with src fab.
Definition AMReX_BaseFab.H:3656
void setValIfNot(T const &val, const Box &bx, const BaseFab< int > &mask, int nstart, int num) noexcept
Definition AMReX_BaseFab.H:1407
BaseFab< T > & shiftHalf(int dir, int n_cell) noexcept
Perform shifts upon the domain of the BaseFab. They are completely analogous to the corresponding Box...
Definition AMReX_BaseFab.H:1365
void prefetchToDevice() const noexcept
Definition AMReX_BaseFab.H:1235
T * dataPtr(int n=0) noexcept
Returns a pointer to an object of type T that is the value of the Nth component associated with the c...
Definition AMReX_BaseFab.H:355
bool ptr_owner
Owner of T*?
Definition AMReX_BaseFab.H:1168
virtual ~BaseFab() noexcept
The destructor deletes the array memory.
Definition AMReX_BaseFab.H:1575
IntVect length() const noexcept
Returns a pointer to an array of SPACEDIM integers giving the length of the domain in each direction.
Definition AMReX_BaseFab.H:300
IntVect maxIndex(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2225
BaseFab< T > & protected_divide(const BaseFab< T > &src) noexcept
Divide wherever "src" is "true" or "non-zero".
Definition AMReX_BaseFab.H:3164
friend class BaseFab
Definition AMReX_BaseFab.H:193
BaseFab< T > & invert(T const &r, const Box &b, int comp=0, int numcomp=1) noexcept
Most general version, specify subbox and which components.
Definition AMReX_BaseFab.H:2761
Array4< T const > array(int start_comp, int num_comps) const noexcept
Definition AMReX_BaseFab.H:394
int nComp() const noexcept
Returns the number of components.
Definition AMReX_BaseFab.H:280
int maskGT(BaseFab< int > &mask, T const &val, int comp=0) const noexcept
Same as above except mark cells with value greater than val.
Definition AMReX_BaseFab.H:2389
T dotmask(const BaseFab< int > &mask, const Box &xbx, int xcomp, const BaseFab< T > &y, const Box &ybx, int ycomp, int numcomp) const noexcept
Definition AMReX_BaseFab.H:2667
BaseFab< T > & plus(const BaseFab< T > &src) noexcept
Definition AMReX_BaseFab.H:3430
BaseFab() noexcept=default
Construct an empty BaseFab, which must be resized (see BaseFab::resize) before use.
BaseFab< T > & divide(const BaseFab< T > &src) noexcept
Definition AMReX_BaseFab.H:3640
void prefetchToHost() const noexcept
Definition AMReX_BaseFab.H:1203
T min(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2003
void setComplement(T const &x, const Box &b, int ns, int num) noexcept
This function is analogous to the fourth form of setVal above, except that instead of setting values ...
Definition AMReX_BaseFab.H:1828
BaseFab< T > & linInterp(const BaseFab< T > &f1, const Box &b1, int comp1, const BaseFab< T > &f2, const Box &b2, int comp2, Real t1, Real t2, Real t, const Box &b, int comp, int numcomp=1) noexcept
Linear interpolation / extrapolation. Result is (t2-t)/(t2-t1)*f1 + (t-t1)/(t2-t1)*f2 Data is taken f...
Definition AMReX_BaseFab.H:3234
IntVect minIndex(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2199
T * dataPtr(const IntVect &p, int n=0) noexcept
Definition AMReX_BaseFab.H:1178
void abs() noexcept
Compute absolute value for all components of this FAB.
Definition AMReX_BaseFab.H:1836
void getVal(T *data, const IntVect &pos) const noexcept
Same as above, except that starts at component 0 and copies all comps.
Definition AMReX_BaseFab.H:1333
BaseFab< T > & plus(T const &val, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3414
Real norm(int p, int scomp=0, int numcomp=1) const
Compute the Lp-norm of this FAB using components (scomp : scomp+ncomp-1). p < 0 -> ERROR p = 0 -> inf...
Definition AMReX_BaseFab.H:1908
A class for managing a List of Boxes that share a common IndexType. This class implements operations ...
Definition AMReX_BoxList.H:52
__host__ __device__ const IntVectND< dim > & bigEnd() const &noexcept
Return the inclusive upper bound of the box.
Definition AMReX_Box.H:123
__host__ __device__ const int * hiVect() const &noexcept
Return a constant pointer the array of high end coordinates. Useful for calls to FORTRAN.
Definition AMReX_Box.H:191
__host__ __device__ Long numPts() const noexcept
Return the number of points contained in the BoxND.
Definition AMReX_Box.H:356
__host__ __device__ IntVectND< dim > length() const noexcept
Return the length of the BoxND.
Definition AMReX_Box.H:154
__host__ __device__ bool contains(const IntVectND< dim > &p) const noexcept
Return true if argument is contained within BoxND.
Definition AMReX_Box.H:212
__host__ __device__ const int * loVect() const &noexcept
Return a constant pointer the array of low end coordinates. Useful for calls to FORTRAN.
Definition AMReX_Box.H:186
__host__ __device__ BoxND & setType(const IndexTypeND< dim > &t) noexcept
Set indexing type.
Definition AMReX_Box.H:505
__host__ __device__ bool ok() const noexcept
Checks if it is a proper BoxND (including a valid type).
Definition AMReX_Box.H:208
__host__ __device__ const IntVectND< dim > & smallEnd() const &noexcept
Return the inclusive lower bound of the box.
Definition AMReX_Box.H:111
GPU-compatible tuple.
Definition AMReX_Tuple.H:98
static int deviceId() noexcept
Definition AMReX_GpuDevice.cpp:692
static int devicePropMajor() noexcept
Definition AMReX_GpuDevice.H:203
Definition AMReX_GpuElixir.H:13
__host__ static __device__ constexpr IntVectND< dim > TheMinVector() noexcept
Definition AMReX_IntVect.H:819
Dynamically allocated vector for trivially copyable data.
Definition AMReX_PODVector.H:308
T * data() noexcept
Definition AMReX_PODVector.H:666
Definition AMReX_Reduce.H:453
Type value()
Definition AMReX_Reduce.H:488
Definition AMReX_Reduce.H:612
std::enable_if_t< IsFabArray< MF >::value > eval(MF const &mf, IntVect const &nghost, D &reduce_data, F &&f)
Definition AMReX_Reduce.H:748
amrex_real Real
Floating Point Type for Fields.
Definition AMReX_REAL.H:79
amrex_long Long
Definition AMReX_INT.H:30
__host__ __device__ Dim3 ubound(Array4< T > const &a) noexcept
Return the inclusive upper bounds of an Array4 in Dim3 form.
Definition AMReX_Array4.H:1331
__host__ __device__ Dim3 length(Array4< T > const &a) noexcept
Return the spatial extents of an Array4 in Dim3 form.
Definition AMReX_Array4.H:1345
__host__ __device__ Dim3 lbound(Array4< T > const &a) noexcept
Return the inclusive lower bounds of an Array4 in Dim3 form.
Definition AMReX_Array4.H:1317
std::array< T, N > Array
Definition AMReX_Array.H:26
__host__ __device__ AMREX_FORCE_INLINE T Exch(T *address, T val) noexcept
Definition AMReX_GpuAtomic.H:487
void streamSynchronize() noexcept
Definition AMReX_GpuDevice.H:310
void dtoh_memcpy_async(void *p_h, const void *p_d, const std::size_t sz) noexcept
Definition AMReX_GpuDevice.H:435
bool inLaunchRegion() noexcept
Definition AMReX_GpuControl.H:88
bool notInLaunchRegion() noexcept
Definition AMReX_GpuControl.H:89
void htod_memcpy_async(void *p_d, const void *p_h, const std::size_t sz) noexcept
Definition AMReX_GpuDevice.H:421
gpuStream_t gpuStream() noexcept
Definition AMReX_GpuDevice.H:291
__host__ __device__ AMREX_FORCE_INLINE void Add(T *const sum, T const value) noexcept
Definition AMReX_GpuAtomic.H:636
Definition AMReX_Amr.cpp:50
MakeType
Definition AMReX_MakeType.H:7
@ make_deep_copy
Definition AMReX_MakeType.H:7
@ make_alias
Definition AMReX_MakeType.H:7
__host__ __device__ void ignore_unused(const Ts &...)
This shuts up the compiler about unused variables.
Definition AMReX.H:139
int nComp(FabArrayBase const &fa)
Definition AMReX_FabArrayBase.cpp:2852
std::enable_if_t< std::is_integral_v< T > > ParallelFor(TypeList< CTOs... > ctos, std::array< int, sizeof...(CTOs)> const &runtime_options, T N, F &&f)
Definition AMReX_CTOParallelForImpl.H:193
__host__ __device__ Array4< T > makeArray4(T *p, Box const &bx, int ncomp) noexcept
Definition AMReX_BaseFab.H:92
RunOn
Definition AMReX_GpuControl.H:65
std::enable_if_t< std::is_arithmetic_v< T > > placementNew(T *const, Long)
Definition AMReX_BaseFab.H:99
cudaStream_t gpuStream_t
Definition AMReX_GpuControl.H:79
bool InitSNaN() noexcept
Definition AMReX.cpp:185
Long TotalBytesAllocatedInFabs() noexcept
Definition AMReX_BaseFab.cpp:66
BoxND< 3 > Box
Box is an alias for amrex::BoxND instantiated with AMREX_SPACEDIM.
Definition AMReX_BaseFwd.H:30
void BaseFab_Initialize()
Definition AMReX_BaseFab.cpp:30
__host__ __device__ constexpr const T & min(const T &a, const T &b) noexcept
Definition AMReX_Algorithm.H:24
void BaseFab_Finalize()
Definition AMReX_BaseFab.cpp:59
__host__ __device__ Dim3 begin(BoxND< dim > const &box) noexcept
Definition AMReX_Box.H:2006
void ResetTotalBytesAllocatedInFabsHWM() noexcept
Definition AMReX_BaseFab.cpp:134
BoxList boxDiff(const Box &b1in, const Box &b2)
Returns BoxList defining the compliment of b2 in b1in.
Definition AMReX_BoxList.cpp:599
__host__ __device__ std::enable_if_t< std::is_floating_point_v< T >, bool > almostEqual(T x, T y, int ulp=2)
Definition AMReX_Algorithm.H:122
IntVectND< 3 > IntVect
IntVect is an alias for amrex::IntVectND instantiated with AMREX_SPACEDIM.
Definition AMReX_BaseFwd.H:33
Long TotalBytesAllocatedInFabsHWM() noexcept
Definition AMReX_BaseFab.cpp:83
__host__ __device__ constexpr const T & max(const T &a, const T &b) noexcept
Definition AMReX_Algorithm.H:44
Long TotalCellsAllocatedInFabsHWM() noexcept
Definition AMReX_BaseFab.cpp:117
Long TotalCellsAllocatedInFabs() noexcept
Definition AMReX_BaseFab.cpp:100
void Error(const std::string &msg)
Print out message to cerr and exit via amrex::Abort().
Definition AMReX.cpp:235
std::enable_if_t< std::is_trivially_destructible_v< T > > placementDelete(T *const, Long)
Definition AMReX_BaseFab.H:124
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition AMReX.cpp:241
void LoopOnCpu(Dim3 lo, Dim3 hi, F const &f) noexcept
Definition AMReX_Loop.H:365
void update_fab_stats(Long n, Long s, size_t szt) noexcept
Definition AMReX_BaseFab.cpp:146
__host__ __device__ Dim3 end(BoxND< dim > const &box) noexcept
Definition AMReX_Box.H:2015
A multidimensional array accessor.
Definition AMReX_Array4.H:283
__host__ __device__ constexpr std::size_t size() const noexcept
Total number of elements in the ArrayND's index region.
Definition AMReX_Array4.H:638
__host__ __device__ T * ptr(idx... i) const noexcept
Multi-index ptr() for accessing pointer to element.
Definition AMReX_Array4.H:556
Definition AMReX_DataAllocator.H:9
void * alloc(std::size_t sz) const noexcept
Definition AMReX_DataAllocator.H:16
Arena * arena() const noexcept
Definition AMReX_DataAllocator.H:24
Definition AMReX_DataAllocator.H:29
Definition AMReX_BaseFab.H:77
int i
Definition AMReX_BaseFab.H:80
__host__ __device__ DestComp(int ai) noexcept
Definition AMReX_BaseFab.H:79
Definition AMReX_Dim3.H:12
int x
Definition AMReX_Dim3.H:12
Definition AMReX_BaseFab.H:83
__host__ __device__ NumComps(int an) noexcept
Definition AMReX_BaseFab.H:85
int n
Definition AMReX_BaseFab.H:86
Definition AMReX_BaseFab.H:71
__host__ __device__ SrcComp(int ai) noexcept
Definition AMReX_BaseFab.H:73
int i
Definition AMReX_BaseFab.H:74