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
44namespace amrex
45{
46
48
49extern std::atomic<Long> atomic_total_bytes_allocated_in_fabs;
50extern std::atomic<Long> atomic_total_bytes_allocated_in_fabs_hwm;
51extern std::atomic<Long> atomic_total_cells_allocated_in_fabs;
52extern std::atomic<Long> atomic_total_cells_allocated_in_fabs_hwm;
53extern Long private_total_bytes_allocated_in_fabs;
54extern Long private_total_bytes_allocated_in_fabs_hwm;
55extern Long private_total_cells_allocated_in_fabs;
56extern Long private_total_cells_allocated_in_fabs_hwm;
57#ifdef AMREX_USE_OMP
58#pragma omp threadprivate(private_total_bytes_allocated_in_fabs)
59#pragma omp threadprivate(private_total_bytes_allocated_in_fabs_hwm)
60#pragma omp threadprivate(private_total_cells_allocated_in_fabs)
61#pragma omp threadprivate(private_total_cells_allocated_in_fabs_hwm)
62#endif
63
65
92void update_fab_stats (Long n, Long s, std::size_t szt) noexcept;
93
95void BaseFab_Initialize ();
97void BaseFab_Finalize ();
98
100struct SrcComp {
102 explicit SrcComp (int ai) noexcept : i(ai) {}
103 int i;
104};
105
107struct DestComp {
109 explicit DestComp (int ai) noexcept : i(ai) {}
110 int i;
111};
112
114struct NumComps {
116 explicit NumComps (int an) noexcept : n(an) {}
117 int n;
118};
119
120template <typename T>
123makeArray4 (T* p, Box const& bx, int ncomp) noexcept
124{
125 return Array4<T>{p, amrex::begin(bx), amrex::end(bx), ncomp};
126}
127
135template <typename T>
136void
137placementNew (T* const ptr, Long n)
138{
140 if constexpr (!std::is_trivially_default_constructible_v<T>) {
142 {
143 new (ptr+i) T;
144 });
145 }
146}
147
156template <typename T>
157void
158placementDelete (T* const ptr, Long n)
159{
161 if constexpr (!std::is_trivially_destructible_v<T>) {
163 {
164 (ptr+i)->~T();
165 });
166 }
167}
168
219template <class T>
221 : public DataAllocator
222{
223public:
224
225 template <class U> friend class BaseFab;
226
227 using value_type = T;
228
230 BaseFab () noexcept = default;
231
233 explicit BaseFab (Arena* ar) noexcept;
234
236 BaseFab (const Box& bx, int n, Arena* ar);
237
239 explicit BaseFab (const Box& bx, int n = 1, bool alloc = true,
240 bool shared = false, Arena* ar = nullptr);
241
251 BaseFab (const BaseFab<T>& rhs, MakeType make_type, int scomp, int ncomp);
252
258 BaseFab (const Box& bx, int ncomp, T* p);
259 BaseFab (const Box& bx, int ncomp, T const* p);
260
262 explicit BaseFab (Array4<T> const& a) noexcept;
263
265 explicit BaseFab (Array4<T> const& a, IndexType t) noexcept;
266
268 explicit BaseFab (Array4<T const> const& a) noexcept;
269
271 explicit BaseFab (Array4<T const> const& a, IndexType t) noexcept;
272
274 virtual ~BaseFab () noexcept;
275
276 BaseFab (const BaseFab<T>& rhs) = delete;
277 BaseFab<T>& operator= (const BaseFab<T>& rhs) = delete;
278
280 BaseFab (BaseFab<T>&& rhs) noexcept;
282 BaseFab<T>& operator= (BaseFab<T>&& rhs) noexcept;
283
289 template <RunOn run_on AMREX_DEFAULT_RUNON>
290 BaseFab& operator= (T const&) noexcept;
291
293 static void Initialize();
295 static void Finalize();
296
310 void resize (const Box& b, int N = 1, Arena* ar = nullptr);
311
318 template <class U=T>
319 requires (std::is_trivially_destructible_v<U>)
320 [[nodiscard]] Elixir elixir () noexcept;
321
326 void clear ();
327
329 [[nodiscard]] std::unique_ptr<T,DataDeleter> release () noexcept;
330
332 [[nodiscard]] std::size_t nBytes () const noexcept { return this->truesize*sizeof(T); }
333
334 [[nodiscard]] std::size_t nBytesOwned () const noexcept {
335 return (this->ptr_owner) ? nBytes() : 0;
336 }
337
339 [[nodiscard]] std::size_t nBytes (const Box& bx, int ncomps) const noexcept
340 { return bx.numPts() * sizeof(T) * ncomps; }
341
343 [[nodiscard]] int nComp () const noexcept { return this->nvar; }
344
346 [[nodiscard]] const int* nCompPtr() const noexcept {
347 return &(this->nvar);
348 }
349
351 [[nodiscard]] Long numPts () const noexcept { return this->domain.numPts(); }
352
354 [[nodiscard]] Long size () const noexcept { return this->nvar*this->domain.numPts(); }
355
357 [[nodiscard]] const Box& box () const noexcept { return this->domain; }
358
363 [[nodiscard]] IntVect length () const noexcept { return this->domain.length(); }
364
369 [[nodiscard]] const IntVect& smallEnd () const noexcept { return this->domain.smallEnd(); }
370
372 [[nodiscard]] const IntVect& bigEnd () const noexcept { return this->domain.bigEnd(); }
373
382 [[nodiscard]] const int* loVect () const noexcept { return this->domain.loVect(); }
383
392 [[nodiscard]] const int* hiVect () const noexcept { return this->domain.hiVect(); }
393
398 [[nodiscard]] bool contains (const BaseFab<T>& fab) const noexcept
399 {
400 return box().contains(fab.box()) && this->nvar >= fab.nvar;
401 }
402
407 [[nodiscard]] bool contains (const Box& bx) const noexcept { return box().contains(bx); }
408
418 [[nodiscard]] T* dataPtr (int n = 0) noexcept {
419 if (this->dptr) {
420 return &(this->dptr[n*this->domain.numPts()]);
421 } else {
422 return nullptr;
423 }
424 }
425
427 [[nodiscard]] const T* dataPtr (int n = 0) const noexcept {
428 if (this->dptr) {
429 return &(this->dptr[n*this->domain.numPts()]);
430 } else {
431 return nullptr;
432 }
433 }
434
441 [[nodiscard]] T* dataPtr (const IntVect& p, int n = 0) noexcept;
442
449 [[nodiscard]] const T* dataPtr (const IntVect& p, int n = 0) const noexcept;
450
457 void setPtr (T* p, Long sz) noexcept { AMREX_ASSERT(this->dptr == nullptr && this->truesize == 0); this->dptr = p; this->truesize = sz; }
458
463 void prefetchToHost () const noexcept;
468 void prefetchToDevice () const noexcept;
469
474 [[nodiscard]] AMREX_FORCE_INLINE
475 Array4<T const> array () const noexcept
476 {
477 return makeArray4<T const>(this->dptr, this->domain, this->nvar);
478 }
479
486 [[nodiscard]] AMREX_FORCE_INLINE
487 Array4<T const> array (int start_comp) const noexcept
488 {
489 return Array4<T const>(makeArray4<T const>(this->dptr, this->domain, this->nvar),start_comp);
490 }
491
499 [[nodiscard]] AMREX_FORCE_INLINE
500 Array4<T const> array (int start_comp, int num_comps) const noexcept
501 {
502 return Array4<T const>(makeArray4<T const>(this->dptr, this->domain, this->nvar), start_comp, num_comps);
503 }
504
509 [[nodiscard]] AMREX_FORCE_INLINE
510 Array4<T> array () noexcept
511 {
512 return makeArray4<T>(this->dptr, this->domain, this->nvar);
513 }
514
521 [[nodiscard]] AMREX_FORCE_INLINE
522 Array4<T> array (int start_comp) noexcept
523 {
524 return Array4<T>(makeArray4<T>(this->dptr, this->domain, this->nvar),start_comp);
525 }
526
534 [[nodiscard]] AMREX_FORCE_INLINE
535 Array4<T> array (int start_comp, int num_comps) noexcept
536 {
537 return Array4<T>(makeArray4<T>(this->dptr, this->domain, this->nvar), start_comp, num_comps);
538 }
539
546 [[nodiscard]] AMREX_FORCE_INLINE
547 Array4<T const> const_array () const noexcept
548 {
549 return makeArray4<T const>(this->dptr, this->domain, this->nvar);
550 }
551
558 [[nodiscard]] AMREX_FORCE_INLINE
559 Array4<T const> const_array (int start_comp) const noexcept
560 {
561 return Array4<T const>(makeArray4<T const>(this->dptr, this->domain, this->nvar),start_comp);
562 }
563
570 [[nodiscard]] AMREX_FORCE_INLINE
571 Array4<T const> const_array (int start_comp, int num_comps) const noexcept
572 {
573 return Array4<T const>(makeArray4<T const>(this->dptr, this->domain, this->nvar), start_comp, num_comps);
574 }
575
577 [[nodiscard]] bool isAllocated () const noexcept { return this->dptr != nullptr; }
578
585 [[nodiscard]] T& operator() (const IntVect& p, int N) noexcept;
586
588 [[nodiscard]] T& operator() (const IntVect& p) noexcept;
589
591 [[nodiscard]] const T& operator() (const IntVect& p, int N) const noexcept;
592
594 [[nodiscard]] const T& operator() (const IntVect& p) const noexcept;
595
601 void getVal (T* data, const IntVect& pos, int N, int numcomp) const noexcept;
603 void getVal (T* data, const IntVect& pos) const noexcept;
604
605 template <RunOn run_on AMREX_DEFAULT_RUNON>
606 requires (std::is_same_v<T,float> || std::is_same_v<T,double>)
607 void fill_snan () noexcept;
608
615 template <RunOn run_on AMREX_DEFAULT_RUNON>
616 void setVal (T const& x, const Box& bx, int dcomp, int ncomp) noexcept;
618 template <RunOn run_on AMREX_DEFAULT_RUNON>
619 void setVal (T const& x, const Box& bx, int N = 0) noexcept;
621 template <RunOn run_on AMREX_DEFAULT_RUNON>
622 void setVal (T const& x, int N) noexcept;
623
624 template <RunOn run_on AMREX_DEFAULT_RUNON>
625 void setValIfNot (T const& val, const Box& bx, const BaseFab<int>& mask, int nstart, int num) noexcept;
626
632 template <RunOn run_on AMREX_DEFAULT_RUNON>
633 void setComplement (T const& x, const Box& b, int ns, int num) noexcept;
634
651 template <RunOn run_on AMREX_DEFAULT_RUNON>
652 BaseFab<T>& copy (const BaseFab<T>& src, const Box& srcbox, int srccomp,
653 const Box& destbox, int destcomp, int numcomp) noexcept;
654
661 template <RunOn run_on AMREX_DEFAULT_RUNON>
662 BaseFab<T>& copy (const BaseFab<T>& src, int srccomp, int destcomp,
663 int numcomp = 1) noexcept;
670 template <RunOn run_on AMREX_DEFAULT_RUNON>
671 BaseFab<T>& copy (const BaseFab<T>& src, const Box& destbox) noexcept;
672
674 template <RunOn run_on AMREX_DEFAULT_RUNON>
675 std::size_t copyToMem (const Box& srcbox, int srccomp,
676 int numcomp, void* dst) const noexcept;
677
679 template <RunOn run_on AMREX_DEFAULT_RUNON, typename BUF = T>
680 std::size_t copyFromMem (const Box& dstbox, int dstcomp,
681 int numcomp, const void* src) noexcept;
682
684 template <RunOn run_on AMREX_DEFAULT_RUNON, typename BUF = T>
685 std::size_t addFromMem (const Box& dstbox, int dstcomp,
686 int numcomp, const void* src) noexcept;
687
693 BaseFab<T>& shift (const IntVect& v) noexcept;
699 BaseFab<T>& shift (int idir, int n_cell) noexcept;
705 BaseFab<T>& shiftHalf (int dir, int n_cell) noexcept;
711 BaseFab<T>& shiftHalf (const IntVect& v) noexcept;
712
713 template <RunOn run_on AMREX_DEFAULT_RUNON>
714 [[nodiscard]] Real norminfmask (const Box& subbox, const BaseFab<int>& mask, int scomp=0, int ncomp=1) const noexcept;
715
722 template <RunOn run_on AMREX_DEFAULT_RUNON>
723 [[nodiscard]] Real norm (int p, int scomp = 0, int numcomp = 1) const;
724
726 template <RunOn run_on AMREX_DEFAULT_RUNON>
727 [[nodiscard]] Real norm (const Box& subbox, int p, int scomp = 0, int numcomp = 1) const;
729 template <RunOn run_on AMREX_DEFAULT_RUNON>
730 void abs () noexcept;
732 template <RunOn run_on AMREX_DEFAULT_RUNON>
733 void abs (int comp, int numcomp=1) noexcept;
737 template <RunOn run_on AMREX_DEFAULT_RUNON>
738 void abs (const Box& subbox, int comp = 0, int numcomp=1) noexcept;
742 template <RunOn run_on AMREX_DEFAULT_RUNON>
743 [[nodiscard]] T min (int comp = 0) const noexcept;
747 template <RunOn run_on AMREX_DEFAULT_RUNON>
748 [[nodiscard]] T min (const Box& subbox, int comp = 0) const noexcept;
752 template <RunOn run_on AMREX_DEFAULT_RUNON>
753 [[nodiscard]] T max (int comp = 0) const noexcept;
757 template <RunOn run_on AMREX_DEFAULT_RUNON>
758 [[nodiscard]] T max (const Box& subbox, int comp = 0) const noexcept;
762 template <RunOn run_on AMREX_DEFAULT_RUNON>
763 [[nodiscard]] std::pair<T,T> minmax (int comp = 0) const noexcept;
767 template <RunOn run_on AMREX_DEFAULT_RUNON>
768 [[nodiscard]] std::pair<T,T> minmax (const Box& subbox, int comp = 0) const noexcept;
772 template <RunOn run_on AMREX_DEFAULT_RUNON>
773 [[nodiscard]] T maxabs (int comp = 0) const noexcept;
777 template <RunOn run_on AMREX_DEFAULT_RUNON>
778 [[nodiscard]] T maxabs (const Box& subbox, int comp = 0) const noexcept;
779
780 /*(
781 * \return location of given value
782 */
783 template <RunOn run_on AMREX_DEFAULT_RUNON>
784 [[nodiscard]] IntVect indexFromValue (const Box& subbox, int comp, T const& value) const noexcept;
785
789 template <RunOn run_on AMREX_DEFAULT_RUNON>
790 [[nodiscard]] IntVect minIndex (int comp = 0) const noexcept;
795 template <RunOn run_on AMREX_DEFAULT_RUNON>
796 [[nodiscard]] IntVect minIndex (const Box& subbox, int comp = 0) const noexcept;
801 template <RunOn run_on AMREX_DEFAULT_RUNON>
802 void minIndex (const Box& subbox, Real& min_val, IntVect& min_idx, int comp = 0) const noexcept;
803
807 template <RunOn run_on AMREX_DEFAULT_RUNON>
808 [[nodiscard]] IntVect maxIndex (int comp = 0) const noexcept;
813 template <RunOn run_on AMREX_DEFAULT_RUNON>
814 [[nodiscard]] IntVect maxIndex (const Box& subbox, int comp = 0) const noexcept;
819 template <RunOn run_on AMREX_DEFAULT_RUNON>
820 void maxIndex (const Box& subbox, Real& max_value, IntVect& max_idx, int comp = 0) const noexcept;
821
828 template <RunOn run_on AMREX_DEFAULT_RUNON>
829 int maskLT (BaseFab<int>& mask, T const& val, int comp = 0) const noexcept;
831 template <RunOn run_on AMREX_DEFAULT_RUNON>
832 int maskLE (BaseFab<int>& mask, T const& val, int comp = 0) const noexcept;
833
835 template <RunOn run_on AMREX_DEFAULT_RUNON>
836 int maskEQ (BaseFab<int>& mask, T const& val, int comp = 0) const noexcept;
838 template <RunOn run_on AMREX_DEFAULT_RUNON>
839 int maskGT (BaseFab<int>& mask, T const& val, int comp = 0) const noexcept;
841 template <RunOn run_on AMREX_DEFAULT_RUNON>
842 int maskGE (BaseFab<int>& mask, T const& val, int comp = 0) const noexcept;
843
845 template <RunOn run_on AMREX_DEFAULT_RUNON>
846 [[nodiscard]] T sum (int comp, int numcomp = 1) const noexcept;
848 template <RunOn run_on AMREX_DEFAULT_RUNON>
849 [[nodiscard]] T sum (const Box& subbox, int comp, int numcomp = 1) const noexcept;
850
852 template <RunOn run_on AMREX_DEFAULT_RUNON>
853 BaseFab<T>& invert (T const& r, const Box& b, int comp=0, int numcomp=1) noexcept;
855 template <RunOn run_on AMREX_DEFAULT_RUNON>
856 BaseFab<T>& invert (T const& r, int comp, int numcomp=1) noexcept;
857
859 template <RunOn run_on AMREX_DEFAULT_RUNON>
860 BaseFab<T>& negate (const Box& b, int comp=0, int numcomp=1) noexcept;
862 template <RunOn run_on AMREX_DEFAULT_RUNON>
863 BaseFab<T>& negate (int comp, int numcomp=1) noexcept;
864
866 template <RunOn run_on AMREX_DEFAULT_RUNON>
867 BaseFab<T>& plus (T const& r, const Box& b, int comp=0, int numcomp=1) noexcept;
868
870 template <RunOn run_on AMREX_DEFAULT_RUNON>
871 BaseFab<T>& plus (T const& r, int comp, int numcomp=1) noexcept;
872
878 template <RunOn run_on AMREX_DEFAULT_RUNON>
879 BaseFab<T>& plus (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp=1) noexcept;
885 template <RunOn run_on AMREX_DEFAULT_RUNON>
886 BaseFab<T>& plus (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp, int numcomp=1) noexcept;
891 template <RunOn run_on AMREX_DEFAULT_RUNON>
892 BaseFab<T>& plus (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
893 int srccomp, int destcomp, int numcomp=1) noexcept;
894
896 template <RunOn run_on AMREX_DEFAULT_RUNON>
897 BaseFab<T>& atomicAdd (const BaseFab<T>& x) noexcept;
898
904 template <RunOn run_on AMREX_DEFAULT_RUNON>
905 BaseFab<T>& atomicAdd (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp=1) noexcept;
911 template <RunOn run_on AMREX_DEFAULT_RUNON>
912 BaseFab<T>& atomicAdd (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
913 int numcomp=1) noexcept;
918 template <RunOn run_on AMREX_DEFAULT_RUNON>
919 BaseFab<T>& atomicAdd (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
920 int srccomp, int destcomp, int numcomp=1) noexcept;
921
927 template <RunOn run_on AMREX_DEFAULT_RUNON>
928 BaseFab<T>& lockAdd (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
929 int srccomp, int destcomp, int numcomp) noexcept;
930
932 template <RunOn run_on AMREX_DEFAULT_RUNON>
933 BaseFab<T>& saxpy (T a, const BaseFab<T>& x, const Box& srcbox, const Box& destbox,
934 int srccomp, int destcomp, int numcomp=1) noexcept;
936 template <RunOn run_on AMREX_DEFAULT_RUNON>
937 BaseFab<T>& saxpy (T a, const BaseFab<T>& x) noexcept;
938
940 template <RunOn run_on AMREX_DEFAULT_RUNON>
941 BaseFab<T>& xpay (T a, const BaseFab<T>& x, const Box& srcbox, const Box& destbox,
942 int srccomp, int destcomp, int numcomp=1) noexcept;
943
945 template <RunOn run_on AMREX_DEFAULT_RUNON>
946 BaseFab<T>& addproduct (const Box& destbox, int destcomp, int numcomp,
947 const BaseFab<T>& src1, int comp1,
948 const BaseFab<T>& src2, int comp2) noexcept;
949
955 template <RunOn run_on AMREX_DEFAULT_RUNON>
956 BaseFab<T>& minus (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp=1) noexcept;
962 template <RunOn run_on AMREX_DEFAULT_RUNON>
963 BaseFab<T>& minus (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
964 int numcomp=1) noexcept;
969 template <RunOn run_on AMREX_DEFAULT_RUNON>
970 BaseFab<T>& minus (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
971 int srccomp, int destcomp, int numcomp=1) noexcept;
972
974 template <RunOn run_on AMREX_DEFAULT_RUNON>
975 BaseFab<T>& mult (T const& r, int comp, int numcomp=1) noexcept;
977 template <RunOn run_on AMREX_DEFAULT_RUNON>
978 BaseFab<T>& mult (T const& r, const Box& b, int comp=0, int numcomp=1) noexcept;
979
985 template <RunOn run_on AMREX_DEFAULT_RUNON>
986 BaseFab<T>& mult (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp=1) noexcept;
987
993 template <RunOn run_on AMREX_DEFAULT_RUNON>
994 BaseFab<T>& mult (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
995 int numcomp=1) noexcept;
996
1001 template <RunOn run_on AMREX_DEFAULT_RUNON>
1002 BaseFab<T>& mult (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
1003 int srccomp, int destcomp, int numcomp=1) noexcept;
1004
1006 template <RunOn run_on AMREX_DEFAULT_RUNON>
1007 BaseFab<T>& divide (T const& r, int comp, int numcomp=1) noexcept;
1008
1010 template <RunOn run_on AMREX_DEFAULT_RUNON>
1011 BaseFab<T>& divide (T const& r, const Box& b, int comp=0, int numcomp=1) noexcept;
1012
1019 template <RunOn run_on AMREX_DEFAULT_RUNON>
1020 BaseFab<T>& divide (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp=1) noexcept;
1026 template <RunOn run_on AMREX_DEFAULT_RUNON>
1027 BaseFab<T>& divide (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
1028 int numcomp=1) noexcept;
1033 template <RunOn run_on AMREX_DEFAULT_RUNON>
1034 BaseFab<T>& divide (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
1035 int srccomp, int destcomp, int numcomp=1) noexcept;
1039 template <RunOn run_on AMREX_DEFAULT_RUNON>
1040 BaseFab<T>& protected_divide (const BaseFab<T>& src) noexcept;
1041
1049 template <RunOn run_on AMREX_DEFAULT_RUNON>
1050 BaseFab<T>& protected_divide (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp=1) noexcept;
1051
1058 template <RunOn run_on AMREX_DEFAULT_RUNON>
1059 BaseFab<T>& protected_divide (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
1060 int numcomp=1) noexcept;
1061
1067 template <RunOn run_on AMREX_DEFAULT_RUNON>
1068 BaseFab<T>& protected_divide (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
1069 int srccomp, int destcomp, int numcomp=1) noexcept;
1070
1081 template <RunOn run_on AMREX_DEFAULT_RUNON>
1082 BaseFab<T>& linInterp (const BaseFab<T>& f1, const Box& b1, int comp1,
1083 const BaseFab<T>& f2, const Box& b2, int comp2,
1084 Real t1, Real t2, Real t,
1085 const Box& b, int comp, int numcomp = 1) noexcept;
1086
1088 template <RunOn run_on AMREX_DEFAULT_RUNON>
1089 BaseFab<T>& linInterp (const BaseFab<T>& f1, int comp1,
1090 const BaseFab<T>& f2, int comp2,
1091 Real t1, Real t2, Real t,
1092 const Box& b, int comp, int numcomp = 1) noexcept;
1093
1103 template <RunOn run_on AMREX_DEFAULT_RUNON>
1104 BaseFab<T>& linComb (const BaseFab<T>& f1, const Box& b1, int comp1,
1105 const BaseFab<T>& f2, const Box& b2, int comp2,
1106 Real alpha, Real beta, const Box& b,
1107 int comp, int numcomp = 1) noexcept;
1108
1110 template <RunOn run_on AMREX_DEFAULT_RUNON>
1111 [[nodiscard]] T dot (const Box& xbx, int xcomp, const BaseFab<T>& y, const Box& ybx, int ycomp,
1112 int numcomp = 1) const noexcept;
1113
1114 template <RunOn run_on AMREX_DEFAULT_RUNON>
1115 [[nodiscard]] T dotmask (const BaseFab<int>& mask, const Box& xbx, int xcomp,
1116 const BaseFab<T>& y, const Box& ybx, int ycomp,
1117 int numcomp) const noexcept;
1118
1120 void SetBoxType (const IndexType& typ) noexcept { this->domain.setType(typ); }
1121
1122 //
1123 // New interfaces
1124 //
1125
1127 template <RunOn run_on AMREX_DEFAULT_RUNON>
1128 void setVal (T const& val) noexcept;
1129 //
1131 template <RunOn run_on AMREX_DEFAULT_RUNON>
1132 void setVal (T const& x, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept;
1133
1134 template <RunOn run_on AMREX_DEFAULT_RUNON>
1135 void setValIf (T const& val, const BaseFab<int>& mask) noexcept;
1136 //
1138 template <RunOn run_on AMREX_DEFAULT_RUNON>
1139 void setValIf (T const& val, Box const& bx, const BaseFab<int>& mask, DestComp dcomp, NumComps ncomp) noexcept;
1140
1141 template <RunOn run_on AMREX_DEFAULT_RUNON>
1142 void setValIfNot (T const& val, const BaseFab<int>& mask) noexcept;
1143 //
1145 template <RunOn run_on AMREX_DEFAULT_RUNON>
1146 void setValIfNot (T const& val, Box const& bx, const BaseFab<int>& mask, DestComp dcomp, NumComps ncomp) noexcept;
1147
1149 template <RunOn run_on AMREX_DEFAULT_RUNON>
1150 void setComplement (T const& x, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept;
1151
1157 template <RunOn run_on AMREX_DEFAULT_RUNON>
1158 BaseFab<T>& copy (const BaseFab<T>& src) noexcept;
1159 //
1161 template <RunOn run_on AMREX_DEFAULT_RUNON>
1162 BaseFab<T>& copy (const BaseFab<T>& src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept;
1163
1165 template <RunOn run_on AMREX_DEFAULT_RUNON>
1166 BaseFab<T>& plus (T const& val) noexcept;
1167 //
1168 template <RunOn run_on AMREX_DEFAULT_RUNON>
1169 BaseFab<T>& operator+= (T const& val) noexcept;
1170 //
1172 template <RunOn run_on AMREX_DEFAULT_RUNON>
1173 BaseFab<T>& plus (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept;
1179 template <RunOn run_on AMREX_DEFAULT_RUNON>
1180 BaseFab<T>& plus (const BaseFab<T>& src) noexcept;
1181 //
1182 template <RunOn run_on AMREX_DEFAULT_RUNON>
1183 BaseFab<T>& operator+= (const BaseFab<T>& src) noexcept;
1184 //
1186 template <RunOn run_on AMREX_DEFAULT_RUNON>
1187 BaseFab<T>& plus (const BaseFab<T>& src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept;
1188
1190 template <RunOn run_on AMREX_DEFAULT_RUNON>
1191 BaseFab<T>& minus (T const& val) noexcept;
1192 //
1193 template <RunOn run_on AMREX_DEFAULT_RUNON>
1194 BaseFab<T>& operator-= (T const& val) noexcept;
1195 //
1197 template <RunOn run_on AMREX_DEFAULT_RUNON>
1198 BaseFab<T>& minus (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept;
1204 template <RunOn run_on AMREX_DEFAULT_RUNON>
1205 BaseFab<T>& minus (const BaseFab<T>& src) noexcept;
1206 //
1207 template <RunOn run_on AMREX_DEFAULT_RUNON>
1208 BaseFab<T>& operator-= (const BaseFab<T>& src) noexcept;
1209 //
1211 template <RunOn run_on AMREX_DEFAULT_RUNON>
1212 BaseFab<T>& minus (const BaseFab<T>& src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept;
1213
1215 template <RunOn run_on AMREX_DEFAULT_RUNON>
1216 BaseFab<T>& mult (T const& val) noexcept;
1217 //
1218 template <RunOn run_on AMREX_DEFAULT_RUNON>
1219 BaseFab<T>& operator*= (T const& val) noexcept;
1220 //
1222 template <RunOn run_on AMREX_DEFAULT_RUNON>
1223 BaseFab<T>& mult (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept;
1229 template <RunOn run_on AMREX_DEFAULT_RUNON>
1230 BaseFab<T>& mult (const BaseFab<T>& src) noexcept;
1231 //
1232 template <RunOn run_on AMREX_DEFAULT_RUNON>
1233 BaseFab<T>& operator*= (const BaseFab<T>& src) noexcept;
1234 //
1236 template <RunOn run_on AMREX_DEFAULT_RUNON>
1237 BaseFab<T>& mult (const BaseFab<T>& src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept;
1238
1240 template <RunOn run_on AMREX_DEFAULT_RUNON>
1241 BaseFab<T>& divide (T const& val) noexcept;
1242 //
1243 template <RunOn run_on AMREX_DEFAULT_RUNON>
1244 BaseFab<T>& operator/= (T const& val) noexcept;
1245 //
1247 template <RunOn run_on AMREX_DEFAULT_RUNON>
1248 BaseFab<T>& divide (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept;
1254 template <RunOn run_on AMREX_DEFAULT_RUNON>
1255 BaseFab<T>& divide (const BaseFab<T>& src) noexcept;
1256 //
1257 template <RunOn run_on AMREX_DEFAULT_RUNON>
1258 BaseFab<T>& operator/= (const BaseFab<T>& src) noexcept;
1259 //
1261 template <RunOn run_on AMREX_DEFAULT_RUNON>
1262 BaseFab<T>& divide (const BaseFab<T>& src, Box bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept;
1263
1265 template <RunOn run_on AMREX_DEFAULT_RUNON>
1266 BaseFab<T>& negate () noexcept;
1267 //
1268 template <RunOn run_on AMREX_DEFAULT_RUNON>
1269 BaseFab<T>& negate (const Box& bx, DestComp dcomp, NumComps ncomp) noexcept;
1270
1272 template <RunOn run_on AMREX_DEFAULT_RUNON>
1273 BaseFab<T>& invert (T const& r) noexcept;
1274 //
1275 template <RunOn run_on AMREX_DEFAULT_RUNON>
1276 BaseFab<T>& invert (T const& r, const Box& bx, DestComp dcomp, NumComps ncomp) noexcept;
1277
1279 template <RunOn run_on AMREX_DEFAULT_RUNON>
1280 [[nodiscard]] T sum (const Box& bx, DestComp dcomp, NumComps ncomp) const noexcept;
1281
1283 template <RunOn run_on AMREX_DEFAULT_RUNON>
1284 [[nodiscard]] T dot (const BaseFab<T>& src, const Box& bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) const noexcept;
1285
1287 template <RunOn run_on AMREX_DEFAULT_RUNON>
1288 [[nodiscard]] T dot (const Box& bx, int destcomp, int numcomp) const noexcept;
1289
1291 template <RunOn run_on AMREX_DEFAULT_RUNON>
1292 [[nodiscard]] T dot (const Box& bx, DestComp dcomp, NumComps ncomp) const noexcept;
1293
1295 template <RunOn run_on AMREX_DEFAULT_RUNON>
1296 [[nodiscard]] T dotmask (const BaseFab<T>& src, const Box& bx, const BaseFab<int>& mask,
1297 SrcComp scomp, DestComp dcomp, NumComps ncomp) const noexcept;
1298
1299protected:
1301 void define ();
1302
1303 T* dptr = nullptr;
1305 int nvar = 0;
1307 bool ptr_owner = false;
1308 bool shared_memory = false;
1309#ifdef AMREX_USE_GPU
1311#endif
1312};
1313
1314template <class T>
1316T*
1317BaseFab<T>::dataPtr (const IntVect& p, int n) noexcept
1318{
1319 AMREX_ASSERT(n >= 0);
1320 AMREX_ASSERT(n < this->nvar);
1321 AMREX_ASSERT(!(this->dptr == nullptr));
1322 AMREX_ASSERT(this->domain.contains(p));
1323
1324 return this->dptr + (this->domain.index(p)+n*this->domain.numPts());
1325}
1326
1327template <class T>
1329const T*
1330BaseFab<T>::dataPtr (const IntVect& p, int n) const noexcept
1331{
1332 AMREX_ASSERT(n >= 0);
1333 AMREX_ASSERT(n < this->nvar);
1334 AMREX_ASSERT(!(this->dptr == nullptr));
1335 AMREX_ASSERT(this->domain.contains(p));
1336
1337 return this->dptr + (this->domain.index(p)+n*this->domain.numPts());
1338}
1339
1340template <class T>
1341void
1343{
1344#ifdef AMREX_USE_GPU
1345 if (this->arena()->isManaged()) {
1346#if defined(AMREX_USE_SYCL)
1347 // xxxxx SYCL todo: prefetchToHost
1348 // std::size_t s = sizeof(T)*this->nvar*this->domain.numPts();
1349 // auto& q = Gpu::Device::streamQueue();
1350 // q.submit([&] (sycl::handler& h) { h.prefetch(this->dptr, s); });
1351#elif defined(AMREX_USE_CUDA) && !defined(_WIN32)
1352 if (Gpu::Device::devicePropMajor() >= 6) {
1353 std::size_t s = sizeof(T)*this->nvar*this->domain.numPts();
1354#if defined(CUDART_VERSION) && (CUDART_VERSION >= 13000)
1355 cudaMemLocation location = {};
1356 location.type = cudaMemLocationTypeHost;
1357 AMREX_CUDA_SAFE_CALL(cudaMemPrefetchAsync(this->dptr, s, location, 0,
1358 Gpu::gpuStream()));
1359#else
1360 AMREX_CUDA_SAFE_CALL(cudaMemPrefetchAsync(this->dptr, s,
1361 cudaCpuDeviceId,
1362 Gpu::gpuStream()));
1363#endif
1364 }
1365#elif defined(AMREX_USE_HIP)
1366 // xxxxx HIP FIX HERE after managed memory is supported
1367#endif
1368 }
1369#endif
1370}
1371
1372template <class T>
1373void
1375{
1376#ifdef AMREX_USE_GPU
1377 if (this->arena()->isManaged()) {
1378#if defined(AMREX_USE_SYCL)
1379 std::size_t s = sizeof(T)*this->nvar*this->domain.numPts();
1380 auto& q = Gpu::Device::streamQueue();
1381 q.submit([&] (sycl::handler& h) { h.prefetch(this->dptr, s); });
1382#elif defined(AMREX_USE_CUDA) && !defined(_WIN32)
1383 if (Gpu::Device::devicePropMajor() >= 6) {
1384 std::size_t s = sizeof(T)*this->nvar*this->domain.numPts();
1385#if defined(CUDART_VERSION) && (CUDART_VERSION >= 13000)
1386 cudaMemLocation location = {};
1387 location.type = cudaMemLocationTypeDevice;
1388 location.id = Gpu::Device::deviceId();
1389 AMREX_CUDA_SAFE_CALL(cudaMemPrefetchAsync(this->dptr, s, location, 0,
1390 Gpu::gpuStream()));
1391#else
1392 AMREX_CUDA_SAFE_CALL(cudaMemPrefetchAsync(this->dptr, s,
1394 Gpu::gpuStream()));
1395#endif
1396 }
1397#elif defined(AMREX_USE_HIP)
1398 // xxxxx HIP FIX HERE after managed memory is supported
1399#endif
1400 }
1401#endif
1402}
1403
1404template <class T>
1406T&
1407BaseFab<T>::operator() (const IntVect& p, int n) noexcept
1408{
1409 AMREX_ASSERT(n >= 0);
1410 AMREX_ASSERT(n < this->nvar);
1411 AMREX_ASSERT(!(this->dptr == nullptr));
1412 AMREX_ASSERT(this->domain.contains(p));
1413
1414 return this->dptr[this->domain.index(p)+n*this->domain.numPts()];
1415}
1416
1417template <class T>
1419T&
1421{
1422 AMREX_ASSERT(!(this->dptr == nullptr));
1423 AMREX_ASSERT(this->domain.contains(p));
1424
1425 return this->dptr[this->domain.index(p)];
1426}
1427
1428template <class T>
1430const T&
1431BaseFab<T>::operator() (const IntVect& p, int n) const noexcept
1432{
1433 AMREX_ASSERT(n >= 0);
1434 AMREX_ASSERT(n < this->nvar);
1435 AMREX_ASSERT(!(this->dptr == nullptr));
1436 AMREX_ASSERT(this->domain.contains(p));
1437
1438 return this->dptr[this->domain.index(p)+n*this->domain.numPts()];
1439}
1440
1441template <class T>
1443const T&
1444BaseFab<T>::operator() (const IntVect& p) const noexcept
1445{
1446 AMREX_ASSERT(!(this->dptr == nullptr));
1447 AMREX_ASSERT(this->domain.contains(p));
1448
1449 return this->dptr[this->domain.index(p)];
1450}
1451
1452template <class T>
1453void
1455 const IntVect& pos,
1456 int n,
1457 int numcomp) const noexcept
1458{
1459 const int loc = this->domain.index(pos);
1460 const Long sz = this->domain.numPts();
1461
1462 AMREX_ASSERT(!(this->dptr == nullptr));
1463 AMREX_ASSERT(n >= 0 && n + numcomp <= this->nvar);
1464
1465 for (int k = 0; k < numcomp; k++) {
1466 data[k] = this->dptr[loc+(n+k)*sz];
1467 }
1468}
1469
1470template <class T>
1471void
1473 const IntVect& pos) const noexcept
1474{
1475 getVal(data,pos,0,this->nvar);
1476}
1477
1478template <class T>
1480BaseFab<T>::shift (const IntVect& v) noexcept
1481{
1482 this->domain += v;
1483 return *this;
1484}
1485
1486template <class T>
1488BaseFab<T>::shift (int idir, int n_cell) noexcept
1489{
1490 this->domain.shift(idir,n_cell);
1491 return *this;
1492}
1493
1494template <class T>
1495BaseFab<T> &
1497{
1498 this->domain.shiftHalf(v);
1499 return *this;
1500}
1501
1502template <class T>
1503BaseFab<T> &
1504BaseFab<T>::shiftHalf (int idir, int n_cell) noexcept
1505{
1506 this->domain.shiftHalf(idir,n_cell);
1507 return *this;
1508}
1509
1510template <class T>
1511template <RunOn run_on>
1512requires (std::is_same_v<T,float> || std::is_same_v<T,double>)
1513void
1515{
1516 amrex::fill_snan<run_on>(this->dptr, this->truesize);
1517}
1518
1519template <class T>
1520template <RunOn run_on>
1521void
1522BaseFab<T>::setVal (T const& x, const Box& bx, int n) noexcept
1523{
1524 this->setVal<run_on>(x, bx, DestComp{n}, NumComps{1});
1525}
1526
1527template <class T>
1528template <RunOn run_on>
1529void
1530BaseFab<T>::setVal (T const& x, int n) noexcept
1531{
1532 this->setVal<run_on>(x, this->domain, DestComp{n}, NumComps{1});
1533}
1534
1535template <class T>
1536template <RunOn run_on>
1537void
1538BaseFab<T>::setVal (T const& x, const Box& bx, int dcomp, int ncomp) noexcept
1539{
1540 this->setVal<run_on>(x, bx, DestComp{dcomp}, NumComps{ncomp});
1541}
1542
1543template <class T>
1544template <RunOn run_on>
1545void
1546BaseFab<T>::setValIfNot (T const& val, const Box& bx, const BaseFab<int>& mask, int ns, int num) noexcept
1547{
1548 this->setValIfNot<run_on>(val, bx, mask, DestComp{ns}, NumComps{num});
1549}
1550
1551template <class T>
1552template <RunOn run_on>
1554BaseFab<T>::copy (const BaseFab<T>& src, const Box& srcbox, int srccomp,
1555 const Box& destbox, int destcomp, int numcomp) noexcept
1556{
1557 AMREX_ASSERT(destbox.ok());
1558 AMREX_ASSERT(srcbox.sameSize(destbox));
1559 AMREX_ASSERT(src.box().contains(srcbox));
1560 AMREX_ASSERT(this->domain.contains(destbox));
1561 AMREX_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
1562 AMREX_ASSERT(destcomp >= 0 && destcomp+numcomp <= this->nvar);
1563
1564 Array4<T> const& d = this->array();
1565 Array4<T const> const& s = src.const_array();
1566 const auto dlo = amrex::lbound(destbox);
1567 const auto slo = amrex::lbound(srcbox);
1568 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
1569
1570 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
1571 {
1572 d(i,j,k,n+destcomp) = s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
1573 });
1574
1575 return *this;
1576}
1577
1578template <class T>
1579template <RunOn run_on>
1581BaseFab<T>::copy (const BaseFab<T>& src, const Box& destbox) noexcept
1582{
1583 return this->copy<run_on>(src, destbox, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
1584}
1585
1586template <class T>
1587template <RunOn run_on>
1589BaseFab<T>::copy (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
1590{
1591 return copy<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
1592}
1593
1594template <class T>
1595void
1597{
1598 AMREX_ASSERT(this->dptr == nullptr);
1599 AMREX_ASSERT(this->domain.numPts() > 0);
1600 AMREX_ASSERT(this->nvar >= 0);
1601 if (this->nvar == 0) { return; }
1602 AMREX_ASSERT(std::numeric_limits<Long>::max()/this->nvar > this->domain.numPts());
1603
1604 this->truesize = this->nvar*this->domain.numPts();
1605 this->ptr_owner = true;
1606 this->dptr = static_cast<T*>(this->alloc(this->truesize*sizeof(T)));
1607#ifdef AMREX_USE_GPU
1608 this->alloc_stream = Gpu::gpuStream();
1609#endif
1610
1611 placementNew(this->dptr, this->truesize);
1612
1613 amrex::update_fab_stats(this->domain.numPts(), this->truesize, sizeof(T));
1614
1615 if constexpr (std::is_same_v<T,float> || std::is_same_v<T,double>) {
1616 if (amrex::InitSNaN() && this->truesize > 0) {
1617#ifdef AMREX_USE_GPU
1618 if (Gpu::inLaunchRegion() && arena()->isDeviceAccessible()) {
1619 this->template fill_snan<RunOn::Device>();
1621 } else
1622#endif
1623 {
1624 this->template fill_snan<RunOn::Host>();
1625 }
1626 }
1627 }
1628}
1629
1630template <class T>
1632 : DataAllocator{ar}
1633{}
1634
1635template <class T>
1636BaseFab<T>::BaseFab (const Box& bx, int n, Arena* ar)
1637 : DataAllocator{ar}, domain(bx), nvar(n)
1638{
1639 define();
1640}
1641
1642template <class T>
1643BaseFab<T>::BaseFab (const Box& bx, int n, bool alloc, bool shared, Arena* ar)
1644 : DataAllocator{ar}, domain(bx), nvar(n), shared_memory(shared)
1645{
1646 if (!this->shared_memory && alloc) { define(); }
1647}
1648
1649template <class T>
1650BaseFab<T>::BaseFab (const BaseFab<T>& rhs, MakeType make_type, int scomp, int ncomp)
1651 : DataAllocator{rhs.arena()},
1652 dptr(const_cast<T*>(rhs.dataPtr(scomp))),
1653 domain(rhs.domain), nvar(ncomp),
1654 truesize(ncomp*rhs.domain.numPts())
1655{
1656 AMREX_ASSERT(scomp+ncomp <= rhs.nComp());
1657 if (make_type == amrex::make_deep_copy)
1658 {
1659 this->dptr = nullptr;
1660 define();
1661 this->copy<RunOn::Device>(rhs, this->domain, scomp, this->domain, 0, ncomp);
1662 } else if (make_type == amrex::make_alias) {
1663 ; // nothing to do
1664 } else {
1665 amrex::Abort("BaseFab: unknown MakeType");
1666 }
1667}
1668
1669template<class T>
1670BaseFab<T>::BaseFab (const Box& bx, int ncomp, T* p)
1671 : dptr(p), domain(bx), nvar(ncomp), truesize(bx.numPts()*ncomp)
1672{
1673}
1674
1675template<class T>
1676BaseFab<T>::BaseFab (const Box& bx, int ncomp, T const* p)
1677 : dptr(const_cast<T*>(p)), domain(bx), nvar(ncomp), truesize(bx.numPts()*ncomp)
1678{
1679}
1680
1681template<class T>
1683 : dptr(a.p),
1684 domain(IntVect(AMREX_D_DECL(a.begin[0],a.begin[1],a.begin[2])),
1685 IntVect(AMREX_D_DECL(a.end[0]-1,a.end[1]-1,a.end[2]-1))),
1686 nvar(a.nComp()), truesize(a.size())
1687{}
1688
1689template<class T>
1691 : dptr(a.p),
1692 domain(IntVect(AMREX_D_DECL(a.begin[0],a.begin[1],a.begin[2])),
1693 IntVect(AMREX_D_DECL(a.end[0]-1,a.end[1]-1,a.end[2]-1)), t),
1694 nvar(a.nComp()), truesize(a.size())
1695{}
1696
1697template<class T>
1699 : dptr(const_cast<T*>(a.p)),
1700 domain(IntVect(AMREX_D_DECL(a.begin[0],a.begin[1],a.begin[2])),
1701 IntVect(AMREX_D_DECL(a.end[0]-1,a.end[1]-1,a.end[2]-1))),
1702 nvar(a.nComp()), truesize(a.size())
1703{}
1704
1705template<class T>
1707 : dptr(const_cast<T*>(a.p)),
1708 domain(IntVect(AMREX_D_DECL(a.begin[0],a.begin[1],a.begin[2])),
1709 IntVect(AMREX_D_DECL(a.end[0]-1,a.end[1]-1,a.end[2]-1)), t),
1710 nvar(a.nComp()), truesize(a.size())
1711{}
1712
1713template <class T>
1715{
1716 clear();
1717}
1718
1719template <class T>
1721 : DataAllocator{rhs.arena()},
1722 dptr(rhs.dptr), domain(rhs.domain),
1723 nvar(rhs.nvar), truesize(rhs.truesize),
1724 ptr_owner(rhs.ptr_owner), shared_memory(rhs.shared_memory)
1725#ifdef AMREX_USE_GPU
1726 , alloc_stream(rhs.alloc_stream)
1727#endif
1728{
1729 rhs.dptr = nullptr;
1730 rhs.ptr_owner = false;
1731}
1732
1733template <class T>
1734BaseFab<T>&
1736{
1737 if (this != &rhs) {
1738 clear();
1739 DataAllocator::operator=(rhs);
1740 dptr = rhs.dptr;
1741 domain = rhs.domain;
1742 nvar = rhs.nvar;
1743 truesize = rhs.truesize;
1744 ptr_owner = rhs.ptr_owner;
1745 shared_memory = rhs.shared_memory;
1746#ifdef AMREX_USE_GPU
1747 alloc_stream = rhs.alloc_stream;
1748#endif
1749
1750 rhs.dptr = nullptr;
1751 rhs.ptr_owner = false;
1752 }
1753 return *this;
1754}
1755
1756template <class T>
1757template <RunOn run_on>
1759BaseFab<T>::operator= (T const& t) noexcept
1760{
1761 setVal<run_on>(t);
1762 return *this;
1763}
1764
1765template <class T>
1766void
1767BaseFab<T>::resize (const Box& b, int n, Arena* ar)
1768{
1769 this->nvar = n;
1770 this->domain = b;
1771
1772 if (ar == nullptr) {
1773 ar = m_arena;
1774 }
1775
1776 if (arena() != DataAllocator(ar).arena()) {
1777 clear();
1778 m_arena = ar;
1779 define();
1780 }
1781 else if (this->dptr == nullptr || !this->ptr_owner)
1782 {
1783 if (this->shared_memory) {
1784 amrex::Abort("BaseFab::resize: BaseFab in shared memory cannot increase size");
1785 }
1786
1787 this->dptr = nullptr;
1788 define();
1789 }
1790 else if (this->nvar*this->domain.numPts() > this->truesize
1791#ifdef AMREX_USE_GPU
1792 || (arena()->isStreamOrderedArena() && alloc_stream != Gpu::gpuStream())
1793#endif
1794 )
1795 {
1796 if (this->shared_memory) {
1797 amrex::Abort("BaseFab::resize: BaseFab in shared memory cannot increase size");
1798 }
1799
1800 clear();
1801
1802 define();
1803 }
1804}
1805
1806template <class T>
1807template <class U>
1808requires (std::is_trivially_destructible_v<U>)
1809Elixir
1811{
1812 bool o;
1813 if (Gpu::inLaunchRegion()) {
1814 o = this->ptr_owner;
1815 this->ptr_owner = false;
1816 if (o && this->dptr) {
1817 if (this->nvar > 1) {
1818 amrex::update_fab_stats(-this->truesize/this->nvar, -this->truesize, sizeof(T));
1819 } else {
1820 amrex::update_fab_stats(0, -this->truesize, sizeof(T));
1821 }
1822 }
1823 } else {
1824 o = false;
1825 }
1826 return Elixir((o ? this->dptr : nullptr), this->arena());
1827}
1828
1829template <class T>
1830void
1832{
1833 if (this->dptr)
1834 {
1835 //
1836 // Call T::~T() on the to-be-destroyed memory.
1837 //
1838 if (this->ptr_owner)
1839 {
1840 if (this->shared_memory)
1841 {
1842 amrex::Abort("BaseFab::clear: BaseFab cannot be owner of shared memory");
1843 }
1844
1845 placementDelete(this->dptr, this->truesize);
1846
1847#ifdef AMREX_USE_GPU
1848 this->arena()->streamOrderedFree(this->dptr, alloc_stream);
1849#else
1850 this->free(this->dptr);
1851#endif
1852
1853 if (this->nvar > 1) {
1854 amrex::update_fab_stats(-this->truesize/this->nvar, -this->truesize, sizeof(T));
1855 } else {
1856 amrex::update_fab_stats(0, -this->truesize, sizeof(T));
1857 }
1858 }
1859
1860 this->dptr = nullptr;
1861 this->truesize = 0;
1862 }
1863}
1864
1865template <class T>
1866std::unique_ptr<T,DataDeleter>
1868{
1869 std::unique_ptr<T,DataDeleter> r(nullptr, DataDeleter{this->arena()});
1870 if (this->dptr && this->ptr_owner) {
1871 r.reset(this->dptr);
1872 this->ptr_owner = false;
1873 if (this->nvar > 1) {
1874 amrex::update_fab_stats(-this->truesize/this->nvar, -this->truesize, sizeof(T));
1875 } else {
1876 amrex::update_fab_stats(0, -this->truesize, sizeof(T));
1877 }
1878 }
1879 return r;
1880}
1881
1882template <class T>
1883template <RunOn run_on>
1884std::size_t
1886 int srccomp,
1887 int numcomp,
1888 void* dst) const noexcept
1889{
1890 BL_ASSERT(box().contains(srcbox));
1891 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= nComp());
1892
1893 if (srcbox.ok())
1894 {
1895 Array4<T> d(static_cast<T*>(dst),amrex::begin(srcbox),amrex::end(srcbox),numcomp);
1896 Array4<T const> const& s = this->const_array();
1897 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, srcbox, numcomp, i, j, k, n,
1898 {
1899 d(i,j,k,n) = s(i,j,k,n+srccomp);
1900 });
1901 return sizeof(T)*d.size();
1902 }
1903 else
1904 {
1905 return 0;
1906 }
1907}
1908
1909template <class T>
1910template <RunOn run_on, typename BUF>
1911std::size_t
1913 int dstcomp,
1914 int numcomp,
1915 const void* src) noexcept
1916{
1917 BL_ASSERT(box().contains(dstbox));
1918 BL_ASSERT(dstcomp >= 0 && dstcomp+numcomp <= nComp());
1919
1920 if (dstbox.ok())
1921 {
1922 Array4<BUF const> s(static_cast<BUF const*>(src), amrex::begin(dstbox),
1923 amrex::end(dstbox), numcomp);
1924 Array4<T> const& d = this->array();
1925 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, dstbox, numcomp, i, j, k, n,
1926 {
1927 d(i,j,k,n+dstcomp) = static_cast<T>(s(i,j,k,n));
1928 });
1929 return sizeof(BUF)*s.size();
1930 }
1931 else
1932 {
1933 return 0;
1934 }
1935}
1936
1937template <class T>
1938template <RunOn run_on, typename BUF>
1939std::size_t
1941 int dstcomp,
1942 int numcomp,
1943 const void* src) noexcept
1944{
1945 BL_ASSERT(box().contains(dstbox));
1946 BL_ASSERT(dstcomp >= 0 && dstcomp+numcomp <= nComp());
1947
1948 if (dstbox.ok())
1949 {
1950 Array4<BUF const> s(static_cast<BUF const*>(src), amrex::begin(dstbox),
1951 amrex::end(dstbox), numcomp);
1952 Array4<T> const& d = this->array();
1953 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, dstbox, numcomp, i, j, k, n,
1954 {
1955 d(i,j,k,n+dstcomp) += static_cast<T>(s(i,j,k,n));
1956 });
1957 return sizeof(BUF)*s.size();
1958 }
1959 else
1960 {
1961 return 0;
1962 }
1963}
1964
1965template <class T>
1966template <RunOn run_on>
1967void
1968BaseFab<T>::setComplement (T const& x, const Box& b, int ns, int num) noexcept
1969{
1970 this->setComplement<run_on>(x, b, DestComp{ns}, NumComps{num});
1971}
1972
1973template <class T>
1974template <RunOn run_on>
1975void
1977{
1978 this->abs<run_on>(this->domain,0,this->nvar);
1979}
1980
1981template <class T>
1982template <RunOn run_on>
1983void
1984BaseFab<T>::abs (int comp, int numcomp) noexcept
1985{
1986 this->abs<run_on>(this->domain,comp,numcomp);
1987}
1988
1989template <class T>
1990template <RunOn run_on>
1991void
1992BaseFab<T>::abs (const Box& subbox, int comp, int numcomp) noexcept
1993{
1994 Array4<T> const& a = this->array();
1995 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG (run_on, subbox, numcomp, i, j, k, n,
1996 {
1997 a(i,j,k,n+comp) = std::abs(a(i,j,k,n+comp));
1998 });
1999}
2000
2001template <class T>
2002template <RunOn run_on>
2003Real
2005 int scomp, int ncomp) const noexcept
2006{
2007 BL_ASSERT(this->domain.contains(subbox));
2008 BL_ASSERT(scomp >= 0 && scomp + ncomp <= this->nvar);
2009
2010 Array4<T const> const& a = this->const_array();
2011 Array4<int const> const& m = mask.const_array();
2012 Real r = 0.0;
2013#ifdef AMREX_USE_GPU
2014 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2015 ReduceOps<ReduceOpMax> reduce_op;
2016 ReduceData<Real> reduce_data(reduce_op);
2017 using ReduceTuple = ReduceData<Real>::Type;
2018 reduce_op.eval(subbox, reduce_data,
2019 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2020 {
2021 Real t = 0.0;
2022 if (m(i,j,k)) {
2023 for (int n = 0; n < ncomp; ++n) {
2024 t = amrex::max(t,static_cast<Real>(std::abs(a(i,j,k,n+scomp))));
2025 }
2026 }
2027 return {t};
2028 });
2029 ReduceTuple hv = reduce_data.value(reduce_op);
2030 r = amrex::get<0>(hv);
2031 } else
2032#endif
2033 {
2034 amrex::LoopOnCpu(subbox, ncomp, [=,&r] (int i, int j, int k, int n) noexcept
2035 {
2036 if (m(i,j,k)) {
2037 Real t = static_cast<Real>(std::abs(a(i,j,k,n+scomp)));
2038 r = amrex::max(r,t);
2039 }
2040 });
2041 }
2042 return r;
2043}
2044
2045template <class T>
2046template <RunOn run_on>
2047Real
2048BaseFab<T>::norm (int p, int comp, int numcomp) const
2049{
2050 return norm<run_on>(this->domain,p,comp,numcomp);
2051}
2052
2053template <class T>
2054template <RunOn run_on>
2055Real
2056BaseFab<T>::norm (const Box& subbox, int p, int comp, int numcomp) const
2057{
2058 BL_ASSERT(this->domain.contains(subbox));
2059 BL_ASSERT(comp >= 0 && comp + numcomp <= this->nvar);
2060
2061 Array4<T const> const& a = this->const_array();
2062 Real nrm = 0.;
2063#ifdef AMREX_USE_GPU
2064 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2065 if (p == 0) {
2066 ReduceOps<ReduceOpMax> reduce_op;
2067 ReduceData<Real> reduce_data(reduce_op);
2068 using ReduceTuple = ReduceData<Real>::Type;
2069 reduce_op.eval(subbox, reduce_data,
2070 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2071 {
2072 Real t = 0.0;
2073 for (int n = 0; n < numcomp; ++n) {
2074 t = amrex::max(t, static_cast<Real>(std::abs(a(i,j,k,n+comp))));
2075 }
2076 return {t};
2077 });
2078 ReduceTuple hv = reduce_data.value(reduce_op);
2079 nrm = amrex::get<0>(hv);
2080 } else if (p == 1) {
2081 ReduceOps<ReduceOpSum> reduce_op;
2082 ReduceData<Real> reduce_data(reduce_op);
2083 using ReduceTuple = ReduceData<Real>::Type;
2084 reduce_op.eval(subbox, reduce_data,
2085 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2086 {
2087 Real t = 0.0;
2088 for (int n = 0; n < numcomp; ++n) {
2089 t += static_cast<Real>(std::abs(a(i,j,k,n+comp)));
2090 }
2091 return {t};
2092 });
2093 ReduceTuple hv = reduce_data.value(reduce_op);
2094 nrm = amrex::get<0>(hv);
2095 } else if (p == 2) {
2096 ReduceOps<ReduceOpSum> reduce_op;
2097 ReduceData<Real> reduce_data(reduce_op);
2098 using ReduceTuple = ReduceData<Real>::Type;
2099 reduce_op.eval(subbox, reduce_data,
2100 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2101 {
2102 Real t = 0.0;
2103 for (int n = 0; n < numcomp; ++n) {
2104 t += static_cast<Real>(a(i,j,k,n+comp)*a(i,j,k,n+comp));
2105 }
2106 return {t};
2107 });
2108 ReduceTuple hv = reduce_data.value(reduce_op);
2109 nrm = amrex::get<0>(hv);
2110 } else {
2111 amrex::Error("BaseFab<T>::norm: wrong p");
2112 }
2113 } else
2114#endif
2115 {
2116 if (p == 0) {
2117 amrex::LoopOnCpu(subbox, numcomp, [=,&nrm] (int i, int j, int k, int n) noexcept
2118 {
2119 Real t = static_cast<Real>(std::abs(a(i,j,k,n+comp)));
2120 nrm = amrex::max(nrm,t);
2121 });
2122 } else if (p == 1) {
2123 amrex::LoopOnCpu(subbox, numcomp, [=,&nrm] (int i, int j, int k, int n) noexcept
2124 {
2125 nrm += std::abs(a(i,j,k,n+comp));
2126 });
2127 } else if (p == 2) {
2128 amrex::LoopOnCpu(subbox, numcomp, [=,&nrm] (int i, int j, int k, int n) noexcept
2129 {
2130 nrm += a(i,j,k,n+comp)*a(i,j,k,n+comp);
2131 });
2132 } else {
2133 amrex::Error("BaseFab<T>::norm: wrong p");
2134 }
2135 }
2136
2137 return nrm;
2138}
2139
2140template <class T>
2141template <RunOn run_on>
2142T
2143BaseFab<T>::min (int comp) const noexcept
2144{
2145 return this->min<run_on>(this->domain,comp);
2146}
2147
2148template <class T>
2149template <RunOn run_on>
2150T
2151BaseFab<T>::min (const Box& subbox, int comp) const noexcept
2152{
2153 Array4<T const> const& a = this->const_array(comp);
2154#ifdef AMREX_USE_GPU
2155 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2156 ReduceOps<ReduceOpMin> reduce_op;
2157 ReduceData<T> reduce_data(reduce_op);
2158 using ReduceTuple = typename decltype(reduce_data)::Type;
2159 reduce_op.eval(subbox, reduce_data,
2160 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2161 {
2162 return { a(i,j,k) };
2163 });
2164 ReduceTuple hv = reduce_data.value(reduce_op);
2165 return amrex::get<0>(hv);
2166 } else
2167#endif
2168 {
2169 T r = std::numeric_limits<T>::max();
2170 amrex::LoopOnCpu(subbox, [=,&r] (int i, int j, int k) noexcept
2171 {
2172 r = amrex::min(r, a(i,j,k));
2173 });
2174 return r;
2175 }
2176}
2177
2178template <class T>
2179template <RunOn run_on>
2180T
2181BaseFab<T>::max (int comp) const noexcept
2182{
2183 return this->max<run_on>(this->domain,comp);
2184}
2185
2186template <class T>
2187template <RunOn run_on>
2188T
2189BaseFab<T>::max (const Box& subbox, int comp) const noexcept
2190{
2191 Array4<T const> const& a = this->const_array(comp);
2192#ifdef AMREX_USE_GPU
2193 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2194 ReduceOps<ReduceOpMax> reduce_op;
2195 ReduceData<T> reduce_data(reduce_op);
2196 using ReduceTuple = typename decltype(reduce_data)::Type;
2197 reduce_op.eval(subbox, reduce_data,
2198 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2199 {
2200 return { a(i,j,k) };
2201 });
2202 ReduceTuple hv = reduce_data.value(reduce_op);
2203 return amrex::get<0>(hv);
2204 } else
2205#endif
2206 {
2207 T r = std::numeric_limits<T>::lowest();
2208 amrex::LoopOnCpu(subbox, [=,&r] (int i, int j, int k) noexcept
2209 {
2210 r = amrex::max(r, a(i,j,k));
2211 });
2212 return r;
2213 }
2214}
2215
2216template <class T>
2217template <RunOn run_on>
2218std::pair<T,T>
2219BaseFab<T>::minmax (int comp) const noexcept
2220{
2221 return this->minmax<run_on>(this->domain,comp);
2222}
2223
2224template <class T>
2225template <RunOn run_on>
2226std::pair<T,T>
2227BaseFab<T>::minmax (const Box& subbox, int comp) const noexcept
2228{
2229 Array4<T const> const& a = this->const_array(comp);
2230#ifdef AMREX_USE_GPU
2231 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2233 ReduceData<T,T> reduce_data(reduce_op);
2234 using ReduceTuple = typename decltype(reduce_data)::Type;
2235 reduce_op.eval(subbox, reduce_data,
2236 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2237 {
2238 auto const x = a(i,j,k);
2239 return { x, x };
2240 });
2241 ReduceTuple hv = reduce_data.value(reduce_op);
2242 return std::make_pair(amrex::get<0>(hv), amrex::get<1>(hv));
2243 } else
2244#endif
2245 {
2246 T rmax = std::numeric_limits<T>::lowest();
2247 T rmin = std::numeric_limits<T>::max();
2248 amrex::LoopOnCpu(subbox, [=,&rmin,&rmax] (int i, int j, int k) noexcept
2249 {
2250 auto const x = a(i,j,k);
2251 rmin = amrex::min(rmin, x);
2252 rmax = amrex::max(rmax, x);
2253 });
2254 return std::make_pair(rmin,rmax);
2255 }
2256}
2257
2258template <class T>
2259template <RunOn run_on>
2260T
2261BaseFab<T>::maxabs (int comp) const noexcept
2262{
2263 return this->maxabs<run_on>(this->domain,comp);
2264}
2265
2266template <class T>
2267template <RunOn run_on>
2268T
2269BaseFab<T>::maxabs (const Box& subbox, int comp) const noexcept
2270{
2271 Array4<T const> const& a = this->const_array(comp);
2272#ifdef AMREX_USE_GPU
2273 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2274 ReduceOps<ReduceOpMax> reduce_op;
2275 ReduceData<T> reduce_data(reduce_op);
2276 using ReduceTuple = typename decltype(reduce_data)::Type;
2277 reduce_op.eval(subbox, reduce_data,
2278 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2279 {
2280 return { std::abs(a(i,j,k)) };
2281 });
2282 ReduceTuple hv = reduce_data.value(reduce_op);
2283 return amrex::get<0>(hv);
2284 } else
2285#endif
2286 {
2287 T r = 0;
2288 amrex::LoopOnCpu(subbox, [=,&r] (int i, int j, int k) noexcept
2289 {
2290 r = amrex::max(r, std::abs(a(i,j,k)));
2291 });
2292 return r;
2293 }
2294}
2295
2296
2297template <class T>
2298template <RunOn run_on>
2299IntVect
2300BaseFab<T>::indexFromValue (Box const& subbox, int comp, T const& value) const noexcept
2301{
2302 Array4<T const> const& a = this->const_array(comp);
2303#ifdef AMREX_USE_GPU
2304 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2305 Array<int,1+AMREX_SPACEDIM> ha{0,AMREX_D_DECL(std::numeric_limits<int>::lowest(),
2306 std::numeric_limits<int>::lowest(),
2307 std::numeric_limits<int>::lowest())};
2308 Gpu::DeviceVector<int> dv(1+AMREX_SPACEDIM);
2309 int* p = dv.data();
2310 Gpu::htod_memcpy_async(p, ha.data(), sizeof(int)*ha.size());
2311 amrex::ParallelFor(subbox, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
2312 {
2313 int* flag = p;
2314 if ((*flag == 0) && (a(i,j,k) == value)) {
2315 if (Gpu::Atomic::Exch(flag,1) == 0) {
2316 AMREX_D_TERM(p[1] = i;,
2317 p[2] = j;,
2318 p[3] = k;);
2319 }
2320 }
2321 });
2322 Gpu::dtoh_memcpy_async(ha.data(), p, sizeof(int)*ha.size());
2324 return IntVect(AMREX_D_DECL(ha[1],ha[2],ha[3]));
2325 } else
2326#endif
2327 {
2328 AMREX_LOOP_3D(subbox, i, j, k,
2329 {
2330 if (a(i,j,k) == value) { return IntVect(AMREX_D_DECL(i,j,k)); }
2331 });
2332 return IntVect::TheMinVector();
2333 }
2334}
2335
2336template <class T>
2337template <RunOn run_on>
2338IntVect
2339BaseFab<T>::minIndex (int comp) const noexcept
2340{
2341 return this->minIndex<run_on>(this->domain,comp);
2342}
2343
2344template <class T>
2345template <RunOn run_on>
2346IntVect
2347BaseFab<T>::minIndex (const Box& subbox, int comp) const noexcept
2348{
2349 T min_val = this->min<run_on>(subbox, comp);
2350 return this->indexFromValue<run_on>(subbox, comp, min_val);
2351}
2352
2353template <class T>
2354template <RunOn run_on>
2355void
2356BaseFab<T>::minIndex (const Box& subbox, Real& min_val, IntVect& min_idx, int comp) const noexcept
2357{
2358 min_val = this->min<run_on>(subbox, comp);
2359 min_idx = this->indexFromValue<run_on>(subbox, comp, min_val);
2360}
2361
2362template <class T>
2363template <RunOn run_on>
2364IntVect
2365BaseFab<T>::maxIndex (int comp) const noexcept
2366{
2367 return this->maxIndex<run_on>(this->domain,comp);
2368}
2369
2370template <class T>
2371template <RunOn run_on>
2372IntVect
2373BaseFab<T>::maxIndex (const Box& subbox, int comp) const noexcept
2374{
2375 T max_val = this->max<run_on>(subbox, comp);
2376 return this->indexFromValue<run_on>(subbox, comp, max_val);
2377}
2378
2379template <class T>
2380template <RunOn run_on>
2381void
2382BaseFab<T>::maxIndex (const Box& subbox, Real& max_val, IntVect& max_idx, int comp) const noexcept
2383{
2384 max_val = this->max<run_on>(subbox, comp);
2385 max_idx = this->indexFromValue<run_on>(subbox, comp, max_val);
2386}
2387
2388template <class T>
2389template <RunOn run_on>
2390int
2391BaseFab<T>::maskLT (BaseFab<int>& mask, T const& val, int comp) const noexcept
2392{
2393 mask.resize(this->domain,1);
2394 int cnt = 0;
2395 Array4<int> const& m = mask.array();
2396 Array4<T const> const& a = this->const_array(comp);
2397#ifdef AMREX_USE_GPU
2398 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2399 ReduceOps<ReduceOpSum> reduce_op;
2400 ReduceData<int> reduce_data(reduce_op);
2401 using ReduceTuple = typename decltype(reduce_data)::Type;
2402 reduce_op.eval(this->domain, reduce_data,
2403 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2404 {
2405 int t;
2406 if (a(i,j,k) < val) {
2407 m(i,j,k) = 1;
2408 t = 1;
2409 } else {
2410 m(i,j,k) = 0;
2411 t = 0;
2412 }
2413 return {t};
2414 });
2415 ReduceTuple hv = reduce_data.value(reduce_op);
2416 cnt = amrex::get<0>(hv);
2417 } else
2418#endif
2419 {
2420 AMREX_LOOP_3D(this->domain, i, j, k,
2421 {
2422 if (a(i,j,k) < val) {
2423 m(i,j,k) = 1;
2424 ++cnt;
2425 } else {
2426 m(i,j,k) = 0;
2427 }
2428 });
2429 }
2430
2431 return cnt;
2432}
2433
2434template <class T>
2435template <RunOn run_on>
2436int
2437BaseFab<T>::maskLE (BaseFab<int>& mask, T const& val, int comp) const noexcept
2438{
2439 mask.resize(this->domain,1);
2440 int cnt = 0;
2441 Array4<int> const& m = mask.array();
2442 Array4<T const> const& a = this->const_array(comp);
2443#ifdef AMREX_USE_GPU
2444 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2445 ReduceOps<ReduceOpSum> reduce_op;
2446 ReduceData<int> reduce_data(reduce_op);
2447 using ReduceTuple = typename decltype(reduce_data)::Type;
2448 reduce_op.eval(this->domain, reduce_data,
2449 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2450 {
2451 int t;
2452 if (a(i,j,k) <= val) {
2453 m(i,j,k) = 1;
2454 t = 1;
2455 } else {
2456 m(i,j,k) = 0;
2457 t = 0;
2458 }
2459 return {t};
2460 });
2461 ReduceTuple hv = reduce_data.value(reduce_op);
2462 cnt = amrex::get<0>(hv);
2463 } else
2464#endif
2465 {
2466 AMREX_LOOP_3D(this->domain, i, j, k,
2467 {
2468 if (a(i,j,k) <= val) {
2469 m(i,j,k) = 1;
2470 ++cnt;
2471 } else {
2472 m(i,j,k) = 0;
2473 }
2474 });
2475 }
2476
2477 return cnt;
2478}
2479
2480template <class T>
2481template <RunOn run_on>
2482int
2483BaseFab<T>::maskEQ (BaseFab<int>& mask, T const& val, int comp) const noexcept
2484{
2485 mask.resize(this->domain,1);
2486 int cnt = 0;
2487 Array4<int> const& m = mask.array();
2488 Array4<T const> const& a = this->const_array(comp);
2489#ifdef AMREX_USE_GPU
2490 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2491 ReduceOps<ReduceOpSum> reduce_op;
2492 ReduceData<int> reduce_data(reduce_op);
2493 using ReduceTuple = typename decltype(reduce_data)::Type;
2494 reduce_op.eval(this->domain, reduce_data,
2495 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2496 {
2497 int t;
2498 if (a(i,j,k) == val) {
2499 m(i,j,k) = 1;
2500 t = 1;
2501 } else {
2502 m(i,j,k) = 0;
2503 t = 0;
2504 }
2505 return {t};
2506 });
2507 ReduceTuple hv = reduce_data.value(reduce_op);
2508 cnt = amrex::get<0>(hv);
2509 } else
2510#endif
2511 {
2512 AMREX_LOOP_3D(this->domain, i, j, k,
2513 {
2514 if (a(i,j,k) == val) {
2515 m(i,j,k) = 1;
2516 ++cnt;
2517 } else {
2518 m(i,j,k) = 0;
2519 }
2520 });
2521 }
2522
2523 return cnt;
2524}
2525
2526template <class T>
2527template <RunOn run_on>
2528int
2529BaseFab<T>::maskGT (BaseFab<int>& mask, T const& val, int comp) const noexcept
2530{
2531 mask.resize(this->domain,1);
2532 int cnt = 0;
2533 Array4<int> const& m = mask.array();
2534 Array4<T const> const& a = this->const_array(comp);
2535#ifdef AMREX_USE_GPU
2536 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2537 ReduceOps<ReduceOpSum> reduce_op;
2538 ReduceData<int> reduce_data(reduce_op);
2539 using ReduceTuple = typename decltype(reduce_data)::Type;
2540 reduce_op.eval(this->domain, reduce_data,
2541 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2542 {
2543 int t;
2544 if (a(i,j,k) > val) {
2545 m(i,j,k) = 1;
2546 t = 1;
2547 } else {
2548 m(i,j,k) = 0;
2549 t = 0;
2550 }
2551 return {t};
2552 });
2553 ReduceTuple hv = reduce_data.value(reduce_op);
2554 cnt = amrex::get<0>(hv);
2555 } else
2556#endif
2557 {
2558 AMREX_LOOP_3D(this->domain, i, j, k,
2559 {
2560 if (a(i,j,k) > val) {
2561 m(i,j,k) = 1;
2562 ++cnt;
2563 } else {
2564 m(i,j,k) = 0;
2565 }
2566 });
2567 }
2568
2569 return cnt;
2570}
2571
2572template <class T>
2573template <RunOn run_on>
2574int
2575BaseFab<T>::maskGE (BaseFab<int>& mask, T const& val, int comp) const noexcept
2576{
2577 mask.resize(this->domain,1);
2578 int cnt = 0;
2579 Array4<int> const& m = mask.array();
2580 Array4<T const> const& a = this->const_array(comp);
2581#ifdef AMREX_USE_GPU
2582 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2583 ReduceOps<ReduceOpSum> reduce_op;
2584 ReduceData<int> reduce_data(reduce_op);
2585 using ReduceTuple = typename decltype(reduce_data)::Type;
2586 reduce_op.eval(this->domain, reduce_data,
2587 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2588 {
2589 int t;
2590 if (a(i,j,k) >= val) {
2591 m(i,j,k) = 1;
2592 t = 1;
2593 } else {
2594 m(i,j,k) = 0;
2595 t = 0;
2596 }
2597 return {t};
2598 });
2599 ReduceTuple hv = reduce_data.value(reduce_op);
2600 cnt = amrex::get<0>(hv);
2601 } else
2602#endif
2603 {
2604 AMREX_LOOP_3D(this->domain, i, j, k,
2605 {
2606 if (a(i,j,k) >= val) {
2607 m(i,j,k) = 1;
2608 ++cnt;
2609 } else {
2610 m(i,j,k) = 0;
2611 }
2612 });
2613 }
2614
2615 return cnt;
2616}
2617
2618template <class T>
2619template <RunOn run_on>
2622{
2623 Box ovlp(this->domain);
2624 ovlp &= x.domain;
2625 return ovlp.ok() ? this->atomicAdd<run_on>(x,ovlp,ovlp,0,0,this->nvar) : *this;
2626}
2627
2628template <class T>
2629template <RunOn run_on>
2631BaseFab<T>::saxpy (T a, const BaseFab<T>& x, const Box& srcbox, const Box& destbox,
2632 int srccomp, int destcomp, int numcomp) noexcept
2633{
2634 BL_ASSERT(srcbox.ok());
2635 BL_ASSERT(x.box().contains(srcbox));
2636 BL_ASSERT(destbox.ok());
2637 BL_ASSERT(box().contains(destbox));
2638 BL_ASSERT(destbox.sameSize(srcbox));
2639 BL_ASSERT( srccomp >= 0 && srccomp+numcomp <= x.nComp());
2640 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2641
2642 Array4<T> const& d = this->array();
2643 Array4<T const> const& s = x.const_array();
2644 const auto dlo = amrex::lbound(destbox);
2645 const auto slo = amrex::lbound(srcbox);
2646 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
2647 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
2648 {
2649 d(i,j,k,n+destcomp) += a * s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
2650 });
2651
2652 return *this;
2653}
2654
2655template <class T>
2656template <RunOn run_on>
2658BaseFab<T>::saxpy (T a, const BaseFab<T>& x) noexcept
2659{
2660 Box ovlp(this->domain);
2661 ovlp &= x.domain;
2662 return ovlp.ok() ? saxpy<run_on>(a,x,ovlp,ovlp,0,0,this->nvar) : *this;
2663}
2664
2665template <class T>
2666template <RunOn run_on>
2669 const Box& srcbox, const Box& destbox,
2670 int srccomp, int destcomp, int numcomp) noexcept
2671{
2672 BL_ASSERT(srcbox.ok());
2673 BL_ASSERT(x.box().contains(srcbox));
2674 BL_ASSERT(destbox.ok());
2675 BL_ASSERT(box().contains(destbox));
2676 BL_ASSERT(destbox.sameSize(srcbox));
2677 BL_ASSERT( srccomp >= 0 && srccomp+numcomp <= x.nComp());
2678 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2679
2680 Array4<T> const& d = this->array();
2681 Array4<T const> const& s = x.const_array();
2682 const auto dlo = amrex::lbound(destbox);
2683 const auto slo = amrex::lbound(srcbox);
2684 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
2685 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
2686 {
2687 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);
2688 });
2689
2690 return *this;
2691}
2692
2693template <class T>
2694template <RunOn run_on>
2696BaseFab<T>::addproduct (const Box& destbox, int destcomp, int numcomp,
2697 const BaseFab<T>& src1, int comp1,
2698 const BaseFab<T>& src2, int comp2) noexcept
2699{
2700 BL_ASSERT(destbox.ok());
2701 BL_ASSERT(box().contains(destbox));
2702 BL_ASSERT( comp1 >= 0 && comp1+numcomp <= src1.nComp());
2703 BL_ASSERT( comp2 >= 0 && comp2+numcomp <= src2.nComp());
2704 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2705
2706 Array4<T> const& d = this->array();
2707 Array4<T const> const& s1 = src1.const_array();
2708 Array4<T const> const& s2 = src2.const_array();
2709 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
2710 {
2711 d(i,j,k,n+destcomp) += s1(i,j,k,n+comp1) * s2(i,j,k,n+comp2);
2712 });
2713
2714 return *this;
2715}
2716
2717template <class T>
2718template <RunOn run_on>
2720BaseFab<T>::linComb (const BaseFab<T>& f1, const Box& b1, int comp1,
2721 const BaseFab<T>& f2, const Box& b2, int comp2,
2722 Real alpha, Real beta, const Box& b,
2723 int comp, int numcomp) noexcept
2724{
2725 BL_ASSERT(b1.ok());
2726 BL_ASSERT(f1.box().contains(b1));
2727 BL_ASSERT(b2.ok());
2728 BL_ASSERT(f2.box().contains(b2));
2729 BL_ASSERT(b.ok());
2730 BL_ASSERT(box().contains(b));
2731 BL_ASSERT(b.sameSize(b1));
2732 BL_ASSERT(b.sameSize(b2));
2733 BL_ASSERT(comp1 >= 0 && comp1+numcomp <= f1.nComp());
2734 BL_ASSERT(comp2 >= 0 && comp2+numcomp <= f2.nComp());
2735 BL_ASSERT(comp >= 0 && comp +numcomp <= nComp());
2736
2737 Array4<T> const& d = this->array();
2738 Array4<T const> const& s1 = f1.const_array();
2739 Array4<T const> const& s2 = f2.const_array();
2740 const auto dlo = amrex::lbound(b);
2741 const auto slo1 = amrex::lbound(b1);
2742 const auto slo2 = amrex::lbound(b2);
2743 const Dim3 off1{.x = slo1.x-dlo.x, .y = slo1.y-dlo.y, .z = slo1.z-dlo.z};
2744 const Dim3 off2{.x = slo2.x-dlo.x, .y = slo2.y-dlo.y, .z = slo2.z-dlo.z};
2745
2746 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, b, numcomp, i, j, k, n,
2747 {
2748 d(i,j,k,n+comp) = alpha*s1(i+off1.x,j+off1.y,k+off1.z,n+comp1)
2749 + beta*s2(i+off2.x,j+off2.y,k+off2.z,n+comp2);
2750 });
2751 return *this;
2752}
2753
2754template <class T>
2755template <RunOn run_on>
2756T
2757BaseFab<T>::dot (const Box& xbx, int xcomp,
2758 const BaseFab<T>& y, const Box& ybx, int ycomp,
2759 int numcomp) const noexcept
2760{
2761 BL_ASSERT(xbx.ok());
2762 BL_ASSERT(box().contains(xbx));
2763 BL_ASSERT(y.box().contains(ybx));
2764 BL_ASSERT(xbx.sameSize(ybx));
2765 BL_ASSERT(xcomp >= 0 && xcomp+numcomp <= nComp());
2766 BL_ASSERT(ycomp >= 0 && ycomp+numcomp <= y.nComp());
2767
2768 T r = 0;
2769
2770 const auto xlo = amrex::lbound(xbx);
2771 const auto ylo = amrex::lbound(ybx);
2772 const Dim3 offset{.x = ylo.x-xlo.x, .y = ylo.y-xlo.y, .z = ylo.z-xlo.z};
2773 Array4<T const> const& xa = this->const_array();
2774 Array4<T const> const& ya = y.const_array();
2775
2776#ifdef AMREX_USE_GPU
2777 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2778 ReduceOps<ReduceOpSum> reduce_op;
2779 ReduceData<T> reduce_data(reduce_op);
2780 using ReduceTuple = typename decltype(reduce_data)::Type;
2781 reduce_op.eval(xbx, reduce_data,
2782 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2783 {
2784 T t = 0;
2785 for (int n = 0; n < numcomp; ++n) {
2786 t += xa(i,j,k,n+xcomp) * ya(i+offset.x,j+offset.y,k+offset.z,n+ycomp);
2787 }
2788 return {t};
2789 });
2790 ReduceTuple hv = reduce_data.value(reduce_op);
2791 r = amrex::get<0>(hv);
2792 } else
2793#endif
2794 {
2795 AMREX_LOOP_4D(xbx, numcomp, i, j, k, n,
2796 {
2797 r += xa(i,j,k,n+xcomp) * ya(i+offset.x,j+offset.y,k+offset.z,n+ycomp);
2798 });
2799 }
2800
2801 return r;
2802}
2803
2804template <class T>
2805template <RunOn run_on>
2806T
2807BaseFab<T>::dotmask (const BaseFab<int>& mask, const Box& xbx, int xcomp,
2808 const BaseFab<T>& y, const Box& ybx, int ycomp,
2809 int numcomp) const noexcept
2810{
2811 BL_ASSERT(xbx.ok());
2812 BL_ASSERT(box().contains(xbx));
2813 BL_ASSERT(y.box().contains(ybx));
2814 BL_ASSERT(xbx.sameSize(ybx));
2815 BL_ASSERT(xcomp >= 0 && xcomp+numcomp <= nComp());
2816 BL_ASSERT(ycomp >= 0 && ycomp+numcomp <= y.nComp());
2817
2818 T r = 0;
2819
2820 const auto xlo = amrex::lbound(xbx);
2821 const auto ylo = amrex::lbound(ybx);
2822 const Dim3 offset{.x = ylo.x-xlo.x, .y = ylo.y-xlo.y, .z = ylo.z-xlo.z};
2823
2824 Array4<T const> const& xa = this->const_array();
2825 Array4<T const> const& ya = y.const_array();
2826 Array4<int const> const& ma = mask.const_array();
2827
2828#ifdef AMREX_USE_GPU
2829 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2830 ReduceOps<ReduceOpSum> reduce_op;
2831 ReduceData<T> reduce_data(reduce_op);
2832 using ReduceTuple = typename decltype(reduce_data)::Type;
2833 reduce_op.eval(xbx, reduce_data,
2834 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2835 {
2836 int m = static_cast<int>(static_cast<bool>(ma(i,j,k)));
2837 T t = 0;
2838 for (int n = 0; n < numcomp; ++n) {
2839 t += xa(i,j,k,n+xcomp) * ya(i+offset.x,j+offset.y,k+offset.z,n+ycomp) * m;
2840 }
2841 return {t};
2842 });
2843 ReduceTuple hv = reduce_data.value(reduce_op);
2844 r = amrex::get<0>(hv);
2845 } else
2846#endif
2847 {
2848 AMREX_LOOP_4D(xbx, numcomp, i, j, k, n,
2849 {
2850 int m = static_cast<int>(static_cast<bool>(ma(i,j,k)));
2851 r += xa(i,j,k,n+xcomp) * ya(i+offset.x,j+offset.y,k+offset.z,n+ycomp) * m;
2852 });
2853 }
2854
2855 return r;
2856}
2857
2858template <class T>
2859template <RunOn run_on>
2860T
2861BaseFab<T>::sum (int comp, int numcomp) const noexcept
2862{
2863 return this->sum<run_on>(this->domain, DestComp{comp}, NumComps{numcomp});
2864}
2865
2866template <class T>
2867template <RunOn run_on>
2868T
2869BaseFab<T>::sum (const Box& subbox, int comp, int numcomp) const noexcept
2870{
2871 return this->sum<run_on>(subbox, DestComp{comp}, NumComps{numcomp});
2872}
2873
2874template <class T>
2875template <RunOn run_on>
2877BaseFab<T>::negate (int comp, int numcomp) noexcept
2878{
2879 return this->negate<run_on>(this->domain, DestComp{comp}, NumComps{numcomp});
2880}
2881
2882template <class T>
2883template <RunOn run_on>
2885BaseFab<T>::negate (const Box& b, int comp, int numcomp) noexcept
2886{
2887 return this->negate<run_on>(b, DestComp{comp}, NumComps{numcomp});
2888}
2889
2890template <class T>
2891template <RunOn run_on>
2893BaseFab<T>::invert (T const& r, int comp, int numcomp) noexcept
2894{
2895 return this->invert<run_on>(r, this->domain, DestComp{comp}, NumComps{numcomp});
2896}
2897
2898template <class T>
2899template <RunOn run_on>
2901BaseFab<T>::invert (T const& r, const Box& b, int comp, int numcomp) noexcept
2902{
2903 return this->invert<run_on>(r, b, DestComp{comp}, NumComps{numcomp});
2904}
2905
2906template <class T>
2907template <RunOn run_on>
2909BaseFab<T>::plus (T const& r, int comp, int numcomp) noexcept
2910{
2911 return this->plus<run_on>(r, this->domain, DestComp{comp}, NumComps{numcomp});
2912}
2913
2914template <class T>
2915template <RunOn run_on>
2917BaseFab<T>::plus (T const& r, const Box& b, int comp, int numcomp) noexcept
2918{
2919 return this->plus<run_on>(r, b, DestComp{comp}, NumComps{numcomp});
2920}
2921
2922template <class T>
2923template <RunOn run_on>
2925BaseFab<T>::plus (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
2926{
2927 return this->plus<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
2928}
2929
2930template <class T>
2931template <RunOn run_on>
2933BaseFab<T>::atomicAdd (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
2934{
2935 Box ovlp(this->domain);
2936 ovlp &= src.domain;
2937 return ovlp.ok() ? this->atomicAdd<run_on>(src,ovlp,ovlp,srccomp,destcomp,numcomp) : *this;
2938}
2939
2940template <class T>
2941template <RunOn run_on>
2943BaseFab<T>::plus (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
2944 int numcomp) noexcept
2945{
2946 return this->plus<run_on>(src, subbox, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
2947}
2948
2949template <class T>
2950template <RunOn run_on>
2952BaseFab<T>::atomicAdd (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
2953 int numcomp) noexcept
2954{
2955 Box ovlp(this->domain);
2956 ovlp &= src.domain;
2957 ovlp &= subbox;
2958 return ovlp.ok() ? this->atomicAdd<run_on>(src,ovlp,ovlp,srccomp,destcomp,numcomp) : *this;
2959}
2960
2961template <class T>
2962template <RunOn run_on>
2964BaseFab<T>::plus (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
2965 int srccomp, int destcomp, int numcomp) noexcept
2966{
2967 BL_ASSERT(destbox.ok());
2968 BL_ASSERT(src.box().contains(srcbox));
2969 BL_ASSERT(box().contains(destbox));
2970 BL_ASSERT(destbox.sameSize(srcbox));
2971 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
2972 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2973
2974 Array4<T> const& d = this->array();
2975 Array4<T const> const& s = src.const_array();
2976 const auto dlo = amrex::lbound(destbox);
2977 const auto slo = amrex::lbound(srcbox);
2978 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
2979 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
2980 {
2981 d(i,j,k,n+destcomp) += s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
2982 });
2983
2984 return *this;
2985}
2986
2987namespace detail
2988{
2989
2990template <RunOn run_on, typename T>
2991requires (HasAtomicAdd<T>::value)
2992void basefab_atomic_add (BaseFab<T>& dfab, const BaseFab<T>& sfab,
2993 const Box& srcbox, const Box& destbox,
2994 int srccomp, int destcomp, int numcomp) noexcept
2995{
2996 Array4<T> const& d = dfab.array();
2997 Array4<T const> const& s = sfab.const_array();
2998 const auto dlo = amrex::lbound(destbox);
2999 const auto slo = amrex::lbound(srcbox);
3000 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3001 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3002 {
3003 T* p = d.ptr(i,j,k,n+destcomp);
3004 HostDevice::Atomic::Add(p, s(i+offset.x,j+offset.y,k+offset.z,n+srccomp));
3005 });
3006}
3007
3008template <RunOn run_on, typename T>
3009requires (!HasAtomicAdd<T>::value)
3010void basefab_atomic_add (BaseFab<T>& dfab, const BaseFab<T>& sfab,
3011 const Box& srcbox, const Box& destbox,
3012 int srccomp, int destcomp, int numcomp)
3013{
3014 amrex::ignore_unused(dfab, sfab, srcbox, destbox, srccomp, destcomp, numcomp);
3015 amrex::Abort("BaseFab: atomicAdd not supported");
3016}
3017
3018}
3019
3020template <class T>
3021template <RunOn run_on>
3022BaseFab<T>&
3023BaseFab<T>::atomicAdd (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3024 int srccomp, int destcomp, int numcomp) noexcept
3025{
3026 BL_ASSERT(destbox.ok());
3027 BL_ASSERT(src.box().contains(srcbox));
3028 BL_ASSERT(box().contains(destbox));
3029 BL_ASSERT(destbox.sameSize(srcbox));
3030 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3031 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3032
3033 detail::basefab_atomic_add<run_on>(*this, src, srcbox, destbox,
3034 srccomp, destcomp, numcomp);
3035
3036 return *this;
3037}
3038
3039template <class T>
3040template <RunOn run_on>
3042BaseFab<T>::lockAdd (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3043 int srccomp, int destcomp, int numcomp) noexcept
3044{
3045#if defined(AMREX_USE_OMP) && (AMREX_SPACEDIM > 1)
3046#if defined(AMREX_USE_GPU)
3047 if (run_on == RunOn::Host || Gpu::notInLaunchRegion()) {
3048#endif
3049 BL_ASSERT(destbox.ok());
3050 BL_ASSERT(src.box().contains(srcbox));
3051 BL_ASSERT(box().contains(destbox));
3052 BL_ASSERT(destbox.sameSize(srcbox));
3053 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3054 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3055
3056 Array4<T> const& d = this->array();
3057 Array4<T const> const& s = src.const_array();
3058 auto const& dlo = amrex::lbound(destbox);
3059 auto const& dhi = amrex::ubound(destbox);
3060 auto const& len = amrex::length(destbox);
3061 auto const& slo = amrex::lbound(srcbox);
3062 Dim3 const offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3063
3064 int planedim;
3065 int nplanes;
3066 int plo;
3067 if (len.z == 1) {
3068 planedim = 1;
3069 nplanes = len.y;
3070 plo = dlo.y;
3071 } else {
3072 planedim = 2;
3073 nplanes = len.z;
3074 plo = dlo.z;
3075 }
3076
3077 auto* mask = (bool*) amrex_mempool_alloc(sizeof(bool)*nplanes);
3078 for (int ip = 0; ip < nplanes; ++ip) {
3079 mask[ip] = false;
3080 }
3081
3082 int mm = 0;
3083 int planes_left = nplanes;
3084 while (planes_left > 0) {
3085 AMREX_ASSERT(mm < nplanes);
3086 auto const m = mm + plo;
3087 auto* lock = OpenMP::get_lock(m);
3088 if (omp_test_lock(lock))
3089 {
3090 auto lo = dlo;
3091 auto hi = dhi;
3092 if (planedim == 1) {
3093 lo.y = m;
3094 hi.y = m;
3095 } else {
3096 lo.z = m;
3097 hi.z = m;
3098 }
3099
3100 for (int n = 0; n < numcomp; ++n) {
3101 for (int k = lo.z; k <= hi.z; ++k) {
3102 for (int j = lo.y; j <= hi.y; ++j) {
3103 auto * pdst = d.ptr(dlo.x,j ,k ,n+destcomp);
3104 auto const* psrc = s.ptr(slo.x,j+offset.y,k+offset.z,n+ srccomp);
3105#pragma omp simd
3106 for (int ii = 0; ii < len.x; ++ii) {
3107 pdst[ii] += psrc[ii];
3108 }
3109 }
3110 }
3111 }
3112
3113 mask[mm] = true;
3114 --planes_left;
3115 omp_unset_lock(lock);
3116 if (planes_left == 0) { break; }
3117 }
3118
3119 ++mm;
3120 for (int ip = 0; ip < nplanes; ++ip) {
3121 int new_mm = (mm+ip) % nplanes;
3122 if ( ! mask[new_mm] ) {
3123 mm = new_mm;
3124 break;
3125 }
3126 }
3127 }
3128
3130
3131 return *this;
3132
3133#if defined(AMREX_USE_GPU)
3134 } else {
3135 return this->template atomicAdd<run_on>(src, srcbox, destbox, srccomp, destcomp, numcomp);
3136 }
3137#endif
3138#else
3139 return this->template atomicAdd<run_on>(src, srcbox, destbox, srccomp, destcomp, numcomp);
3140#endif
3141}
3142
3143template <class T>
3144template <RunOn run_on>
3146BaseFab<T>::minus (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
3147{
3148 return this->minus<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3149}
3150
3151template <class T>
3152template <RunOn run_on>
3154BaseFab<T>::minus (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp, int numcomp) noexcept
3155{
3156 return this->minus<run_on>(src, subbox, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3157}
3158
3159template <class T>
3160template <RunOn run_on>
3162BaseFab<T>::minus (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3163 int srccomp, int destcomp, int numcomp) noexcept
3164{
3165 BL_ASSERT(destbox.ok());
3166 BL_ASSERT(src.box().contains(srcbox));
3167 BL_ASSERT(box().contains(destbox));
3168 BL_ASSERT(destbox.sameSize(srcbox));
3169 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3170 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3171
3172 Array4<T> const& d = this->array();
3173 Array4<T const> const& s = src.const_array();
3174 const auto dlo = amrex::lbound(destbox);
3175 const auto slo = amrex::lbound(srcbox);
3176 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3177 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3178 {
3179 d(i,j,k,n+destcomp) -= s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
3180 });
3181
3182 return *this;
3183}
3184
3185template <class T>
3186template <RunOn run_on>
3188BaseFab<T>::mult (T const& r, int comp, int numcomp) noexcept
3189{
3190 return this->mult<run_on>(r, this->domain, DestComp{comp}, NumComps{numcomp});
3191}
3192
3193template <class T>
3194template <RunOn run_on>
3196BaseFab<T>::mult (T const& r, const Box& b, int comp, int numcomp) noexcept
3197{
3198 return this->mult<run_on>(r, b, DestComp{comp}, NumComps{numcomp});
3199}
3200
3201template <class T>
3202template <RunOn run_on>
3204BaseFab<T>::mult (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
3205{
3206 return this->mult<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3207}
3208
3209template <class T>
3210template <RunOn run_on>
3212BaseFab<T>::mult (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp, int numcomp) noexcept
3213{
3214 return this->mult<run_on>(src, subbox, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3215}
3216
3217template <class T>
3218template <RunOn run_on>
3220BaseFab<T>::mult (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3221 int srccomp, int destcomp, int numcomp) noexcept
3222{
3223 BL_ASSERT(destbox.ok());
3224 BL_ASSERT(src.box().contains(srcbox));
3225 BL_ASSERT(box().contains(destbox));
3226 BL_ASSERT(destbox.sameSize(srcbox));
3227 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3228 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3229
3230 Array4<T> const& d = this->array();
3231 Array4<T const> const& s = src.const_array();
3232 const auto dlo = amrex::lbound(destbox);
3233 const auto slo = amrex::lbound(srcbox);
3234 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3235 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3236 {
3237 d(i,j,k,n+destcomp) *= s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
3238 });
3239
3240 return *this;
3241}
3242
3243template <class T>
3244template <RunOn run_on>
3246BaseFab<T>::divide (T const& r, int comp, int numcomp) noexcept
3247{
3248 return this->divide<run_on>(r, this->domain, DestComp{comp}, NumComps{numcomp});
3249}
3250
3251template <class T>
3252template <RunOn run_on>
3254BaseFab<T>::divide (T const& r, const Box& b, int comp, int numcomp) noexcept
3255{
3256 return this->divide<run_on>(r, b, DestComp{comp}, NumComps{numcomp});
3257}
3258
3259template <class T>
3260template <RunOn run_on>
3262BaseFab<T>::divide (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
3263{
3264 return this->divide<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3265}
3266
3267template <class T>
3268template <RunOn run_on>
3270BaseFab<T>::divide (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp, int numcomp) noexcept
3271{
3272 return this->divide<run_on>(src, subbox, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3273}
3274
3275template <class T>
3276template <RunOn run_on>
3278BaseFab<T>::divide (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3279 int srccomp, int destcomp, int numcomp) noexcept
3280{
3281 BL_ASSERT(destbox.ok());
3282 BL_ASSERT(src.box().contains(srcbox));
3283 BL_ASSERT(box().contains(destbox));
3284 BL_ASSERT(destbox.sameSize(srcbox));
3285 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3286 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3287
3288 Array4<T> const& d = this->array();
3289 Array4<T const> const& s = src.const_array();
3290 const auto dlo = amrex::lbound(destbox);
3291 const auto slo = amrex::lbound(srcbox);
3292 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3293 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3294 {
3295 d(i,j,k,n+destcomp) /= s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
3296 });
3297
3298 return *this;
3299}
3300
3301template <class T>
3302template <RunOn run_on>
3305{
3306 Box ovlp(this->domain);
3307 ovlp &= src.domain;
3308 return ovlp.ok() ? this->protected_divide<run_on>(src,ovlp,ovlp,0,0,this->nvar) : *this;
3309}
3310
3311template <class T>
3312template <RunOn run_on>
3314BaseFab<T>::protected_divide (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
3315{
3316 Box ovlp(this->domain);
3317 ovlp &= src.domain;
3318 return ovlp.ok() ? this->protected_divide<run_on>(src,ovlp,ovlp,srccomp,destcomp,numcomp) : *this;
3319}
3320
3321template <class T>
3322template <RunOn run_on>
3324BaseFab<T>::protected_divide (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
3325 int numcomp) noexcept
3326{
3327 Box ovlp(this->domain);
3328 ovlp &= src.domain;
3329 ovlp &= subbox;
3330 return ovlp.ok() ? this->protected_divide<run_on>(src,ovlp,ovlp,srccomp,destcomp,numcomp) : *this;
3331}
3332
3333template <class T>
3334template <RunOn run_on>
3336BaseFab<T>::protected_divide (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3337 int srccomp, int destcomp, int numcomp) noexcept
3338{
3339 BL_ASSERT(destbox.ok());
3340 BL_ASSERT(src.box().contains(srcbox));
3341 BL_ASSERT(box().contains(destbox));
3342 BL_ASSERT(destbox.sameSize(srcbox));
3343 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3344 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3345
3346 Array4<T> const& d = this->array();
3347 Array4<T const> const& s = src.const_array();
3348 const auto dlo = amrex::lbound(destbox);
3349 const auto slo = amrex::lbound(srcbox);
3350 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3351 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3352 {
3353 if (s(i+offset.x,j+offset.y,k+offset.z,n+srccomp) != 0) {
3354 d(i,j,k,n+destcomp) /= s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
3355 }
3356 });
3357
3358 return *this;
3359}
3360
3371template <class T>
3372template <RunOn run_on>
3374BaseFab<T>::linInterp (const BaseFab<T>& f1, const Box& b1, int comp1,
3375 const BaseFab<T>& f2, const Box& b2, int comp2,
3376 Real t1, Real t2, Real t,
3377 const Box& b, int comp, int numcomp) noexcept
3378{
3379 if (amrex::almostEqual(t1,t2) || amrex::almostEqual(t1,t)) {
3380 return copy<run_on>(f1,b1,comp1,b,comp,numcomp);
3381 } else if (amrex::almostEqual(t2,t)) {
3382 return copy<run_on>(f2,b2,comp2,b,comp,numcomp);
3383 } else {
3384 Real alpha = (t2-t)/(t2-t1);
3385 Real beta = (t-t1)/(t2-t1);
3386 return linComb<run_on>(f1,b1,comp1,f2,b2,comp2,alpha,beta,b,comp,numcomp);
3387 }
3388}
3389
3390template <class T>
3391template <RunOn run_on>
3393BaseFab<T>::linInterp (const BaseFab<T>& f1, int comp1,
3394 const BaseFab<T>& f2, int comp2,
3395 Real t1, Real t2, Real t,
3396 const Box& b, int comp, int numcomp) noexcept
3397{
3398 if (amrex::almostEqual(t1,t2) || amrex::almostEqual(t1,t)) {
3399 return copy<run_on>(f1,b,comp1,b,comp,numcomp);
3400 } else if (amrex::almostEqual(t2,t)) {
3401 return copy<run_on>(f2,b,comp2,b,comp,numcomp);
3402 } else {
3403 Real alpha = (t2-t)/(t2-t1);
3404 Real beta = (t-t1)/(t2-t1);
3405 return linComb<run_on>(f1,b,comp1,f2,b,comp2,alpha,beta,b,comp,numcomp);
3406 }
3407}
3408
3409//
3410// New interfaces
3411//
3412
3413template <class T>
3414template <RunOn run_on>
3415void
3416BaseFab<T>::setVal (T const& val) noexcept
3417{
3418 this->setVal<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3419}
3420
3421template <class T>
3422template <RunOn run_on>
3423void
3424BaseFab<T>::setVal (T const& x, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept
3425{
3426 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3427 Array4<T> const& a = this->array();
3428 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG (run_on, bx, ncomp.n, i, j, k, n,
3429 {
3430 a(i,j,k,n+dcomp.i) = x;
3431 });
3432}
3433
3434template <class T>
3435template <RunOn run_on>
3436void
3437BaseFab<T>::setValIf (T const& val, const BaseFab<int>& mask) noexcept
3438{
3439 this->setValIf<run_on>(val, this->domain, mask, DestComp{0}, NumComps{this->nvar});
3440}
3441
3442template <class T>
3443template <RunOn run_on>
3444void
3445BaseFab<T>::setValIf (T const& val, Box const& bx, const BaseFab<int>& mask, DestComp dcomp, NumComps ncomp) noexcept
3446{
3447 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3448 Array4<T> const& a = this->array();
3449 Array4<int const> const& m = mask.const_array();
3450 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG (run_on, bx, ncomp.n, i, j, k, n,
3451 {
3452 if (m(i,j,k)) { a(i,j,k,n+dcomp.i) = val; }
3453 });
3454}
3455
3456template <class T>
3457template <RunOn run_on>
3458void
3459BaseFab<T>::setValIfNot (T const& val, const BaseFab<int>& mask) noexcept
3460{
3461 this->setValIfNot<run_on>(val, this->domain, mask, DestComp{0}, NumComps{this->nvar});
3462}
3463
3464template <class T>
3465template <RunOn run_on>
3466void
3467BaseFab<T>::setValIfNot (T const& val, Box const& bx, const BaseFab<int>& mask, DestComp dcomp, NumComps ncomp) noexcept
3468{
3469 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3470 Array4<T> const& a = this->array();
3471 Array4<int const> const& m = mask.const_array();
3472 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG (run_on, bx, ncomp.n, i, j, k, n,
3473 {
3474 if (!m(i,j,k)) { a(i,j,k,n+dcomp.i) = val; }
3475 });
3476}
3477
3478template <class T>
3479template <RunOn run_on>
3480void
3481BaseFab<T>::setComplement (T const& x, const Box& bx, DestComp dcomp, NumComps ncomp) noexcept
3482{
3483#ifdef AMREX_USE_GPU
3484 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
3485 Array4<T> const& a = this->array();
3486 amrex::ParallelFor(this->domain, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
3487 {
3488 if (! bx.contains(IntVect(AMREX_D_DECL(i,j,k)))) {
3489 for (int n = dcomp.i; n < ncomp.n+dcomp.i; ++n) {
3490 a(i,j,k,n) = x;
3491 }
3492 }
3493 });
3494 } else
3495#endif
3496 {
3497 const BoxList b_lst = amrex::boxDiff(this->domain,bx);
3498 for (auto const& b : b_lst) {
3499 this->setVal<RunOn::Host>(x, b, dcomp, ncomp);
3500 }
3501 }
3502}
3503
3504template <class T>
3505template <RunOn run_on>
3507BaseFab<T>::copy (const BaseFab<T>& src) noexcept
3508{
3509 this->copy<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3510 return *this;
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>::plus (T const& val) noexcept
3539{
3540 return this->plus<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->plus<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3549}
3550
3551template <class T>
3552template <RunOn run_on>
3554BaseFab<T>::plus (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>::plus (const BaseFab<T>& src) noexcept
3571{
3572 return this->plus<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->plus<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>::minus (T const& val) noexcept
3609{
3610 return this->minus<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->minus<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3619}
3620
3621template <class T>
3622template <RunOn run_on>
3624BaseFab<T>::minus (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>::minus (const BaseFab<T>& src) noexcept
3641{
3642 return this->minus<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->minus<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>
3678BaseFab<T>::mult (T const& val) noexcept
3679{
3680 return this->mult<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3681}
3682
3683template <class T>
3684template <RunOn run_on>
3686BaseFab<T>::operator*= (T const& val) noexcept
3687{
3688 return this->mult<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3689}
3690
3691template <class T>
3692template <RunOn run_on>
3694BaseFab<T>::mult (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept
3695{
3696 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3697
3698 Array4<T> const& a = this->array();
3699 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3700 {
3701 a(i,j,k,n+dcomp.i) *= val;
3702 });
3703
3704 return *this;
3705}
3706
3707template <class T>
3708template <RunOn run_on>
3710BaseFab<T>::mult (const BaseFab<T>& src) noexcept
3711{
3712 return this->mult<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3713}
3714
3715template <class T>
3716template <RunOn run_on>
3719{
3720 return this->mult<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3721}
3722
3723template <class T>
3724template <RunOn run_on>
3727 SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
3728{
3729 AMREX_ASSERT(this->domain.sameType(src.domain));
3730 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3731 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3732
3733 bx &= src.domain;
3734
3735 Array4<T> const& d = this->array();
3736 Array4<T const> const& s = src.const_array();
3737 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3738 {
3739 d(i,j,k,n+dcomp.i) *= s(i,j,k,n+scomp.i);
3740 });
3741
3742 return *this;
3743}
3744
3745template <class T>
3746template <RunOn run_on>
3748BaseFab<T>::divide (T const& val) noexcept
3749{
3750 return this->divide<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3751}
3752
3753template <class T>
3754template <RunOn run_on>
3756BaseFab<T>::operator/= (T const& val) noexcept
3757{
3758 return this->divide<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3759}
3760
3761template <class T>
3762template <RunOn run_on>
3764BaseFab<T>::divide (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept
3765{
3766 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3767
3768 Array4<T> const& a = this->array();
3769 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3770 {
3771 a(i,j,k,n+dcomp.i) /= val;
3772 });
3773
3774 return *this;
3775}
3776
3777template <class T>
3778template <RunOn run_on>
3780BaseFab<T>::divide (const BaseFab<T>& src) noexcept
3781{
3782 return this->divide<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3783}
3784
3785template <class T>
3786template <RunOn run_on>
3789{
3790 return this->divide<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3791}
3792
3793template <class T>
3794template <RunOn run_on>
3797 SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
3798{
3799 AMREX_ASSERT(this->domain.sameType(src.domain));
3800 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3801 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3802
3803 bx &= src.domain;
3804
3805 Array4<T> const& d = this->array();
3806 Array4<T const> const& s = src.const_array();
3807 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3808 {
3809 d(i,j,k,n+dcomp.i) /= s(i,j,k,n+scomp.i);
3810 });
3811
3812 return *this;
3813}
3814
3815template <class T>
3816template <RunOn run_on>
3819{
3820 return this->negate<run_on>(this->domain, DestComp{0}, NumComps{this->nvar});
3821}
3822
3823template <class T>
3824template <RunOn run_on>
3826BaseFab<T>::negate (const Box& bx, DestComp dcomp, NumComps ncomp) noexcept
3827{
3828 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3829
3830 Array4<T> const& a = this->array();
3831 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3832 {
3833 a(i,j,k,n+dcomp.i) = -a(i,j,k,n+dcomp.i);
3834 });
3835
3836 return *this;
3837}
3838
3839template <class T>
3840template <RunOn run_on>
3842BaseFab<T>::invert (T const& r) noexcept
3843{
3844 return this->invert<run_on>(r, this->domain, DestComp{0}, NumComps{this->nvar});
3845}
3846
3847template <class T>
3848template <RunOn run_on>
3850BaseFab<T>::invert (T const& r, const Box& bx, DestComp dcomp, NumComps ncomp) noexcept
3851{
3852 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3853
3854 Array4<T> const& a = this->array();
3855 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3856 {
3857 a(i,j,k,n+dcomp.i) = r / a(i,j,k,n+dcomp.i);
3858 });
3859
3860 return *this;
3861}
3862
3863template <class T>
3864template <RunOn run_on>
3865T
3866BaseFab<T>::sum (const Box& bx, DestComp dcomp, NumComps ncomp) const noexcept
3867{
3868 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3869
3870 T r = 0;
3871 Array4<T const> const& a = this->const_array();
3872#ifdef AMREX_USE_GPU
3873 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
3874 ReduceOps<ReduceOpSum> reduce_op;
3875 ReduceData<T> reduce_data(reduce_op);
3876 using ReduceTuple = typename decltype(reduce_data)::Type;
3877 reduce_op.eval(bx, reduce_data,
3878 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
3879 {
3880 T t = 0;
3881 for (int n = 0; n < ncomp.n; ++n) {
3882 t += a(i,j,k,n+dcomp.i);
3883 }
3884 return { t };
3885 });
3886 ReduceTuple hv = reduce_data.value(reduce_op);
3887 r = amrex::get<0>(hv);
3888 } else
3889#endif
3890 {
3891 amrex::LoopOnCpu(bx, ncomp.n, [=,&r] (int i, int j, int k, int n) noexcept
3892 {
3893 r += a(i,j,k,n+dcomp.i);
3894 });
3895 }
3896
3897 return r;
3898}
3899
3900template <class T>
3901template <RunOn run_on>
3902T
3903BaseFab<T>::dot (const BaseFab<T>& src, const Box& bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) const noexcept
3904{
3905 AMREX_ASSERT(this->domain.sameType(src.domain));
3906 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3907 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3908
3909 T r = 0;
3910 Array4<T const> const& d = this->const_array();
3911 Array4<T const> const& s = src.const_array();
3912#ifdef AMREX_USE_GPU
3913 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
3914 ReduceOps<ReduceOpSum> reduce_op;
3915 ReduceData<T> reduce_data(reduce_op);
3916 using ReduceTuple = typename decltype(reduce_data)::Type;
3917 reduce_op.eval(bx, reduce_data,
3918 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
3919 {
3920 T t = 0;
3921 for (int n = 0; n < ncomp.n; ++n) {
3922 t += d(i,j,k,n+dcomp.i) * s(i,j,k,n+scomp.i);
3923 }
3924 return { t };
3925 });
3926 ReduceTuple hv = reduce_data.value(reduce_op);
3927 r = amrex::get<0>(hv);
3928 } else
3929#endif
3930 {
3931 amrex::LoopOnCpu(bx, ncomp.n, [=,&r] (int i, int j, int k, int n) noexcept
3932 {
3933 r += d(i,j,k,n+dcomp.i) * s(i,j,k,n+scomp.i);
3934 });
3935 }
3936
3937 return r;
3938}
3939
3940template <class T>
3941template <RunOn run_on>
3942T
3943BaseFab<T>::dot (const Box& bx, int destcomp, int numcomp) const noexcept
3944{
3945 return dot<run_on>(bx, DestComp{destcomp}, NumComps{numcomp});
3946}
3947
3948
3949template <class T>
3950template <RunOn run_on>
3951T
3952BaseFab<T>::dot (const Box& bx, DestComp dcomp, NumComps ncomp) const noexcept
3953{
3954 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3955
3956 T r = 0;
3957 Array4<T const> const& a = this->const_array();
3958#ifdef AMREX_USE_GPU
3959 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
3960 ReduceOps<ReduceOpSum> reduce_op;
3961 ReduceData<T> reduce_data(reduce_op);
3962 using ReduceTuple = typename decltype(reduce_data)::Type;
3963 reduce_op.eval(bx, reduce_data,
3964 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
3965 {
3966 T t = 0;
3967 for (int n = 0; n < ncomp.n; ++n) {
3968 t += a(i,j,k,n+dcomp.i)*a(i,j,k,n+dcomp.i);
3969 }
3970 return { t };
3971 });
3972 ReduceTuple hv = reduce_data.value(reduce_op);
3973 r = amrex::get<0>(hv);
3974 } else
3975#endif
3976 {
3977 amrex::LoopOnCpu(bx, ncomp.n, [=,&r] (int i, int j, int k, int n) noexcept
3978 {
3979 r += a(i,j,k,n+dcomp.i)*a(i,j,k,n+dcomp.i);
3980 });
3981 }
3982
3983 return r;
3984}
3985
3986template <class T>
3987template <RunOn run_on>
3988T
3989BaseFab<T>::dotmask (const BaseFab<T>& src, const Box& bx, const BaseFab<int>& mask,
3990 SrcComp scomp, DestComp dcomp, NumComps ncomp) const noexcept
3991{
3992 AMREX_ASSERT(this->domain.sameType(src.domain));
3993 AMREX_ASSERT(this->domain.sameType(mask.domain));
3994 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3995 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3996
3997 T r = 0;
3998 Array4<T const> const& d = this->const_array();
3999 Array4<T const> const& s = src.const_array();
4000 Array4<int const> const& m = mask.const_array();
4001#ifdef AMREX_USE_GPU
4002 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
4003 ReduceOps<ReduceOpSum> reduce_op;
4004 ReduceData<T> reduce_data(reduce_op);
4005 using ReduceTuple = typename decltype(reduce_data)::Type;
4006 reduce_op.eval(bx, reduce_data,
4007 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
4008 {
4009 T t = 0;
4010 T mi = static_cast<T>(static_cast<int>(static_cast<bool>(m(i,j,k))));
4011 for (int n = 0; n < ncomp.n; ++n) {
4012 t += d(i,j,k,n+dcomp.i)*s(i,j,k,n+scomp.i)*mi;
4013 }
4014 return { t };
4015 });
4016 ReduceTuple hv = reduce_data.value(reduce_op);
4017 r = amrex::get<0>(hv);
4018 } else
4019#endif
4020 {
4021 amrex::LoopOnCpu(bx, ncomp.n, [=,&r] (int i, int j, int k, int n) noexcept
4022 {
4023 int mi = static_cast<int>(static_cast<bool>(m(i,j,k)));
4024 r += d(i,j,k,n+dcomp.i)*s(i,j,k,n+scomp.i)*mi;
4025 });
4026 }
4027
4028 return r;
4029}
4030
4031}
4032
4033#endif /*BL_BASEFAB_H*/
General-purpose algorithm utilities available on both host and device.
Fixed-size array types for use on GPU and CPU.
Aligned CPU heap-backed implementation of the Arena interface.
Assertion macros used across AMReX for runtime consistency checks.
#define BL_ASSERT(EX)
Definition AMReX_BLassert.H:39
#define AMREX_ASSERT(EX)
Definition AMReX_BLassert.H:38
Iterator over all IntVects in a Box.
A list of Boxes sharing a common IndexType.
Integer-lattice boxes and helpers for defining index-space regions.
Coalescing first-fit dynamic memory arena.
Lightweight Arena-backed allocator and deleter helpers.
Compiler- and backend-specific extension macros (e.g., restrict, SIMD, inline).
#define AMREX_FORCE_INLINE
Definition AMReX_Extension.H:124
#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:1133
Array4< int const > offset
Definition AMReX_HypreMLABecLap.cpp:1129
Real * pdst
Definition AMReX_HypreMLABecLap.cpp:1130
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
Abstract base class for memory arenas.
Definition AMReX_Arena.H:153
A FortranArrayBox(FAB)-like object.
Definition AMReX_BaseFab.H:222
int maskGE(BaseFab< int > &mask, T const &val, int comp=0) const noexcept
Mark cells with value greater than or equal to val.
Definition AMReX_BaseFab.H:2575
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:2631
Array4< T const > const_array() const noexcept
Return a const-qualified Array4 view over all components.
Definition AMReX_BaseFab.H:547
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:3516
T sum(int comp, int numcomp=1) const noexcept
Returns sum of given component of FAB state vector.
Definition AMReX_BaseFab.H:2861
gpuStream_t alloc_stream
Definition AMReX_BaseFab.H:1310
Real norminfmask(const Box &subbox, const BaseFab< int > &mask, int scomp=0, int ncomp=1) const noexcept
Definition AMReX_BaseFab.H:2004
BaseFab< T > & divide(T const &val) noexcept
Scalar division on the whole domain and all components.
Definition AMReX_BaseFab.H:3748
const int * hiVect() const noexcept
Returns the upper corner of the domain.
Definition AMReX_BaseFab.H:392
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:3656
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:3042
static void Initialize()
Initialize static tracking/profiling metadata used by all BaseFabs.
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:1885
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:1940
std::size_t nBytesOwned() const noexcept
Definition AMReX_BaseFab.H:334
BaseFab< T > & copy(const BaseFab< T > &src) noexcept
Definition AMReX_BaseFab.H:3507
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:2696
BaseFab< T > & minus(T const &val) noexcept
Scalar subtraction on the whole domain and all components.
Definition AMReX_BaseFab.H:3608
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:2391
BaseFab< T > & plus(T const &val) noexcept
Scalar addition on the whole domain and all components.
Definition AMReX_BaseFab.H:3538
BaseFab< T > & mult(T const &val, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3694
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:3726
std::size_t nBytes(const Box &bx, int ncomps) const noexcept
Returns bytes used in the Box for those components.
Definition AMReX_BaseFab.H:339
void setPtr(T *p, Long sz) noexcept
Attach externally allocated storage without taking ownership.
Definition AMReX_BaseFab.H:457
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:2720
void define()
Allocates memory for the BaseFab<T>.
Definition AMReX_BaseFab.H:1596
BaseFab< T > & operator*=(T const &val) noexcept
Definition AMReX_BaseFab.H:3686
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:1767
void clear()
The function returns the BaseFab to the invalid state. The memory is freed.
Definition AMReX_BaseFab.H:1831
const IntVect & smallEnd() const noexcept
Returns the lower corner of the domain See class Box for analogue.
Definition AMReX_BaseFab.H:369
BaseFab< T > & mult(T const &r, int comp, int numcomp=1) noexcept
Scalar multiplication applied to selected components over the entire domain.
Definition AMReX_BaseFab.H:3188
BaseFab< T > & atomicAdd(const BaseFab< T > &x) noexcept
Atomic FAB addition (a[i] <- a[i] + b[i]).
Definition AMReX_BaseFab.H:2621
int maskEQ(BaseFab< int > &mask, T const &val, int comp=0) const noexcept
Mark cells with value equal to val.
Definition AMReX_BaseFab.H:2483
const int * loVect() const noexcept
Returns the lower corner of the domain.
Definition AMReX_BaseFab.H:382
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:398
bool isAllocated() const noexcept
Returns true if the data for the FAB has been allocated.
Definition AMReX_BaseFab.H:577
std::unique_ptr< T, DataDeleter > release() noexcept
Release ownership of memory.
Definition AMReX_BaseFab.H:1867
void setVal(T const &x, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3424
BaseFab< T > & operator-=(T const &val) noexcept
Definition AMReX_BaseFab.H:3616
const Box & box() const noexcept
Returns the domain (box) where the array is defined.
Definition AMReX_BaseFab.H:357
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:3445
static void Finalize()
Release static tracking/profiling metadata allocated by Initialize().
Array4< T > array() noexcept
Create a mutable Array4 view over all components.
Definition AMReX_BaseFab.H:510
IntVect indexFromValue(const Box &subbox, int comp, T const &value) const noexcept
Definition AMReX_BaseFab.H:2300
BaseFab< T > & mult(const BaseFab< T > &src) noexcept
Definition AMReX_BaseFab.H:3710
bool shared_memory
Is the memory allocated in shared memory?
Definition AMReX_BaseFab.H:1308
int maskLE(BaseFab< int > &mask, T const &val, int comp=0) const noexcept
Mark cells with value <= val.
Definition AMReX_BaseFab.H:2437
void setValIf(T const &val, const BaseFab< int > &mask) noexcept
Definition AMReX_BaseFab.H:3437
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:3586
void setValIfNot(T const &val, const BaseFab< int > &mask) noexcept
Definition AMReX_BaseFab.H:3459
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:2668
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:1407
std::size_t nBytes() const noexcept
Returns how many bytes used.
Definition AMReX_BaseFab.H:332
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:1912
BaseFab< T > & negate() noexcept
on the whole domain and all components
Definition AMReX_BaseFab.H:3818
BaseFab< T > & minus(const BaseFab< T > &src) noexcept
Definition AMReX_BaseFab.H:3640
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:3146
T value_type
Definition AMReX_BaseFab.H:227
void SetBoxType(const IndexType &typ) noexcept
Change the Box type without change the length.
Definition AMReX_BaseFab.H:1120
Array4< T const > array() const noexcept
Create an Array4 view over all components.
Definition AMReX_BaseFab.H:475
T maxabs(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2261
Elixir elixir() noexcept
Transfer owned storage to an Elixir when running in a GPU launch region.
Definition AMReX_BaseFab.H:1810
BaseFab< T > & operator+=(T const &val) noexcept
Definition AMReX_BaseFab.H:3546
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:3481
BaseFab< T > & minus(T const &val, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3624
Long truesize
nvar*numpts that was allocated on heap.
Definition AMReX_BaseFab.H:1306
void setVal(T const &val) noexcept
Set value on the whole domain and all components.
Definition AMReX_BaseFab.H:3416
const int * nCompPtr() const noexcept
for calls to fortran.
Definition AMReX_BaseFab.H:346
Array4< T const > const_array(int start_comp, int num_comps) const noexcept
Return a const-qualified Array4 view spanning num_comps components.
Definition AMReX_BaseFab.H:571
Box domain
My index space.
Definition AMReX_BaseFab.H:1304
void fill_snan() noexcept
Definition AMReX_BaseFab.H:1514
bool contains(const Box &bx) const noexcept
Returns true if bx is totally contained within the domain of this BaseFab.
Definition AMReX_BaseFab.H:407
T * dptr
The data pointer.
Definition AMReX_BaseFab.H:1303
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:1480
BaseFab< T > & divide(T const &val, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3764
T max(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2181
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:1554
Array4< T > array(int start_comp, int num_comps) noexcept
Mutable overload spanning num_comps components.
Definition AMReX_BaseFab.H:535
int nvar
Number components.
Definition AMReX_BaseFab.H:1305
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:2757
BaseFab< T > & divide(T const &r, int comp, int numcomp=1) noexcept
Divide selected components by a scalar across the entire domain.
Definition AMReX_BaseFab.H:3246
Array4< T const > array(int start_comp) const noexcept
Create a sliced Array4 view beginning at component start_comp.
Definition AMReX_BaseFab.H:487
BaseFab< T > & operator/=(T const &val) noexcept
Definition AMReX_BaseFab.H:3756
std::pair< T, T > minmax(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2219
Array4< T const > const_array(int start_comp) const noexcept
Return a const-qualified Array4 view beginning at component start_comp.
Definition AMReX_BaseFab.H:559
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:1538
Long size() const noexcept
Returns the total number of points of all components.
Definition AMReX_BaseFab.H:354
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:2917
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:1454
const IntVect & bigEnd() const noexcept
Returns the upper corner of the domain. See class Box for analogue.
Definition AMReX_BaseFab.H:372
Array4< T > array(int start_comp) noexcept
Mutable overload selecting components starting at start_comp.
Definition AMReX_BaseFab.H:522
Long numPts() const noexcept
Returns the number of points.
Definition AMReX_BaseFab.H:351
const T * dataPtr(int n=0) const noexcept
Const-qualified variant of dataPtr().
Definition AMReX_BaseFab.H:427
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:3467
BaseFab< T > & mult(T const &val) noexcept
Scalar multiplication on the whole domain and all components.
Definition AMReX_BaseFab.H:3678
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:3796
void setValIfNot(T const &val, const Box &bx, const BaseFab< int > &mask, int nstart, int num) noexcept
Definition AMReX_BaseFab.H:1546
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:1504
void setVal(T const &x, const Box &bx, int N=0) noexcept
Set a single component N over bx to x.
Definition AMReX_BaseFab.H:1522
void prefetchToDevice() const noexcept
Advise the runtime to migrate storage to the active device.
Definition AMReX_BaseFab.H:1374
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:418
bool ptr_owner
Owner of T*?
Definition AMReX_BaseFab.H:1307
virtual ~BaseFab() noexcept
The destructor deletes the array memory.
Definition AMReX_BaseFab.H:1714
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:363
IntVect maxIndex(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2365
BaseFab< T > & protected_divide(const BaseFab< T > &src) noexcept
Divide wherever "src" is "true" or "non-zero".
Definition AMReX_BaseFab.H:3304
friend class BaseFab
Definition AMReX_BaseFab.H:225
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:2901
Array4< T const > array(int start_comp, int num_comps) const noexcept
Create a sliced Array4 view spanning num_comps components.
Definition AMReX_BaseFab.H:500
int nComp() const noexcept
Returns the number of components.
Definition AMReX_BaseFab.H:343
int maskGT(BaseFab< int > &mask, T const &val, int comp=0) const noexcept
Mark cells with value greater than val.
Definition AMReX_BaseFab.H:2529
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:2807
BaseFab< T > & plus(const BaseFab< T > &src) noexcept
Definition AMReX_BaseFab.H:3570
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:3780
void prefetchToHost() const noexcept
Advise the runtime to migrate storage to the host.
Definition AMReX_BaseFab.H:1342
T min(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2143
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:1968
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:3374
IntVect minIndex(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2339
T * dataPtr(const IntVect &p, int n=0) noexcept
Return a pointer to the element at spatial index p and component n.
Definition AMReX_BaseFab.H:1317
void abs() noexcept
Compute absolute value for all components of this FAB.
Definition AMReX_BaseFab.H:1976
void getVal(T *data, const IntVect &pos) const noexcept
Copy all components beginning at component 0.
Definition AMReX_BaseFab.H:1472
BaseFab< T > & plus(T const &val, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3554
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:2048
A list of Boxes sharing a common IndexType.
Definition AMReX_BoxList.H:109
__host__ __device__ const IntVectND< dim > & bigEnd() const &noexcept
Return the inclusive upper bound of the box.
Definition AMReX_Box.H:136
__host__ __device__ const int * hiVect() const &noexcept
Return a pointer to the high-end coordinates (useful for Fortran calls).
Definition AMReX_Box.H:210
__host__ __device__ Long numPts() const noexcept
Return the number of points contained in the BoxND.
Definition AMReX_Box.H:385
__host__ __device__ bool sameSize(const BoxND &b) const noexcept
Return true is Boxes same size, ie translates of each other,. It is an error if they have different t...
Definition AMReX_Box.H:312
__host__ __device__ IntVectND< dim > length() const noexcept
Return the length of the BoxND.
Definition AMReX_Box.H:167
__host__ __device__ bool contains(const IntVectND< dim > &p) const noexcept
Return true if argument is contained within BoxND.
Definition AMReX_Box.H:233
__host__ __device__ const int * loVect() const &noexcept
Return a pointer to the low-end coordinates (useful for Fortran calls).
Definition AMReX_Box.H:205
__host__ __device__ BoxND & setType(const IndexTypeND< dim > &t) noexcept
Set indexing type.
Definition AMReX_Box.H:534
__host__ __device__ bool ok() const noexcept
Return true if high bounds are >= low bounds and the index type is valid.
Definition AMReX_Box.H:229
__host__ __device__ const IntVectND< dim > & smallEnd() const &noexcept
Return the inclusive lower bound of the box.
Definition AMReX_Box.H:124
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:438
Type value()
Definition AMReX_Reduce.H:473
Definition AMReX_Reduce.H:597
void eval(MF const &mf, IntVect const &nghost, D &reduce_data, F &&f)
Definition AMReX_Reduce.H:731
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:1364
__host__ __device__ Dim3 length(Array4< T > const &a) noexcept
Return the spatial extents of an Array4 in Dim3 form.
Definition AMReX_Array4.H:1378
__host__ __device__ Dim3 lbound(Array4< T > const &a) noexcept
Return the inclusive lower bounds of an Array4 in Dim3 form.
Definition AMReX_Array4.H:1350
std::array< T, N > Array
Definition AMReX_Array.H:31
__host__ __device__ constexpr const T & min(const T &a, const T &b) noexcept
Definition AMReX_Algorithm.H:31
__host__ __device__ constexpr const T & max(const T &a, const T &b) noexcept
Definition AMReX_Algorithm.H:53
__host__ __device__ bool almostEqual(T x, T y, int ulp=2)
Definition AMReX_Algorithm.H:139
__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:640
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 &...)
No-op helper that marks variables as intentionally unused.
Definition AMReX.H:259
int nComp(FabArrayBase const &fa)
Definition AMReX_FabArrayBase.cpp:2852
__host__ __device__ Array4< T > makeArray4(T *p, Box const &bx, int ncomp) noexcept
Definition AMReX_BaseFab.H:123
RunOn
Definition AMReX_GpuControl.H:65
__host__ __device__ Dim3 begin(BoxND< dim > const &box) noexcept
Return the iterator begin coordinate of box as Dim3.
Definition AMReX_Box.H:2239
cudaStream_t gpuStream_t
Definition AMReX_GpuControl.H:79
bool InitSNaN() noexcept
Whether new allocations are filled with signaling NaNs when enabled.
Definition AMReX.cpp:185
void ParallelFor(TypeList< CTOs... > ctos, std::array< int, sizeof...(CTOs)> const &runtime_options, T N, F &&f)
Definition AMReX_CTOParallelForImpl.H:202
Long TotalBytesAllocatedInFabs() noexcept
Total bytes currently held by all BaseFab instances on this rank.
Definition AMReX_BaseFab.cpp:66
BoxND< 3 > Box
Box is an alias for amrex::BoxND instantiated with AMREX_SPACEDIM.
Definition AMReX_BaseFwd.H:35
void BaseFab_Initialize()
Initialize BaseFab subsystem resources (profilers, stats counters, etc.).
Definition AMReX_BaseFab.cpp:30
void BaseFab_Finalize()
Release any BaseFab subsystem resources.
Definition AMReX_BaseFab.cpp:59
void ResetTotalBytesAllocatedInFabsHWM() noexcept
Reset the byte high-water mark reported by TotalBytesAllocatedInFabsHWM().
Definition AMReX_BaseFab.cpp:134
BoxList boxDiff(const Box &b1in, const Box &b2)
Returns a BoxList defining the complement of b2 in b1in.
Definition AMReX_BoxList.cpp:597
void placementNew(T *const ptr, Long n)
Start the lifetime of objects in a preallocated buffer.
Definition AMReX_BaseFab.H:137
IntVectND< 3 > IntVect
IntVect is an alias for amrex::IntVectND instantiated with AMREX_SPACEDIM.
Definition AMReX_BaseFwd.H:38
Long TotalBytesAllocatedInFabsHWM() noexcept
High-water mark of BaseFab bytes since the last reset.
Definition AMReX_BaseFab.cpp:83
void placementDelete(T *const ptr, Long n)
Destroy objects previously constructed via placementNew.
Definition AMReX_BaseFab.H:158
Long TotalCellsAllocatedInFabsHWM() noexcept
High-water mark of Real BaseFab grid-point count since the last reset.
Definition AMReX_BaseFab.cpp:117
Long TotalCellsAllocatedInFabs() noexcept
Total spatial grid points (numPts) tracked for Real BaseFabs on this rank.
Definition AMReX_BaseFab.cpp:100
void Error(const std::string &msg)
Print a message to stderr and abort the program.
Definition AMReX.cpp:235
void Abort(const std::string &msg)
Print a fatal-error message to stderr and abort execution.
Definition AMReX.cpp:241
void LoopOnCpu(Dim3 lo, Dim3 hi, F const &f) noexcept
Definition AMReX_Loop.H:365
__host__ __device__ Dim3 end(BoxND< dim > const &box) noexcept
Return the iterator end coordinate of box as Dim3.
Definition AMReX_Box.H:2257
void update_fab_stats(Long n, Long s, size_t szt) noexcept
Definition AMReX_BaseFab.cpp:146
A multidimensional array accessor.
Definition AMReX_Array4.H:288
__host__ __device__ constexpr std::size_t size() const noexcept
Total number of elements in the ArrayND's index region.
Definition AMReX_Array4.H:642
__host__ __device__ T * ptr(idx... i) const noexcept
Multi-index ptr() for accessing pointer to element.
Definition AMReX_Array4.H:560
Lightweight allocator/deallocator backed by an Arena.
Definition AMReX_DataAllocator.H:21
void * alloc(std::size_t sz) const noexcept
Allocate sz bytes from the arena.
Definition AMReX_DataAllocator.H:31
Arena * arena() const noexcept
Return the stored arena, or The_Arena() if none was set.
Definition AMReX_DataAllocator.H:41
Custom deleter that returns memory to an Arena.
Definition AMReX_DataAllocator.H:53
Destination-component descriptor.
Definition AMReX_BaseFab.H:107
int i
Definition AMReX_BaseFab.H:110
__host__ __device__ DestComp(int ai) noexcept
Definition AMReX_BaseFab.H:109
A simple struct holding 3 int values for a 3D index.
Definition AMReX_Dim3.H:24
int x
Definition AMReX_Dim3.H:24
Definition AMReX_TypeTraits.H:51
Number-of-components descriptor.
Definition AMReX_BaseFab.H:114
__host__ __device__ NumComps(int an) noexcept
Definition AMReX_BaseFab.H:116
int n
Definition AMReX_BaseFab.H:117
Source-component descriptor.
Definition AMReX_BaseFab.H:100
__host__ __device__ SrcComp(int ai) noexcept
Definition AMReX_BaseFab.H:102
int i
Definition AMReX_BaseFab.H:103