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 Long 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#ifdef AMREX_USE_GPU
1662 if (Gpu::inLaunchRegion() && arena()->isDeviceAccessible()) {
1663 this->copy<RunOn::Device>(rhs, this->domain, scomp, this->domain, 0, ncomp);
1664 } else
1665#endif
1666 {
1667 this->copy<RunOn::Host>(rhs, this->domain, scomp, this->domain, 0, ncomp);
1668 }
1669 } else if (make_type == amrex::make_alias) {
1670 ; // nothing to do
1671 } else {
1672 amrex::Abort("BaseFab: unknown MakeType");
1673 }
1674}
1675
1676template<class T>
1677BaseFab<T>::BaseFab (const Box& bx, int ncomp, T* p)
1678 : dptr(p), domain(bx), nvar(ncomp), truesize(bx.numPts()*ncomp)
1679{
1680}
1681
1682template<class T>
1683BaseFab<T>::BaseFab (const Box& bx, int ncomp, T const* p)
1684 : dptr(const_cast<T*>(p)), domain(bx), nvar(ncomp), truesize(bx.numPts()*ncomp)
1685{
1686}
1687
1688template<class T>
1690 : dptr(a.p),
1691 domain(IntVect(AMREX_D_DECL(a.begin[0],a.begin[1],a.begin[2])),
1692 IntVect(AMREX_D_DECL(a.end[0]-1,a.end[1]-1,a.end[2]-1))),
1693 nvar(a.nComp()), truesize(a.size())
1694{}
1695
1696template<class T>
1698 : dptr(a.p),
1699 domain(IntVect(AMREX_D_DECL(a.begin[0],a.begin[1],a.begin[2])),
1700 IntVect(AMREX_D_DECL(a.end[0]-1,a.end[1]-1,a.end[2]-1)), t),
1701 nvar(a.nComp()), truesize(a.size())
1702{}
1703
1704template<class T>
1706 : dptr(const_cast<T*>(a.p)),
1707 domain(IntVect(AMREX_D_DECL(a.begin[0],a.begin[1],a.begin[2])),
1708 IntVect(AMREX_D_DECL(a.end[0]-1,a.end[1]-1,a.end[2]-1))),
1709 nvar(a.nComp()), truesize(a.size())
1710{}
1711
1712template<class T>
1714 : dptr(const_cast<T*>(a.p)),
1715 domain(IntVect(AMREX_D_DECL(a.begin[0],a.begin[1],a.begin[2])),
1716 IntVect(AMREX_D_DECL(a.end[0]-1,a.end[1]-1,a.end[2]-1)), t),
1717 nvar(a.nComp()), truesize(a.size())
1718{}
1719
1720template <class T>
1722{
1723 clear();
1724}
1725
1726template <class T>
1728 : DataAllocator{rhs.arena()},
1729 dptr(rhs.dptr), domain(rhs.domain),
1730 nvar(rhs.nvar), truesize(rhs.truesize),
1731 ptr_owner(rhs.ptr_owner), shared_memory(rhs.shared_memory)
1732#ifdef AMREX_USE_GPU
1733 , alloc_stream(rhs.alloc_stream)
1734#endif
1735{
1736 rhs.dptr = nullptr;
1737 rhs.ptr_owner = false;
1738}
1739
1740template <class T>
1741BaseFab<T>&
1743{
1744 if (this != &rhs) {
1745 clear();
1746 DataAllocator::operator=(rhs);
1747 dptr = rhs.dptr;
1748 domain = rhs.domain;
1749 nvar = rhs.nvar;
1750 truesize = rhs.truesize;
1751 ptr_owner = rhs.ptr_owner;
1752 shared_memory = rhs.shared_memory;
1753#ifdef AMREX_USE_GPU
1754 alloc_stream = rhs.alloc_stream;
1755#endif
1756
1757 rhs.dptr = nullptr;
1758 rhs.ptr_owner = false;
1759 }
1760 return *this;
1761}
1762
1763template <class T>
1764template <RunOn run_on>
1766BaseFab<T>::operator= (T const& t) noexcept
1767{
1768 setVal<run_on>(t);
1769 return *this;
1770}
1771
1772template <class T>
1773void
1774BaseFab<T>::resize (const Box& b, int n, Arena* ar)
1775{
1776 this->nvar = n;
1777 this->domain = b;
1778
1779 if (ar == nullptr) {
1780 ar = m_arena;
1781 }
1782
1783 if (arena() != DataAllocator(ar).arena()) {
1784 clear();
1785 m_arena = ar;
1786 define();
1787 }
1788 else if (this->dptr == nullptr || !this->ptr_owner)
1789 {
1790 if (this->shared_memory) {
1791 amrex::Abort("BaseFab::resize: BaseFab in shared memory cannot increase size");
1792 }
1793
1794 this->dptr = nullptr;
1795 define();
1796 }
1797 else if (this->nvar*this->domain.numPts() > this->truesize
1798#ifdef AMREX_USE_GPU
1799 || (arena()->isStreamOrderedArena() && alloc_stream != Gpu::gpuStream())
1800#endif
1801 )
1802 {
1803 if (this->shared_memory) {
1804 amrex::Abort("BaseFab::resize: BaseFab in shared memory cannot increase size");
1805 }
1806
1807 clear();
1808
1809 define();
1810 }
1811}
1812
1813template <class T>
1814template <class U>
1815requires (std::is_trivially_destructible_v<U>)
1816Elixir
1818{
1819 bool o;
1820 if (Gpu::inLaunchRegion()) {
1821 o = this->ptr_owner;
1822 this->ptr_owner = false;
1823 if (o && this->dptr) {
1824 if (this->nvar > 0) {
1825 amrex::update_fab_stats(-this->truesize/this->nvar, -this->truesize, sizeof(T));
1826 } else {
1827 amrex::update_fab_stats(0, -this->truesize, sizeof(T));
1828 }
1829 }
1830 } else {
1831 o = false;
1832 }
1833 return Elixir((o ? this->dptr : nullptr), this->arena());
1834}
1835
1836template <class T>
1837void
1839{
1840 if (this->dptr)
1841 {
1842 //
1843 // Call T::~T() on the to-be-destroyed memory.
1844 //
1845 if (this->ptr_owner)
1846 {
1847 if (this->shared_memory)
1848 {
1849 amrex::Abort("BaseFab::clear: BaseFab cannot be owner of shared memory");
1850 }
1851
1852 placementDelete(this->dptr, this->truesize);
1853
1854#ifdef AMREX_USE_GPU
1855 this->arena()->streamOrderedFree(this->dptr, alloc_stream);
1856#else
1857 this->free(this->dptr);
1858#endif
1859
1860 if (this->nvar > 0) {
1861 amrex::update_fab_stats(-this->truesize/this->nvar, -this->truesize, sizeof(T));
1862 } else {
1863 amrex::update_fab_stats(0, -this->truesize, sizeof(T));
1864 }
1865 }
1866
1867 this->dptr = nullptr;
1868 this->truesize = 0;
1869 }
1870}
1871
1872template <class T>
1873std::unique_ptr<T,DataDeleter>
1875{
1876 std::unique_ptr<T,DataDeleter> r(nullptr, DataDeleter{this->arena()});
1877 if (this->dptr && this->ptr_owner) {
1878 r.reset(this->dptr);
1879 this->ptr_owner = false;
1880 if (this->nvar > 0) {
1881 amrex::update_fab_stats(-this->truesize/this->nvar, -this->truesize, sizeof(T));
1882 } else {
1883 amrex::update_fab_stats(0, -this->truesize, sizeof(T));
1884 }
1885 }
1886 return r;
1887}
1888
1889template <class T>
1890template <RunOn run_on>
1891std::size_t
1893 int srccomp,
1894 int numcomp,
1895 void* dst) const noexcept
1896{
1897 BL_ASSERT(box().contains(srcbox));
1898 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= nComp());
1899
1900 if (srcbox.ok())
1901 {
1902 Array4<T> d(static_cast<T*>(dst),amrex::begin(srcbox),amrex::end(srcbox),numcomp);
1903 Array4<T const> const& s = this->const_array();
1904 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, srcbox, numcomp, i, j, k, n,
1905 {
1906 d(i,j,k,n) = s(i,j,k,n+srccomp);
1907 });
1908 return sizeof(T)*d.size();
1909 }
1910 else
1911 {
1912 return 0;
1913 }
1914}
1915
1916template <class T>
1917template <RunOn run_on, typename BUF>
1918std::size_t
1920 int dstcomp,
1921 int numcomp,
1922 const void* src) noexcept
1923{
1924 BL_ASSERT(box().contains(dstbox));
1925 BL_ASSERT(dstcomp >= 0 && dstcomp+numcomp <= nComp());
1926
1927 if (dstbox.ok())
1928 {
1929 Array4<BUF const> s(static_cast<BUF const*>(src), amrex::begin(dstbox),
1930 amrex::end(dstbox), numcomp);
1931 Array4<T> const& d = this->array();
1932 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, dstbox, numcomp, i, j, k, n,
1933 {
1934 d(i,j,k,n+dstcomp) = static_cast<T>(s(i,j,k,n));
1935 });
1936 return sizeof(BUF)*s.size();
1937 }
1938 else
1939 {
1940 return 0;
1941 }
1942}
1943
1944template <class T>
1945template <RunOn run_on, typename BUF>
1946std::size_t
1948 int dstcomp,
1949 int numcomp,
1950 const void* src) noexcept
1951{
1952 BL_ASSERT(box().contains(dstbox));
1953 BL_ASSERT(dstcomp >= 0 && dstcomp+numcomp <= nComp());
1954
1955 if (dstbox.ok())
1956 {
1957 Array4<BUF const> s(static_cast<BUF const*>(src), amrex::begin(dstbox),
1958 amrex::end(dstbox), numcomp);
1959 Array4<T> const& d = this->array();
1960 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, dstbox, numcomp, i, j, k, n,
1961 {
1962 d(i,j,k,n+dstcomp) += static_cast<T>(s(i,j,k,n));
1963 });
1964 return sizeof(BUF)*s.size();
1965 }
1966 else
1967 {
1968 return 0;
1969 }
1970}
1971
1972template <class T>
1973template <RunOn run_on>
1974void
1975BaseFab<T>::setComplement (T const& x, const Box& b, int ns, int num) noexcept
1976{
1977 this->setComplement<run_on>(x, b, DestComp{ns}, NumComps{num});
1978}
1979
1980template <class T>
1981template <RunOn run_on>
1982void
1984{
1985 this->abs<run_on>(this->domain,0,this->nvar);
1986}
1987
1988template <class T>
1989template <RunOn run_on>
1990void
1991BaseFab<T>::abs (int comp, int numcomp) noexcept
1992{
1993 this->abs<run_on>(this->domain,comp,numcomp);
1994}
1995
1996template <class T>
1997template <RunOn run_on>
1998void
1999BaseFab<T>::abs (const Box& subbox, int comp, int numcomp) noexcept
2000{
2001 Array4<T> const& a = this->array();
2002 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG (run_on, subbox, numcomp, i, j, k, n,
2003 {
2004 a(i,j,k,n+comp) = std::abs(a(i,j,k,n+comp));
2005 });
2006}
2007
2008template <class T>
2009template <RunOn run_on>
2010Real
2012 int scomp, int ncomp) const noexcept
2013{
2014 BL_ASSERT(this->domain.contains(subbox));
2015 BL_ASSERT(scomp >= 0 && scomp + ncomp <= this->nvar);
2016
2017 Array4<T const> const& a = this->const_array();
2018 Array4<int const> const& m = mask.const_array();
2019 Real r = 0.0;
2020#ifdef AMREX_USE_GPU
2021 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2022 ReduceOps<ReduceOpMax> reduce_op;
2023 ReduceData<Real> reduce_data(reduce_op);
2024 using ReduceTuple = ReduceData<Real>::Type;
2025 reduce_op.eval(subbox, reduce_data,
2026 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2027 {
2028 Real t = 0.0;
2029 if (m(i,j,k)) {
2030 for (int n = 0; n < ncomp; ++n) {
2031 t = amrex::max(t,static_cast<Real>(std::abs(a(i,j,k,n+scomp))));
2032 }
2033 }
2034 return {t};
2035 });
2036 ReduceTuple hv = reduce_data.value(reduce_op);
2037 r = amrex::get<0>(hv);
2038 } else
2039#endif
2040 {
2041 amrex::LoopOnCpu(subbox, ncomp, [=,&r] (int i, int j, int k, int n) noexcept
2042 {
2043 if (m(i,j,k)) {
2044 Real t = static_cast<Real>(std::abs(a(i,j,k,n+scomp)));
2045 r = amrex::max(r,t);
2046 }
2047 });
2048 }
2049 return r;
2050}
2051
2052template <class T>
2053template <RunOn run_on>
2054Real
2055BaseFab<T>::norm (int p, int comp, int numcomp) const
2056{
2057 return norm<run_on>(this->domain,p,comp,numcomp);
2058}
2059
2060template <class T>
2061template <RunOn run_on>
2062Real
2063BaseFab<T>::norm (const Box& subbox, int p, int comp, int numcomp) const
2064{
2065 BL_ASSERT(this->domain.contains(subbox));
2066 BL_ASSERT(comp >= 0 && comp + numcomp <= this->nvar);
2067
2068 Array4<T const> const& a = this->const_array();
2069 Real nrm = 0.;
2070#ifdef AMREX_USE_GPU
2071 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2072 if (p == 0) {
2073 ReduceOps<ReduceOpMax> reduce_op;
2074 ReduceData<Real> reduce_data(reduce_op);
2075 using ReduceTuple = ReduceData<Real>::Type;
2076 reduce_op.eval(subbox, reduce_data,
2077 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2078 {
2079 Real t = 0.0;
2080 for (int n = 0; n < numcomp; ++n) {
2081 t = amrex::max(t, static_cast<Real>(std::abs(a(i,j,k,n+comp))));
2082 }
2083 return {t};
2084 });
2085 ReduceTuple hv = reduce_data.value(reduce_op);
2086 nrm = amrex::get<0>(hv);
2087 } else if (p == 1) {
2088 ReduceOps<ReduceOpSum> reduce_op;
2089 ReduceData<Real> reduce_data(reduce_op);
2090 using ReduceTuple = ReduceData<Real>::Type;
2091 reduce_op.eval(subbox, reduce_data,
2092 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2093 {
2094 Real t = 0.0;
2095 for (int n = 0; n < numcomp; ++n) {
2096 t += static_cast<Real>(std::abs(a(i,j,k,n+comp)));
2097 }
2098 return {t};
2099 });
2100 ReduceTuple hv = reduce_data.value(reduce_op);
2101 nrm = amrex::get<0>(hv);
2102 } else if (p == 2) {
2103 ReduceOps<ReduceOpSum> reduce_op;
2104 ReduceData<Real> reduce_data(reduce_op);
2105 using ReduceTuple = ReduceData<Real>::Type;
2106 reduce_op.eval(subbox, reduce_data,
2107 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2108 {
2109 Real t = 0.0;
2110 for (int n = 0; n < numcomp; ++n) {
2111 t += static_cast<Real>(a(i,j,k,n+comp)*a(i,j,k,n+comp));
2112 }
2113 return {t};
2114 });
2115 ReduceTuple hv = reduce_data.value(reduce_op);
2116 nrm = amrex::get<0>(hv);
2117 } else {
2118 amrex::Error("BaseFab<T>::norm: wrong p");
2119 }
2120 } else
2121#endif
2122 {
2123 if (p == 0) {
2124 amrex::LoopOnCpu(subbox, numcomp, [=,&nrm] (int i, int j, int k, int n) noexcept
2125 {
2126 Real t = static_cast<Real>(std::abs(a(i,j,k,n+comp)));
2127 nrm = amrex::max(nrm,t);
2128 });
2129 } else if (p == 1) {
2130 amrex::LoopOnCpu(subbox, numcomp, [=,&nrm] (int i, int j, int k, int n) noexcept
2131 {
2132 nrm += std::abs(a(i,j,k,n+comp));
2133 });
2134 } else if (p == 2) {
2135 amrex::LoopOnCpu(subbox, numcomp, [=,&nrm] (int i, int j, int k, int n) noexcept
2136 {
2137 nrm += a(i,j,k,n+comp)*a(i,j,k,n+comp);
2138 });
2139 } else {
2140 amrex::Error("BaseFab<T>::norm: wrong p");
2141 }
2142 }
2143
2144 return nrm;
2145}
2146
2147template <class T>
2148template <RunOn run_on>
2149T
2150BaseFab<T>::min (int comp) const noexcept
2151{
2152 return this->min<run_on>(this->domain,comp);
2153}
2154
2155template <class T>
2156template <RunOn run_on>
2157T
2158BaseFab<T>::min (const Box& subbox, int comp) const noexcept
2159{
2160 Array4<T const> const& a = this->const_array(comp);
2161#ifdef AMREX_USE_GPU
2162 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2163 ReduceOps<ReduceOpMin> reduce_op;
2164 ReduceData<T> reduce_data(reduce_op);
2165 using ReduceTuple = typename decltype(reduce_data)::Type;
2166 reduce_op.eval(subbox, reduce_data,
2167 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2168 {
2169 return { a(i,j,k) };
2170 });
2171 ReduceTuple hv = reduce_data.value(reduce_op);
2172 return amrex::get<0>(hv);
2173 } else
2174#endif
2175 {
2176 T r = std::numeric_limits<T>::max();
2177 amrex::LoopOnCpu(subbox, [=,&r] (int i, int j, int k) noexcept
2178 {
2179 r = amrex::min(r, a(i,j,k));
2180 });
2181 return r;
2182 }
2183}
2184
2185template <class T>
2186template <RunOn run_on>
2187T
2188BaseFab<T>::max (int comp) const noexcept
2189{
2190 return this->max<run_on>(this->domain,comp);
2191}
2192
2193template <class T>
2194template <RunOn run_on>
2195T
2196BaseFab<T>::max (const Box& subbox, int comp) const noexcept
2197{
2198 Array4<T const> const& a = this->const_array(comp);
2199#ifdef AMREX_USE_GPU
2200 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2201 ReduceOps<ReduceOpMax> reduce_op;
2202 ReduceData<T> reduce_data(reduce_op);
2203 using ReduceTuple = typename decltype(reduce_data)::Type;
2204 reduce_op.eval(subbox, reduce_data,
2205 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2206 {
2207 return { a(i,j,k) };
2208 });
2209 ReduceTuple hv = reduce_data.value(reduce_op);
2210 return amrex::get<0>(hv);
2211 } else
2212#endif
2213 {
2214 T r = std::numeric_limits<T>::lowest();
2215 amrex::LoopOnCpu(subbox, [=,&r] (int i, int j, int k) noexcept
2216 {
2217 r = amrex::max(r, a(i,j,k));
2218 });
2219 return r;
2220 }
2221}
2222
2223template <class T>
2224template <RunOn run_on>
2225std::pair<T,T>
2226BaseFab<T>::minmax (int comp) const noexcept
2227{
2228 return this->minmax<run_on>(this->domain,comp);
2229}
2230
2231template <class T>
2232template <RunOn run_on>
2233std::pair<T,T>
2234BaseFab<T>::minmax (const Box& subbox, int comp) const noexcept
2235{
2236 Array4<T const> const& a = this->const_array(comp);
2237#ifdef AMREX_USE_GPU
2238 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2240 ReduceData<T,T> reduce_data(reduce_op);
2241 using ReduceTuple = typename decltype(reduce_data)::Type;
2242 reduce_op.eval(subbox, reduce_data,
2243 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2244 {
2245 auto const x = a(i,j,k);
2246 return { x, x };
2247 });
2248 ReduceTuple hv = reduce_data.value(reduce_op);
2249 return std::make_pair(amrex::get<0>(hv), amrex::get<1>(hv));
2250 } else
2251#endif
2252 {
2253 T rmax = std::numeric_limits<T>::lowest();
2254 T rmin = std::numeric_limits<T>::max();
2255 amrex::LoopOnCpu(subbox, [=,&rmin,&rmax] (int i, int j, int k) noexcept
2256 {
2257 auto const x = a(i,j,k);
2258 rmin = amrex::min(rmin, x);
2259 rmax = amrex::max(rmax, x);
2260 });
2261 return std::make_pair(rmin,rmax);
2262 }
2263}
2264
2265template <class T>
2266template <RunOn run_on>
2267T
2268BaseFab<T>::maxabs (int comp) const noexcept
2269{
2270 return this->maxabs<run_on>(this->domain,comp);
2271}
2272
2273template <class T>
2274template <RunOn run_on>
2275T
2276BaseFab<T>::maxabs (const Box& subbox, int comp) const noexcept
2277{
2278 Array4<T const> const& a = this->const_array(comp);
2279#ifdef AMREX_USE_GPU
2280 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2281 ReduceOps<ReduceOpMax> reduce_op;
2282 ReduceData<T> reduce_data(reduce_op);
2283 using ReduceTuple = typename decltype(reduce_data)::Type;
2284 reduce_op.eval(subbox, reduce_data,
2285 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2286 {
2287 return { std::abs(a(i,j,k)) };
2288 });
2289 ReduceTuple hv = reduce_data.value(reduce_op);
2290 return amrex::get<0>(hv);
2291 } else
2292#endif
2293 {
2294 T r = 0;
2295 amrex::LoopOnCpu(subbox, [=,&r] (int i, int j, int k) noexcept
2296 {
2297 r = amrex::max(r, std::abs(a(i,j,k)));
2298 });
2299 return r;
2300 }
2301}
2302
2303
2304template <class T>
2305template <RunOn run_on>
2306IntVect
2307BaseFab<T>::indexFromValue (Box const& subbox, int comp, T const& value) const noexcept
2308{
2309 Array4<T const> const& a = this->const_array(comp);
2310#ifdef AMREX_USE_GPU
2311 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2312 Array<int,1+AMREX_SPACEDIM> ha{0,AMREX_D_DECL(std::numeric_limits<int>::lowest(),
2313 std::numeric_limits<int>::lowest(),
2314 std::numeric_limits<int>::lowest())};
2315 Gpu::DeviceVector<int> dv(1+AMREX_SPACEDIM);
2316 int* p = dv.data();
2317 Gpu::htod_memcpy_async(p, ha.data(), sizeof(int)*ha.size());
2318 amrex::ParallelFor(subbox, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
2319 {
2320 int* flag = p;
2321 if ((*flag == 0) && (a(i,j,k) == value)) {
2322 if (Gpu::Atomic::Exch(flag,1) == 0) {
2323 AMREX_D_TERM(p[1] = i;,
2324 p[2] = j;,
2325 p[3] = k;);
2326 }
2327 }
2328 });
2329 Gpu::dtoh_memcpy_async(ha.data(), p, sizeof(int)*ha.size());
2331 return IntVect(AMREX_D_DECL(ha[1],ha[2],ha[3]));
2332 } else
2333#endif
2334 {
2335 AMREX_LOOP_3D(subbox, i, j, k,
2336 {
2337 if (a(i,j,k) == value) { return IntVect(AMREX_D_DECL(i,j,k)); }
2338 });
2339 return IntVect::TheMinVector();
2340 }
2341}
2342
2343template <class T>
2344template <RunOn run_on>
2345IntVect
2346BaseFab<T>::minIndex (int comp) const noexcept
2347{
2348 return this->minIndex<run_on>(this->domain,comp);
2349}
2350
2351template <class T>
2352template <RunOn run_on>
2353IntVect
2354BaseFab<T>::minIndex (const Box& subbox, int comp) const noexcept
2355{
2356 T min_val = this->min<run_on>(subbox, comp);
2357 return this->indexFromValue<run_on>(subbox, comp, min_val);
2358}
2359
2360template <class T>
2361template <RunOn run_on>
2362void
2363BaseFab<T>::minIndex (const Box& subbox, Real& min_val, IntVect& min_idx, int comp) const noexcept
2364{
2365 min_val = this->min<run_on>(subbox, comp);
2366 min_idx = this->indexFromValue<run_on>(subbox, comp, min_val);
2367}
2368
2369template <class T>
2370template <RunOn run_on>
2371IntVect
2372BaseFab<T>::maxIndex (int comp) const noexcept
2373{
2374 return this->maxIndex<run_on>(this->domain,comp);
2375}
2376
2377template <class T>
2378template <RunOn run_on>
2379IntVect
2380BaseFab<T>::maxIndex (const Box& subbox, int comp) const noexcept
2381{
2382 T max_val = this->max<run_on>(subbox, comp);
2383 return this->indexFromValue<run_on>(subbox, comp, max_val);
2384}
2385
2386template <class T>
2387template <RunOn run_on>
2388void
2389BaseFab<T>::maxIndex (const Box& subbox, Real& max_val, IntVect& max_idx, int comp) const noexcept
2390{
2391 max_val = this->max<run_on>(subbox, comp);
2392 max_idx = this->indexFromValue<run_on>(subbox, comp, max_val);
2393}
2394
2395template <class T>
2396template <RunOn run_on>
2397int
2398BaseFab<T>::maskLT (BaseFab<int>& mask, T const& val, int comp) const noexcept
2399{
2400 mask.resize(this->domain,1);
2401 int cnt = 0;
2402 Array4<int> const& m = mask.array();
2403 Array4<T const> const& a = this->const_array(comp);
2404#ifdef AMREX_USE_GPU
2405 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2406 ReduceOps<ReduceOpSum> reduce_op;
2407 ReduceData<int> reduce_data(reduce_op);
2408 using ReduceTuple = typename decltype(reduce_data)::Type;
2409 reduce_op.eval(this->domain, reduce_data,
2410 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2411 {
2412 int t;
2413 if (a(i,j,k) < val) {
2414 m(i,j,k) = 1;
2415 t = 1;
2416 } else {
2417 m(i,j,k) = 0;
2418 t = 0;
2419 }
2420 return {t};
2421 });
2422 ReduceTuple hv = reduce_data.value(reduce_op);
2423 cnt = amrex::get<0>(hv);
2424 } else
2425#endif
2426 {
2427 AMREX_LOOP_3D(this->domain, i, j, k,
2428 {
2429 if (a(i,j,k) < val) {
2430 m(i,j,k) = 1;
2431 ++cnt;
2432 } else {
2433 m(i,j,k) = 0;
2434 }
2435 });
2436 }
2437
2438 return cnt;
2439}
2440
2441template <class T>
2442template <RunOn run_on>
2443int
2444BaseFab<T>::maskLE (BaseFab<int>& mask, T const& val, int comp) const noexcept
2445{
2446 mask.resize(this->domain,1);
2447 int cnt = 0;
2448 Array4<int> const& m = mask.array();
2449 Array4<T const> const& a = this->const_array(comp);
2450#ifdef AMREX_USE_GPU
2451 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2452 ReduceOps<ReduceOpSum> reduce_op;
2453 ReduceData<int> reduce_data(reduce_op);
2454 using ReduceTuple = typename decltype(reduce_data)::Type;
2455 reduce_op.eval(this->domain, reduce_data,
2456 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2457 {
2458 int t;
2459 if (a(i,j,k) <= val) {
2460 m(i,j,k) = 1;
2461 t = 1;
2462 } else {
2463 m(i,j,k) = 0;
2464 t = 0;
2465 }
2466 return {t};
2467 });
2468 ReduceTuple hv = reduce_data.value(reduce_op);
2469 cnt = amrex::get<0>(hv);
2470 } else
2471#endif
2472 {
2473 AMREX_LOOP_3D(this->domain, i, j, k,
2474 {
2475 if (a(i,j,k) <= val) {
2476 m(i,j,k) = 1;
2477 ++cnt;
2478 } else {
2479 m(i,j,k) = 0;
2480 }
2481 });
2482 }
2483
2484 return cnt;
2485}
2486
2487template <class T>
2488template <RunOn run_on>
2489int
2490BaseFab<T>::maskEQ (BaseFab<int>& mask, T const& val, int comp) const noexcept
2491{
2492 mask.resize(this->domain,1);
2493 int cnt = 0;
2494 Array4<int> const& m = mask.array();
2495 Array4<T const> const& a = this->const_array(comp);
2496#ifdef AMREX_USE_GPU
2497 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2498 ReduceOps<ReduceOpSum> reduce_op;
2499 ReduceData<int> reduce_data(reduce_op);
2500 using ReduceTuple = typename decltype(reduce_data)::Type;
2501 reduce_op.eval(this->domain, reduce_data,
2502 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2503 {
2504 int t;
2505 if (a(i,j,k) == val) {
2506 m(i,j,k) = 1;
2507 t = 1;
2508 } else {
2509 m(i,j,k) = 0;
2510 t = 0;
2511 }
2512 return {t};
2513 });
2514 ReduceTuple hv = reduce_data.value(reduce_op);
2515 cnt = amrex::get<0>(hv);
2516 } else
2517#endif
2518 {
2519 AMREX_LOOP_3D(this->domain, i, j, k,
2520 {
2521 if (a(i,j,k) == val) {
2522 m(i,j,k) = 1;
2523 ++cnt;
2524 } else {
2525 m(i,j,k) = 0;
2526 }
2527 });
2528 }
2529
2530 return cnt;
2531}
2532
2533template <class T>
2534template <RunOn run_on>
2535int
2536BaseFab<T>::maskGT (BaseFab<int>& mask, T const& val, int comp) const noexcept
2537{
2538 mask.resize(this->domain,1);
2539 int cnt = 0;
2540 Array4<int> const& m = mask.array();
2541 Array4<T const> const& a = this->const_array(comp);
2542#ifdef AMREX_USE_GPU
2543 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2544 ReduceOps<ReduceOpSum> reduce_op;
2545 ReduceData<int> reduce_data(reduce_op);
2546 using ReduceTuple = typename decltype(reduce_data)::Type;
2547 reduce_op.eval(this->domain, reduce_data,
2548 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2549 {
2550 int t;
2551 if (a(i,j,k) > val) {
2552 m(i,j,k) = 1;
2553 t = 1;
2554 } else {
2555 m(i,j,k) = 0;
2556 t = 0;
2557 }
2558 return {t};
2559 });
2560 ReduceTuple hv = reduce_data.value(reduce_op);
2561 cnt = amrex::get<0>(hv);
2562 } else
2563#endif
2564 {
2565 AMREX_LOOP_3D(this->domain, i, j, k,
2566 {
2567 if (a(i,j,k) > val) {
2568 m(i,j,k) = 1;
2569 ++cnt;
2570 } else {
2571 m(i,j,k) = 0;
2572 }
2573 });
2574 }
2575
2576 return cnt;
2577}
2578
2579template <class T>
2580template <RunOn run_on>
2581int
2582BaseFab<T>::maskGE (BaseFab<int>& mask, T const& val, int comp) const noexcept
2583{
2584 mask.resize(this->domain,1);
2585 int cnt = 0;
2586 Array4<int> const& m = mask.array();
2587 Array4<T const> const& a = this->const_array(comp);
2588#ifdef AMREX_USE_GPU
2589 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2590 ReduceOps<ReduceOpSum> reduce_op;
2591 ReduceData<int> reduce_data(reduce_op);
2592 using ReduceTuple = typename decltype(reduce_data)::Type;
2593 reduce_op.eval(this->domain, reduce_data,
2594 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2595 {
2596 int t;
2597 if (a(i,j,k) >= val) {
2598 m(i,j,k) = 1;
2599 t = 1;
2600 } else {
2601 m(i,j,k) = 0;
2602 t = 0;
2603 }
2604 return {t};
2605 });
2606 ReduceTuple hv = reduce_data.value(reduce_op);
2607 cnt = amrex::get<0>(hv);
2608 } else
2609#endif
2610 {
2611 AMREX_LOOP_3D(this->domain, i, j, k,
2612 {
2613 if (a(i,j,k) >= val) {
2614 m(i,j,k) = 1;
2615 ++cnt;
2616 } else {
2617 m(i,j,k) = 0;
2618 }
2619 });
2620 }
2621
2622 return cnt;
2623}
2624
2625template <class T>
2626template <RunOn run_on>
2629{
2630 Box ovlp(this->domain);
2631 ovlp &= x.domain;
2632 return ovlp.ok() ? this->atomicAdd<run_on>(x,ovlp,ovlp,0,0,this->nvar) : *this;
2633}
2634
2635template <class T>
2636template <RunOn run_on>
2638BaseFab<T>::saxpy (T a, const BaseFab<T>& x, const Box& srcbox, const Box& destbox,
2639 int srccomp, int destcomp, int numcomp) noexcept
2640{
2641 BL_ASSERT(srcbox.ok());
2642 BL_ASSERT(x.box().contains(srcbox));
2643 BL_ASSERT(destbox.ok());
2644 BL_ASSERT(box().contains(destbox));
2645 BL_ASSERT(destbox.sameSize(srcbox));
2646 BL_ASSERT( srccomp >= 0 && srccomp+numcomp <= x.nComp());
2647 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2648
2649 Array4<T> const& d = this->array();
2650 Array4<T const> const& s = x.const_array();
2651 const auto dlo = amrex::lbound(destbox);
2652 const auto slo = amrex::lbound(srcbox);
2653 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
2654 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
2655 {
2656 d(i,j,k,n+destcomp) += a * s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
2657 });
2658
2659 return *this;
2660}
2661
2662template <class T>
2663template <RunOn run_on>
2665BaseFab<T>::saxpy (T a, const BaseFab<T>& x) noexcept
2666{
2667 Box ovlp(this->domain);
2668 ovlp &= x.domain;
2669 return ovlp.ok() ? saxpy<run_on>(a,x,ovlp,ovlp,0,0,this->nvar) : *this;
2670}
2671
2672template <class T>
2673template <RunOn run_on>
2676 const Box& srcbox, const Box& destbox,
2677 int srccomp, int destcomp, int numcomp) noexcept
2678{
2679 BL_ASSERT(srcbox.ok());
2680 BL_ASSERT(x.box().contains(srcbox));
2681 BL_ASSERT(destbox.ok());
2682 BL_ASSERT(box().contains(destbox));
2683 BL_ASSERT(destbox.sameSize(srcbox));
2684 BL_ASSERT( srccomp >= 0 && srccomp+numcomp <= x.nComp());
2685 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2686
2687 Array4<T> const& d = this->array();
2688 Array4<T const> const& s = x.const_array();
2689 const auto dlo = amrex::lbound(destbox);
2690 const auto slo = amrex::lbound(srcbox);
2691 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
2692 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
2693 {
2694 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);
2695 });
2696
2697 return *this;
2698}
2699
2700template <class T>
2701template <RunOn run_on>
2703BaseFab<T>::addproduct (const Box& destbox, int destcomp, int numcomp,
2704 const BaseFab<T>& src1, int comp1,
2705 const BaseFab<T>& src2, int comp2) noexcept
2706{
2707 BL_ASSERT(destbox.ok());
2708 BL_ASSERT(box().contains(destbox));
2709 BL_ASSERT( comp1 >= 0 && comp1+numcomp <= src1.nComp());
2710 BL_ASSERT( comp2 >= 0 && comp2+numcomp <= src2.nComp());
2711 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2712
2713 Array4<T> const& d = this->array();
2714 Array4<T const> const& s1 = src1.const_array();
2715 Array4<T const> const& s2 = src2.const_array();
2716 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
2717 {
2718 d(i,j,k,n+destcomp) += s1(i,j,k,n+comp1) * s2(i,j,k,n+comp2);
2719 });
2720
2721 return *this;
2722}
2723
2724template <class T>
2725template <RunOn run_on>
2727BaseFab<T>::linComb (const BaseFab<T>& f1, const Box& b1, int comp1,
2728 const BaseFab<T>& f2, const Box& b2, int comp2,
2729 Real alpha, Real beta, const Box& b,
2730 int comp, int numcomp) noexcept
2731{
2732 BL_ASSERT(b1.ok());
2733 BL_ASSERT(f1.box().contains(b1));
2734 BL_ASSERT(b2.ok());
2735 BL_ASSERT(f2.box().contains(b2));
2736 BL_ASSERT(b.ok());
2737 BL_ASSERT(box().contains(b));
2738 BL_ASSERT(b.sameSize(b1));
2739 BL_ASSERT(b.sameSize(b2));
2740 BL_ASSERT(comp1 >= 0 && comp1+numcomp <= f1.nComp());
2741 BL_ASSERT(comp2 >= 0 && comp2+numcomp <= f2.nComp());
2742 BL_ASSERT(comp >= 0 && comp +numcomp <= nComp());
2743
2744 Array4<T> const& d = this->array();
2745 Array4<T const> const& s1 = f1.const_array();
2746 Array4<T const> const& s2 = f2.const_array();
2747 const auto dlo = amrex::lbound(b);
2748 const auto slo1 = amrex::lbound(b1);
2749 const auto slo2 = amrex::lbound(b2);
2750 const Dim3 off1{.x = slo1.x-dlo.x, .y = slo1.y-dlo.y, .z = slo1.z-dlo.z};
2751 const Dim3 off2{.x = slo2.x-dlo.x, .y = slo2.y-dlo.y, .z = slo2.z-dlo.z};
2752
2753 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, b, numcomp, i, j, k, n,
2754 {
2755 d(i,j,k,n+comp) = alpha*s1(i+off1.x,j+off1.y,k+off1.z,n+comp1)
2756 + beta*s2(i+off2.x,j+off2.y,k+off2.z,n+comp2);
2757 });
2758 return *this;
2759}
2760
2761template <class T>
2762template <RunOn run_on>
2763T
2764BaseFab<T>::dot (const Box& xbx, int xcomp,
2765 const BaseFab<T>& y, const Box& ybx, int ycomp,
2766 int numcomp) const noexcept
2767{
2768 BL_ASSERT(xbx.ok());
2769 BL_ASSERT(box().contains(xbx));
2770 BL_ASSERT(y.box().contains(ybx));
2771 BL_ASSERT(xbx.sameSize(ybx));
2772 BL_ASSERT(xcomp >= 0 && xcomp+numcomp <= nComp());
2773 BL_ASSERT(ycomp >= 0 && ycomp+numcomp <= y.nComp());
2774
2775 T r = 0;
2776
2777 const auto xlo = amrex::lbound(xbx);
2778 const auto ylo = amrex::lbound(ybx);
2779 const Dim3 offset{.x = ylo.x-xlo.x, .y = ylo.y-xlo.y, .z = ylo.z-xlo.z};
2780 Array4<T const> const& xa = this->const_array();
2781 Array4<T const> const& ya = y.const_array();
2782
2783#ifdef AMREX_USE_GPU
2784 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2785 ReduceOps<ReduceOpSum> reduce_op;
2786 ReduceData<T> reduce_data(reduce_op);
2787 using ReduceTuple = typename decltype(reduce_data)::Type;
2788 reduce_op.eval(xbx, reduce_data,
2789 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2790 {
2791 T t = 0;
2792 for (int n = 0; n < numcomp; ++n) {
2793 t += xa(i,j,k,n+xcomp) * ya(i+offset.x,j+offset.y,k+offset.z,n+ycomp);
2794 }
2795 return {t};
2796 });
2797 ReduceTuple hv = reduce_data.value(reduce_op);
2798 r = amrex::get<0>(hv);
2799 } else
2800#endif
2801 {
2802 AMREX_LOOP_4D(xbx, numcomp, i, j, k, n,
2803 {
2804 r += xa(i,j,k,n+xcomp) * ya(i+offset.x,j+offset.y,k+offset.z,n+ycomp);
2805 });
2806 }
2807
2808 return r;
2809}
2810
2811template <class T>
2812template <RunOn run_on>
2813T
2814BaseFab<T>::dotmask (const BaseFab<int>& mask, const Box& xbx, int xcomp,
2815 const BaseFab<T>& y, const Box& ybx, int ycomp,
2816 int numcomp) const noexcept
2817{
2818 BL_ASSERT(xbx.ok());
2819 BL_ASSERT(box().contains(xbx));
2820 BL_ASSERT(y.box().contains(ybx));
2821 BL_ASSERT(xbx.sameSize(ybx));
2822 BL_ASSERT(xcomp >= 0 && xcomp+numcomp <= nComp());
2823 BL_ASSERT(ycomp >= 0 && ycomp+numcomp <= y.nComp());
2824
2825 T r = 0;
2826
2827 const auto xlo = amrex::lbound(xbx);
2828 const auto ylo = amrex::lbound(ybx);
2829 const Dim3 offset{.x = ylo.x-xlo.x, .y = ylo.y-xlo.y, .z = ylo.z-xlo.z};
2830
2831 Array4<T const> const& xa = this->const_array();
2832 Array4<T const> const& ya = y.const_array();
2833 Array4<int const> const& ma = mask.const_array();
2834
2835#ifdef AMREX_USE_GPU
2836 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
2837 ReduceOps<ReduceOpSum> reduce_op;
2838 ReduceData<T> reduce_data(reduce_op);
2839 using ReduceTuple = typename decltype(reduce_data)::Type;
2840 reduce_op.eval(xbx, reduce_data,
2841 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
2842 {
2843 int m = static_cast<int>(static_cast<bool>(ma(i,j,k)));
2844 T t = 0;
2845 for (int n = 0; n < numcomp; ++n) {
2846 t += xa(i,j,k,n+xcomp) * ya(i+offset.x,j+offset.y,k+offset.z,n+ycomp) * m;
2847 }
2848 return {t};
2849 });
2850 ReduceTuple hv = reduce_data.value(reduce_op);
2851 r = amrex::get<0>(hv);
2852 } else
2853#endif
2854 {
2855 AMREX_LOOP_4D(xbx, numcomp, i, j, k, n,
2856 {
2857 int m = static_cast<int>(static_cast<bool>(ma(i,j,k)));
2858 r += xa(i,j,k,n+xcomp) * ya(i+offset.x,j+offset.y,k+offset.z,n+ycomp) * m;
2859 });
2860 }
2861
2862 return r;
2863}
2864
2865template <class T>
2866template <RunOn run_on>
2867T
2868BaseFab<T>::sum (int comp, int numcomp) const noexcept
2869{
2870 return this->sum<run_on>(this->domain, DestComp{comp}, NumComps{numcomp});
2871}
2872
2873template <class T>
2874template <RunOn run_on>
2875T
2876BaseFab<T>::sum (const Box& subbox, int comp, int numcomp) const noexcept
2877{
2878 return this->sum<run_on>(subbox, DestComp{comp}, NumComps{numcomp});
2879}
2880
2881template <class T>
2882template <RunOn run_on>
2884BaseFab<T>::negate (int comp, int numcomp) noexcept
2885{
2886 return this->negate<run_on>(this->domain, DestComp{comp}, NumComps{numcomp});
2887}
2888
2889template <class T>
2890template <RunOn run_on>
2892BaseFab<T>::negate (const Box& b, int comp, int numcomp) noexcept
2893{
2894 return this->negate<run_on>(b, DestComp{comp}, NumComps{numcomp});
2895}
2896
2897template <class T>
2898template <RunOn run_on>
2900BaseFab<T>::invert (T const& r, int comp, int numcomp) noexcept
2901{
2902 return this->invert<run_on>(r, this->domain, DestComp{comp}, NumComps{numcomp});
2903}
2904
2905template <class T>
2906template <RunOn run_on>
2908BaseFab<T>::invert (T const& r, const Box& b, int comp, int numcomp) noexcept
2909{
2910 return this->invert<run_on>(r, b, DestComp{comp}, NumComps{numcomp});
2911}
2912
2913template <class T>
2914template <RunOn run_on>
2916BaseFab<T>::plus (T const& r, int comp, int numcomp) noexcept
2917{
2918 return this->plus<run_on>(r, this->domain, DestComp{comp}, NumComps{numcomp});
2919}
2920
2921template <class T>
2922template <RunOn run_on>
2924BaseFab<T>::plus (T const& r, const Box& b, int comp, int numcomp) noexcept
2925{
2926 return this->plus<run_on>(r, b, DestComp{comp}, NumComps{numcomp});
2927}
2928
2929template <class T>
2930template <RunOn run_on>
2932BaseFab<T>::plus (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
2933{
2934 return this->plus<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
2935}
2936
2937template <class T>
2938template <RunOn run_on>
2940BaseFab<T>::atomicAdd (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
2941{
2942 Box ovlp(this->domain);
2943 ovlp &= src.domain;
2944 return ovlp.ok() ? this->atomicAdd<run_on>(src,ovlp,ovlp,srccomp,destcomp,numcomp) : *this;
2945}
2946
2947template <class T>
2948template <RunOn run_on>
2950BaseFab<T>::plus (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
2951 int numcomp) noexcept
2952{
2953 return this->plus<run_on>(src, subbox, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
2954}
2955
2956template <class T>
2957template <RunOn run_on>
2959BaseFab<T>::atomicAdd (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
2960 int numcomp) noexcept
2961{
2962 Box ovlp(this->domain);
2963 ovlp &= src.domain;
2964 ovlp &= subbox;
2965 return ovlp.ok() ? this->atomicAdd<run_on>(src,ovlp,ovlp,srccomp,destcomp,numcomp) : *this;
2966}
2967
2968template <class T>
2969template <RunOn run_on>
2971BaseFab<T>::plus (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
2972 int srccomp, int destcomp, int numcomp) noexcept
2973{
2974 BL_ASSERT(destbox.ok());
2975 BL_ASSERT(src.box().contains(srcbox));
2976 BL_ASSERT(box().contains(destbox));
2977 BL_ASSERT(destbox.sameSize(srcbox));
2978 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
2979 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
2980
2981 Array4<T> const& d = this->array();
2982 Array4<T const> const& s = src.const_array();
2983 const auto dlo = amrex::lbound(destbox);
2984 const auto slo = amrex::lbound(srcbox);
2985 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
2986 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
2987 {
2988 d(i,j,k,n+destcomp) += s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
2989 });
2990
2991 return *this;
2992}
2993
2994namespace detail
2995{
2996
2997template <RunOn run_on, typename T>
2998requires (HasAtomicAdd<T>::value)
2999void basefab_atomic_add (BaseFab<T>& dfab, const BaseFab<T>& sfab,
3000 const Box& srcbox, const Box& destbox,
3001 int srccomp, int destcomp, int numcomp) noexcept
3002{
3003 Array4<T> const& d = dfab.array();
3004 Array4<T const> const& s = sfab.const_array();
3005 const auto dlo = amrex::lbound(destbox);
3006 const auto slo = amrex::lbound(srcbox);
3007 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3008 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3009 {
3010 T* p = d.ptr(i,j,k,n+destcomp);
3011 HostDevice::Atomic::Add(p, s(i+offset.x,j+offset.y,k+offset.z,n+srccomp));
3012 });
3013}
3014
3015template <RunOn run_on, typename T>
3016requires (!HasAtomicAdd<T>::value)
3017void basefab_atomic_add (BaseFab<T>& dfab, const BaseFab<T>& sfab,
3018 const Box& srcbox, const Box& destbox,
3019 int srccomp, int destcomp, int numcomp)
3020{
3021 amrex::ignore_unused(dfab, sfab, srcbox, destbox, srccomp, destcomp, numcomp);
3022 amrex::Abort("BaseFab: atomicAdd not supported");
3023}
3024
3025}
3026
3027template <class T>
3028template <RunOn run_on>
3029BaseFab<T>&
3030BaseFab<T>::atomicAdd (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3031 int srccomp, int destcomp, int numcomp) noexcept
3032{
3033 BL_ASSERT(destbox.ok());
3034 BL_ASSERT(src.box().contains(srcbox));
3035 BL_ASSERT(box().contains(destbox));
3036 BL_ASSERT(destbox.sameSize(srcbox));
3037 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3038 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3039
3040 detail::basefab_atomic_add<run_on>(*this, src, srcbox, destbox,
3041 srccomp, destcomp, numcomp);
3042
3043 return *this;
3044}
3045
3046template <class T>
3047template <RunOn run_on>
3049BaseFab<T>::lockAdd (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3050 int srccomp, int destcomp, int numcomp) noexcept
3051{
3052#if defined(AMREX_USE_OMP) && (AMREX_SPACEDIM > 1)
3053#if defined(AMREX_USE_GPU)
3054 if (run_on == RunOn::Host || Gpu::notInLaunchRegion()) {
3055#endif
3056 BL_ASSERT(destbox.ok());
3057 BL_ASSERT(src.box().contains(srcbox));
3058 BL_ASSERT(box().contains(destbox));
3059 BL_ASSERT(destbox.sameSize(srcbox));
3060 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3061 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3062
3063 Array4<T> const& d = this->array();
3064 Array4<T const> const& s = src.const_array();
3065 auto const& dlo = amrex::lbound(destbox);
3066 auto const& dhi = amrex::ubound(destbox);
3067 auto const& len = amrex::length(destbox);
3068 auto const& slo = amrex::lbound(srcbox);
3069 Dim3 const offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3070
3071 // The lock dimension must be a property of the destination fab, not of
3072 // destbox, so that concurrent calls into the same fab use the same locks.
3073 int planedim;
3074 int nplanes;
3075 int plo;
3076 if (amrex::length(this->domain).z == 1) {
3077 planedim = 1;
3078 nplanes = len.y;
3079 plo = dlo.y;
3080 } else {
3081 planedim = 2;
3082 nplanes = len.z;
3083 plo = dlo.z;
3084 }
3085
3086 auto* mask = (bool*) amrex_mempool_alloc(sizeof(bool)*nplanes);
3087 for (int ip = 0; ip < nplanes; ++ip) {
3088 mask[ip] = false;
3089 }
3090
3091 int mm = 0;
3092 int planes_left = nplanes;
3093 while (planes_left > 0) {
3094 AMREX_ASSERT(mm < nplanes);
3095 auto const m = mm + plo;
3096 auto* lock = OpenMP::get_lock(m);
3097 if (omp_test_lock(lock))
3098 {
3099 auto lo = dlo;
3100 auto hi = dhi;
3101 if (planedim == 1) {
3102 lo.y = m;
3103 hi.y = m;
3104 } else {
3105 lo.z = m;
3106 hi.z = m;
3107 }
3108
3109 for (int n = 0; n < numcomp; ++n) {
3110 for (int k = lo.z; k <= hi.z; ++k) {
3111 for (int j = lo.y; j <= hi.y; ++j) {
3112 auto * pdst = d.ptr(dlo.x,j ,k ,n+destcomp);
3113 auto const* psrc = s.ptr(slo.x,j+offset.y,k+offset.z,n+ srccomp);
3114#pragma omp simd
3115 for (int ii = 0; ii < len.x; ++ii) {
3116 pdst[ii] += psrc[ii];
3117 }
3118 }
3119 }
3120 }
3121
3122 mask[mm] = true;
3123 --planes_left;
3124 omp_unset_lock(lock);
3125 if (planes_left == 0) { break; }
3126 }
3127
3128 ++mm;
3129 for (int ip = 0; ip < nplanes; ++ip) {
3130 int new_mm = (mm+ip) % nplanes;
3131 if ( ! mask[new_mm] ) {
3132 mm = new_mm;
3133 break;
3134 }
3135 }
3136 }
3137
3139
3140 return *this;
3141
3142#if defined(AMREX_USE_GPU)
3143 } else {
3144 return this->template atomicAdd<run_on>(src, srcbox, destbox, srccomp, destcomp, numcomp);
3145 }
3146#endif
3147#else
3148 return this->template atomicAdd<run_on>(src, srcbox, destbox, srccomp, destcomp, numcomp);
3149#endif
3150}
3151
3152template <class T>
3153template <RunOn run_on>
3155BaseFab<T>::minus (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
3156{
3157 return this->minus<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3158}
3159
3160template <class T>
3161template <RunOn run_on>
3163BaseFab<T>::minus (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp, int numcomp) noexcept
3164{
3165 return this->minus<run_on>(src, subbox, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3166}
3167
3168template <class T>
3169template <RunOn run_on>
3171BaseFab<T>::minus (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3172 int srccomp, int destcomp, int numcomp) noexcept
3173{
3174 BL_ASSERT(destbox.ok());
3175 BL_ASSERT(src.box().contains(srcbox));
3176 BL_ASSERT(box().contains(destbox));
3177 BL_ASSERT(destbox.sameSize(srcbox));
3178 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3179 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3180
3181 Array4<T> const& d = this->array();
3182 Array4<T const> const& s = src.const_array();
3183 const auto dlo = amrex::lbound(destbox);
3184 const auto slo = amrex::lbound(srcbox);
3185 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3186 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3187 {
3188 d(i,j,k,n+destcomp) -= s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
3189 });
3190
3191 return *this;
3192}
3193
3194template <class T>
3195template <RunOn run_on>
3197BaseFab<T>::mult (T const& r, int comp, int numcomp) noexcept
3198{
3199 return this->mult<run_on>(r, this->domain, DestComp{comp}, NumComps{numcomp});
3200}
3201
3202template <class T>
3203template <RunOn run_on>
3205BaseFab<T>::mult (T const& r, const Box& b, int comp, int numcomp) noexcept
3206{
3207 return this->mult<run_on>(r, b, DestComp{comp}, NumComps{numcomp});
3208}
3209
3210template <class T>
3211template <RunOn run_on>
3213BaseFab<T>::mult (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
3214{
3215 return this->mult<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3216}
3217
3218template <class T>
3219template <RunOn run_on>
3221BaseFab<T>::mult (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp, int numcomp) noexcept
3222{
3223 return this->mult<run_on>(src, subbox, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3224}
3225
3226template <class T>
3227template <RunOn run_on>
3229BaseFab<T>::mult (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3230 int srccomp, int destcomp, int numcomp) noexcept
3231{
3232 BL_ASSERT(destbox.ok());
3233 BL_ASSERT(src.box().contains(srcbox));
3234 BL_ASSERT(box().contains(destbox));
3235 BL_ASSERT(destbox.sameSize(srcbox));
3236 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3237 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3238
3239 Array4<T> const& d = this->array();
3240 Array4<T const> const& s = src.const_array();
3241 const auto dlo = amrex::lbound(destbox);
3242 const auto slo = amrex::lbound(srcbox);
3243 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3244 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3245 {
3246 d(i,j,k,n+destcomp) *= s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
3247 });
3248
3249 return *this;
3250}
3251
3252template <class T>
3253template <RunOn run_on>
3255BaseFab<T>::divide (T const& r, int comp, int numcomp) noexcept
3256{
3257 return this->divide<run_on>(r, this->domain, DestComp{comp}, NumComps{numcomp});
3258}
3259
3260template <class T>
3261template <RunOn run_on>
3263BaseFab<T>::divide (T const& r, const Box& b, int comp, int numcomp) noexcept
3264{
3265 return this->divide<run_on>(r, b, DestComp{comp}, NumComps{numcomp});
3266}
3267
3268template <class T>
3269template <RunOn run_on>
3271BaseFab<T>::divide (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
3272{
3273 return this->divide<run_on>(src, this->domain, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3274}
3275
3276template <class T>
3277template <RunOn run_on>
3279BaseFab<T>::divide (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp, int numcomp) noexcept
3280{
3281 return this->divide<run_on>(src, subbox, SrcComp{srccomp}, DestComp{destcomp}, NumComps{numcomp});
3282}
3283
3284template <class T>
3285template <RunOn run_on>
3287BaseFab<T>::divide (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3288 int srccomp, int destcomp, int numcomp) noexcept
3289{
3290 BL_ASSERT(destbox.ok());
3291 BL_ASSERT(src.box().contains(srcbox));
3292 BL_ASSERT(box().contains(destbox));
3293 BL_ASSERT(destbox.sameSize(srcbox));
3294 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3295 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3296
3297 Array4<T> const& d = this->array();
3298 Array4<T const> const& s = src.const_array();
3299 const auto dlo = amrex::lbound(destbox);
3300 const auto slo = amrex::lbound(srcbox);
3301 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3302 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3303 {
3304 d(i,j,k,n+destcomp) /= s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
3305 });
3306
3307 return *this;
3308}
3309
3310template <class T>
3311template <RunOn run_on>
3314{
3315 Box ovlp(this->domain);
3316 ovlp &= src.domain;
3317 return ovlp.ok() ? this->protected_divide<run_on>(src,ovlp,ovlp,0,0,this->nvar) : *this;
3318}
3319
3320template <class T>
3321template <RunOn run_on>
3323BaseFab<T>::protected_divide (const BaseFab<T>& src, int srccomp, int destcomp, int numcomp) noexcept
3324{
3325 Box ovlp(this->domain);
3326 ovlp &= src.domain;
3327 return ovlp.ok() ? this->protected_divide<run_on>(src,ovlp,ovlp,srccomp,destcomp,numcomp) : *this;
3328}
3329
3330template <class T>
3331template <RunOn run_on>
3333BaseFab<T>::protected_divide (const BaseFab<T>& src, const Box& subbox, int srccomp, int destcomp,
3334 int numcomp) noexcept
3335{
3336 Box ovlp(this->domain);
3337 ovlp &= src.domain;
3338 ovlp &= subbox;
3339 return ovlp.ok() ? this->protected_divide<run_on>(src,ovlp,ovlp,srccomp,destcomp,numcomp) : *this;
3340}
3341
3342template <class T>
3343template <RunOn run_on>
3345BaseFab<T>::protected_divide (const BaseFab<T>& src, const Box& srcbox, const Box& destbox,
3346 int srccomp, int destcomp, int numcomp) noexcept
3347{
3348 BL_ASSERT(destbox.ok());
3349 BL_ASSERT(src.box().contains(srcbox));
3350 BL_ASSERT(box().contains(destbox));
3351 BL_ASSERT(destbox.sameSize(srcbox));
3352 BL_ASSERT(srccomp >= 0 && srccomp+numcomp <= src.nComp());
3353 BL_ASSERT(destcomp >= 0 && destcomp+numcomp <= nComp());
3354
3355 Array4<T> const& d = this->array();
3356 Array4<T const> const& s = src.const_array();
3357 const auto dlo = amrex::lbound(destbox);
3358 const auto slo = amrex::lbound(srcbox);
3359 const Dim3 offset{.x = slo.x-dlo.x, .y = slo.y-dlo.y, .z = slo.z-dlo.z};
3360 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, destbox, numcomp, i, j, k, n,
3361 {
3362 if (s(i+offset.x,j+offset.y,k+offset.z,n+srccomp) != 0) {
3363 d(i,j,k,n+destcomp) /= s(i+offset.x,j+offset.y,k+offset.z,n+srccomp);
3364 }
3365 });
3366
3367 return *this;
3368}
3369
3380template <class T>
3381template <RunOn run_on>
3383BaseFab<T>::linInterp (const BaseFab<T>& f1, const Box& b1, int comp1,
3384 const BaseFab<T>& f2, const Box& b2, int comp2,
3385 Real t1, Real t2, Real t,
3386 const Box& b, int comp, int numcomp) noexcept
3387{
3388 if (amrex::almostEqual(t1,t2) || amrex::almostEqual(t1,t)) {
3389 return copy<run_on>(f1,b1,comp1,b,comp,numcomp);
3390 } else if (amrex::almostEqual(t2,t)) {
3391 return copy<run_on>(f2,b2,comp2,b,comp,numcomp);
3392 } else {
3393 Real alpha = (t2-t)/(t2-t1);
3394 Real beta = (t-t1)/(t2-t1);
3395 return linComb<run_on>(f1,b1,comp1,f2,b2,comp2,alpha,beta,b,comp,numcomp);
3396 }
3397}
3398
3399template <class T>
3400template <RunOn run_on>
3402BaseFab<T>::linInterp (const BaseFab<T>& f1, int comp1,
3403 const BaseFab<T>& f2, int comp2,
3404 Real t1, Real t2, Real t,
3405 const Box& b, int comp, int numcomp) noexcept
3406{
3407 if (amrex::almostEqual(t1,t2) || amrex::almostEqual(t1,t)) {
3408 return copy<run_on>(f1,b,comp1,b,comp,numcomp);
3409 } else if (amrex::almostEqual(t2,t)) {
3410 return copy<run_on>(f2,b,comp2,b,comp,numcomp);
3411 } else {
3412 Real alpha = (t2-t)/(t2-t1);
3413 Real beta = (t-t1)/(t2-t1);
3414 return linComb<run_on>(f1,b,comp1,f2,b,comp2,alpha,beta,b,comp,numcomp);
3415 }
3416}
3417
3418//
3419// New interfaces
3420//
3421
3422template <class T>
3423template <RunOn run_on>
3424void
3425BaseFab<T>::setVal (T const& val) noexcept
3426{
3427 this->setVal<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3428}
3429
3430template <class T>
3431template <RunOn run_on>
3432void
3433BaseFab<T>::setVal (T const& x, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept
3434{
3435 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3436 Array4<T> const& a = this->array();
3437 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG (run_on, bx, ncomp.n, i, j, k, n,
3438 {
3439 a(i,j,k,n+dcomp.i) = x;
3440 });
3441}
3442
3443template <class T>
3444template <RunOn run_on>
3445void
3446BaseFab<T>::setValIf (T const& val, const BaseFab<int>& mask) noexcept
3447{
3448 this->setValIf<run_on>(val, this->domain, mask, DestComp{0}, NumComps{this->nvar});
3449}
3450
3451template <class T>
3452template <RunOn run_on>
3453void
3454BaseFab<T>::setValIf (T const& val, Box const& bx, const BaseFab<int>& mask, DestComp dcomp, NumComps ncomp) noexcept
3455{
3456 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3457 Array4<T> const& a = this->array();
3458 Array4<int const> const& m = mask.const_array();
3459 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG (run_on, bx, ncomp.n, i, j, k, n,
3460 {
3461 if (m(i,j,k)) { a(i,j,k,n+dcomp.i) = val; }
3462 });
3463}
3464
3465template <class T>
3466template <RunOn run_on>
3467void
3468BaseFab<T>::setValIfNot (T const& val, const BaseFab<int>& mask) noexcept
3469{
3470 this->setValIfNot<run_on>(val, this->domain, mask, DestComp{0}, NumComps{this->nvar});
3471}
3472
3473template <class T>
3474template <RunOn run_on>
3475void
3476BaseFab<T>::setValIfNot (T const& val, Box const& bx, const BaseFab<int>& mask, DestComp dcomp, NumComps ncomp) noexcept
3477{
3478 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3479 Array4<T> const& a = this->array();
3480 Array4<int const> const& m = mask.const_array();
3481 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG (run_on, bx, ncomp.n, i, j, k, n,
3482 {
3483 if (!m(i,j,k)) { a(i,j,k,n+dcomp.i) = val; }
3484 });
3485}
3486
3487template <class T>
3488template <RunOn run_on>
3489void
3490BaseFab<T>::setComplement (T const& x, const Box& bx, DestComp dcomp, NumComps ncomp) noexcept
3491{
3492#ifdef AMREX_USE_GPU
3493 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
3494 Array4<T> const& a = this->array();
3495 amrex::ParallelFor(this->domain, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
3496 {
3497 if (! bx.contains(IntVect(AMREX_D_DECL(i,j,k)))) {
3498 for (int n = dcomp.i; n < ncomp.n+dcomp.i; ++n) {
3499 a(i,j,k,n) = x;
3500 }
3501 }
3502 });
3503 } else
3504#endif
3505 {
3506 const BoxList b_lst = amrex::boxDiff(this->domain,bx);
3507 for (auto const& b : b_lst) {
3508 this->setVal<RunOn::Host>(x, b, dcomp, ncomp);
3509 }
3510 }
3511}
3512
3513template <class T>
3514template <RunOn run_on>
3516BaseFab<T>::copy (const BaseFab<T>& src) noexcept
3517{
3518 this->copy<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3519 return *this;
3520}
3521
3522template <class T>
3523template <RunOn run_on>
3526 SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
3527{
3528 AMREX_ASSERT(this->domain.sameType(src.domain));
3529 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3530 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3531
3532 bx &= src.domain;
3533
3534 Array4<T> const& d = this->array();
3535 Array4<T const> const& s = src.const_array();
3536 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3537 {
3538 d(i,j,k,n+dcomp.i) = s(i,j,k,n+scomp.i);
3539 });
3540
3541 return *this;
3542}
3543
3544template <class T>
3545template <RunOn run_on>
3547BaseFab<T>::plus (T const& val) noexcept
3548{
3549 return this->plus<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3550}
3551
3552template <class T>
3553template <RunOn run_on>
3555BaseFab<T>::operator+= (T const& val) noexcept
3556{
3557 return this->plus<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3558}
3559
3560template <class T>
3561template <RunOn run_on>
3563BaseFab<T>::plus (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept
3564{
3565 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3566
3567 Array4<T> const& a = this->array();
3568 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3569 {
3570 a(i,j,k,n+dcomp.i) += val;
3571 });
3572
3573 return *this;
3574}
3575
3576template <class T>
3577template <RunOn run_on>
3579BaseFab<T>::plus (const BaseFab<T>& src) noexcept
3580{
3581 return this->plus<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3582}
3583
3584template <class T>
3585template <RunOn run_on>
3588{
3589 return this->plus<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3590}
3591
3592template <class T>
3593template <RunOn run_on>
3596 SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
3597{
3598 AMREX_ASSERT(this->domain.sameType(src.domain));
3599 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3600 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3601
3602 bx &= src.domain;
3603
3604 Array4<T> const& d = this->array();
3605 Array4<T const> const& s = src.const_array();
3606 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3607 {
3608 d(i,j,k,n+dcomp.i) += s(i,j,k,n+scomp.i);
3609 });
3610
3611 return *this;
3612}
3613
3614template <class T>
3615template <RunOn run_on>
3617BaseFab<T>::minus (T const& val) noexcept
3618{
3619 return this->minus<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3620}
3621
3622template <class T>
3623template <RunOn run_on>
3625BaseFab<T>::operator-= (T const& val) noexcept
3626{
3627 return this->minus<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3628}
3629
3630template <class T>
3631template <RunOn run_on>
3633BaseFab<T>::minus (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept
3634{
3635 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3636
3637 Array4<T> const& a = this->array();
3638 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3639 {
3640 a(i,j,k,n+dcomp.i) -= val;
3641 });
3642
3643 return *this;
3644}
3645
3646template <class T>
3647template <RunOn run_on>
3649BaseFab<T>::minus (const BaseFab<T>& src) noexcept
3650{
3651 return this->minus<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3652}
3653
3654template <class T>
3655template <RunOn run_on>
3658{
3659 return this->minus<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3660}
3661
3662template <class T>
3663template <RunOn run_on>
3666 SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
3667{
3668 AMREX_ASSERT(this->domain.sameType(src.domain));
3669 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3670 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3671
3672 bx &= src.domain;
3673
3674 Array4<T> const& d = this->array();
3675 Array4<T const> const& s = src.const_array();
3676 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3677 {
3678 d(i,j,k,n+dcomp.i) -= s(i,j,k,n+scomp.i);
3679 });
3680
3681 return *this;
3682}
3683
3684template <class T>
3685template <RunOn run_on>
3687BaseFab<T>::mult (T const& val) noexcept
3688{
3689 return this->mult<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3690}
3691
3692template <class T>
3693template <RunOn run_on>
3695BaseFab<T>::operator*= (T const& val) noexcept
3696{
3697 return this->mult<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3698}
3699
3700template <class T>
3701template <RunOn run_on>
3703BaseFab<T>::mult (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept
3704{
3705 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3706
3707 Array4<T> const& a = this->array();
3708 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3709 {
3710 a(i,j,k,n+dcomp.i) *= val;
3711 });
3712
3713 return *this;
3714}
3715
3716template <class T>
3717template <RunOn run_on>
3719BaseFab<T>::mult (const BaseFab<T>& src) noexcept
3720{
3721 return this->mult<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3722}
3723
3724template <class T>
3725template <RunOn run_on>
3728{
3729 return this->mult<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3730}
3731
3732template <class T>
3733template <RunOn run_on>
3736 SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
3737{
3738 AMREX_ASSERT(this->domain.sameType(src.domain));
3739 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3740 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3741
3742 bx &= src.domain;
3743
3744 Array4<T> const& d = this->array();
3745 Array4<T const> const& s = src.const_array();
3746 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3747 {
3748 d(i,j,k,n+dcomp.i) *= s(i,j,k,n+scomp.i);
3749 });
3750
3751 return *this;
3752}
3753
3754template <class T>
3755template <RunOn run_on>
3757BaseFab<T>::divide (T const& val) noexcept
3758{
3759 return this->divide<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3760}
3761
3762template <class T>
3763template <RunOn run_on>
3765BaseFab<T>::operator/= (T const& val) noexcept
3766{
3767 return this->divide<run_on>(val, this->domain, DestComp{0}, NumComps{this->nvar});
3768}
3769
3770template <class T>
3771template <RunOn run_on>
3773BaseFab<T>::divide (T const& val, Box const& bx, DestComp dcomp, NumComps ncomp) noexcept
3774{
3775 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3776
3777 Array4<T> const& a = this->array();
3778 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3779 {
3780 a(i,j,k,n+dcomp.i) /= val;
3781 });
3782
3783 return *this;
3784}
3785
3786template <class T>
3787template <RunOn run_on>
3789BaseFab<T>::divide (const BaseFab<T>& src) noexcept
3790{
3791 return this->divide<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3792}
3793
3794template <class T>
3795template <RunOn run_on>
3798{
3799 return this->divide<run_on>(src, this->domain, SrcComp{0}, DestComp{0}, NumComps{this->nvar});
3800}
3801
3802template <class T>
3803template <RunOn run_on>
3806 SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
3807{
3808 AMREX_ASSERT(this->domain.sameType(src.domain));
3809 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3810 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3811
3812 bx &= src.domain;
3813
3814 Array4<T> const& d = this->array();
3815 Array4<T const> const& s = src.const_array();
3816 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3817 {
3818 d(i,j,k,n+dcomp.i) /= s(i,j,k,n+scomp.i);
3819 });
3820
3821 return *this;
3822}
3823
3824template <class T>
3825template <RunOn run_on>
3828{
3829 return this->negate<run_on>(this->domain, DestComp{0}, NumComps{this->nvar});
3830}
3831
3832template <class T>
3833template <RunOn run_on>
3835BaseFab<T>::negate (const Box& bx, DestComp dcomp, NumComps ncomp) noexcept
3836{
3837 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3838
3839 Array4<T> const& a = this->array();
3840 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3841 {
3842 a(i,j,k,n+dcomp.i) = -a(i,j,k,n+dcomp.i);
3843 });
3844
3845 return *this;
3846}
3847
3848template <class T>
3849template <RunOn run_on>
3851BaseFab<T>::invert (T const& r) noexcept
3852{
3853 return this->invert<run_on>(r, this->domain, DestComp{0}, NumComps{this->nvar});
3854}
3855
3856template <class T>
3857template <RunOn run_on>
3859BaseFab<T>::invert (T const& r, const Box& bx, DestComp dcomp, NumComps ncomp) noexcept
3860{
3861 BL_ASSERT(dcomp.i >= 0 && dcomp.i + ncomp.n <= this->nvar);
3862
3863 Array4<T> const& a = this->array();
3864 AMREX_HOST_DEVICE_PARALLEL_FOR_4D_FLAG(run_on, bx, ncomp.n, i, j, k, n,
3865 {
3866 a(i,j,k,n+dcomp.i) = r / a(i,j,k,n+dcomp.i);
3867 });
3868
3869 return *this;
3870}
3871
3872template <class T>
3873template <RunOn run_on>
3874T
3875BaseFab<T>::sum (const Box& bx, DestComp dcomp, NumComps ncomp) const noexcept
3876{
3877 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3878
3879 T r = 0;
3880 Array4<T const> const& a = this->const_array();
3881#ifdef AMREX_USE_GPU
3882 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
3883 ReduceOps<ReduceOpSum> reduce_op;
3884 ReduceData<T> reduce_data(reduce_op);
3885 using ReduceTuple = typename decltype(reduce_data)::Type;
3886 reduce_op.eval(bx, reduce_data,
3887 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
3888 {
3889 T t = 0;
3890 for (int n = 0; n < ncomp.n; ++n) {
3891 t += a(i,j,k,n+dcomp.i);
3892 }
3893 return { t };
3894 });
3895 ReduceTuple hv = reduce_data.value(reduce_op);
3896 r = amrex::get<0>(hv);
3897 } else
3898#endif
3899 {
3900 amrex::LoopOnCpu(bx, ncomp.n, [=,&r] (int i, int j, int k, int n) noexcept
3901 {
3902 r += a(i,j,k,n+dcomp.i);
3903 });
3904 }
3905
3906 return r;
3907}
3908
3909template <class T>
3910template <RunOn run_on>
3911T
3912BaseFab<T>::dot (const BaseFab<T>& src, const Box& bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) const noexcept
3913{
3914 AMREX_ASSERT(this->domain.sameType(src.domain));
3915 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
3916 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3917
3918 T r = 0;
3919 Array4<T const> const& d = this->const_array();
3920 Array4<T const> const& s = src.const_array();
3921#ifdef AMREX_USE_GPU
3922 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
3923 ReduceOps<ReduceOpSum> reduce_op;
3924 ReduceData<T> reduce_data(reduce_op);
3925 using ReduceTuple = typename decltype(reduce_data)::Type;
3926 reduce_op.eval(bx, reduce_data,
3927 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
3928 {
3929 T t = 0;
3930 for (int n = 0; n < ncomp.n; ++n) {
3931 t += d(i,j,k,n+dcomp.i) * s(i,j,k,n+scomp.i);
3932 }
3933 return { t };
3934 });
3935 ReduceTuple hv = reduce_data.value(reduce_op);
3936 r = amrex::get<0>(hv);
3937 } else
3938#endif
3939 {
3940 amrex::LoopOnCpu(bx, ncomp.n, [=,&r] (int i, int j, int k, int n) noexcept
3941 {
3942 r += d(i,j,k,n+dcomp.i) * s(i,j,k,n+scomp.i);
3943 });
3944 }
3945
3946 return r;
3947}
3948
3949template <class T>
3950template <RunOn run_on>
3951T
3952BaseFab<T>::dot (const Box& bx, int destcomp, int numcomp) const noexcept
3953{
3954 return dot<run_on>(bx, DestComp{destcomp}, NumComps{numcomp});
3955}
3956
3957
3958template <class T>
3959template <RunOn run_on>
3960T
3961BaseFab<T>::dot (const Box& bx, DestComp dcomp, NumComps ncomp) const noexcept
3962{
3963 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
3964
3965 T r = 0;
3966 Array4<T const> const& a = this->const_array();
3967#ifdef AMREX_USE_GPU
3968 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
3969 ReduceOps<ReduceOpSum> reduce_op;
3970 ReduceData<T> reduce_data(reduce_op);
3971 using ReduceTuple = typename decltype(reduce_data)::Type;
3972 reduce_op.eval(bx, reduce_data,
3973 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
3974 {
3975 T t = 0;
3976 for (int n = 0; n < ncomp.n; ++n) {
3977 t += a(i,j,k,n+dcomp.i)*a(i,j,k,n+dcomp.i);
3978 }
3979 return { t };
3980 });
3981 ReduceTuple hv = reduce_data.value(reduce_op);
3982 r = amrex::get<0>(hv);
3983 } else
3984#endif
3985 {
3986 amrex::LoopOnCpu(bx, ncomp.n, [=,&r] (int i, int j, int k, int n) noexcept
3987 {
3988 r += a(i,j,k,n+dcomp.i)*a(i,j,k,n+dcomp.i);
3989 });
3990 }
3991
3992 return r;
3993}
3994
3995template <class T>
3996template <RunOn run_on>
3997T
3998BaseFab<T>::dotmask (const BaseFab<T>& src, const Box& bx, const BaseFab<int>& mask,
3999 SrcComp scomp, DestComp dcomp, NumComps ncomp) const noexcept
4000{
4001 AMREX_ASSERT(this->domain.sameType(src.domain));
4002 AMREX_ASSERT(this->domain.sameType(mask.domain));
4003 AMREX_ASSERT(scomp.i >= 0 && scomp.i+ncomp.n <= src.nvar);
4004 AMREX_ASSERT(dcomp.i >= 0 && dcomp.i+ncomp.n <= this->nvar);
4005
4006 T r = 0;
4007 Array4<T const> const& d = this->const_array();
4008 Array4<T const> const& s = src.const_array();
4009 Array4<int const> const& m = mask.const_array();
4010#ifdef AMREX_USE_GPU
4011 if (run_on == RunOn::Device && Gpu::inLaunchRegion()) {
4012 ReduceOps<ReduceOpSum> reduce_op;
4013 ReduceData<T> reduce_data(reduce_op);
4014 using ReduceTuple = typename decltype(reduce_data)::Type;
4015 reduce_op.eval(bx, reduce_data,
4016 [=] AMREX_GPU_DEVICE (int i, int j, int k) -> ReduceTuple
4017 {
4018 T t = 0;
4019 T mi = static_cast<T>(static_cast<int>(static_cast<bool>(m(i,j,k))));
4020 for (int n = 0; n < ncomp.n; ++n) {
4021 t += d(i,j,k,n+dcomp.i)*s(i,j,k,n+scomp.i)*mi;
4022 }
4023 return { t };
4024 });
4025 ReduceTuple hv = reduce_data.value(reduce_op);
4026 r = amrex::get<0>(hv);
4027 } else
4028#endif
4029 {
4030 amrex::LoopOnCpu(bx, ncomp.n, [=,&r] (int i, int j, int k, int n) noexcept
4031 {
4032 int mi = static_cast<int>(static_cast<bool>(m(i,j,k)));
4033 r += d(i,j,k,n+dcomp.i)*s(i,j,k,n+scomp.i)*mi;
4034 });
4035 }
4036
4037 return r;
4038}
4039
4040}
4041
4042#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
Convenience header for the core AMReX GPU facilities.
int idir
Definition AMReX_HypreMLABecLap.cpp:1135
Array4< int const > offset
Definition AMReX_HypreMLABecLap.cpp:1131
Real * pdst
Definition AMReX_HypreMLABecLap.cpp:1132
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
GpuArray< Real, 3 > beta
Definition AMReX_MLEBNodeFDLaplacian.cpp:1099
GpuArray< int, 3 > dhi
Definition AMReX_MLEBNodeFDLaplacian.cpp:1100
GpuArray< int, 3 > dlo
Definition AMReX_MLEBNodeFDLaplacian.cpp:1100
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:2582
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:2638
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:3525
T sum(int comp, int numcomp=1) const noexcept
Returns sum of given component of FAB state vector.
Definition AMReX_BaseFab.H:2868
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:2011
BaseFab< T > & divide(T const &val) noexcept
Scalar division on the whole domain and all components.
Definition AMReX_BaseFab.H:3757
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:3665
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:3049
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:1892
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:1947
std::size_t nBytesOwned() const noexcept
Definition AMReX_BaseFab.H:334
BaseFab< T > & copy(const BaseFab< T > &src) noexcept
Definition AMReX_BaseFab.H:3516
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:2703
BaseFab< T > & minus(T const &val) noexcept
Scalar subtraction on the whole domain and all components.
Definition AMReX_BaseFab.H:3617
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:2398
BaseFab< T > & plus(T const &val) noexcept
Scalar addition on the whole domain and all components.
Definition AMReX_BaseFab.H:3547
BaseFab< T > & mult(T const &val, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3703
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:3735
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:2727
void define()
Allocates memory for the BaseFab<T>.
Definition AMReX_BaseFab.H:1596
BaseFab< T > & operator*=(T const &val) noexcept
Definition AMReX_BaseFab.H:3695
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:1774
void clear()
The function returns the BaseFab to the invalid state. The memory is freed.
Definition AMReX_BaseFab.H:1838
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:3197
BaseFab< T > & atomicAdd(const BaseFab< T > &x) noexcept
Atomic FAB addition (a[i] <- a[i] + b[i]).
Definition AMReX_BaseFab.H:2628
int maskEQ(BaseFab< int > &mask, T const &val, int comp=0) const noexcept
Mark cells with value equal to val.
Definition AMReX_BaseFab.H:2490
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:1874
void setVal(T const &x, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3433
BaseFab< T > & operator-=(T const &val) noexcept
Definition AMReX_BaseFab.H:3625
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:3454
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:2307
BaseFab< T > & mult(const BaseFab< T > &src) noexcept
Definition AMReX_BaseFab.H:3719
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:2444
void setValIf(T const &val, const BaseFab< int > &mask) noexcept
Definition AMReX_BaseFab.H:3446
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:3595
void setValIfNot(T const &val, const BaseFab< int > &mask) noexcept
Definition AMReX_BaseFab.H:3468
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:2675
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:1919
BaseFab< T > & negate() noexcept
on the whole domain and all components
Definition AMReX_BaseFab.H:3827
BaseFab< T > & minus(const BaseFab< T > &src) noexcept
Definition AMReX_BaseFab.H:3649
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:3155
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:2268
Elixir elixir() noexcept
Transfer owned storage to an Elixir when running in a GPU launch region.
Definition AMReX_BaseFab.H:1817
BaseFab< T > & operator+=(T const &val) noexcept
Definition AMReX_BaseFab.H:3555
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:3490
BaseFab< T > & minus(T const &val, Box const &bx, DestComp dcomp, NumComps ncomp) noexcept
Do nothing if bx is empty.
Definition AMReX_BaseFab.H:3633
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:3425
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:3773
T max(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2188
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:2764
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:3255
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:3765
std::pair< T, T > minmax(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2226
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:2924
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:3476
BaseFab< T > & mult(T const &val) noexcept
Scalar multiplication on the whole domain and all components.
Definition AMReX_BaseFab.H:3687
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:3805
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:1721
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:2372
BaseFab< T > & protected_divide(const BaseFab< T > &src) noexcept
Divide wherever "src" is "true" or "non-zero".
Definition AMReX_BaseFab.H:3313
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:2908
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:2536
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:2814
BaseFab< T > & plus(const BaseFab< T > &src) noexcept
Definition AMReX_BaseFab.H:3579
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:3789
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:2150
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:1975
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:3383
IntVect minIndex(int comp=0) const noexcept
Definition AMReX_BaseFab.H:2346
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:1983
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:3563
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:2055
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:104
static int deviceId() noexcept
Definition AMReX_GpuDevice.cpp:697
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:672
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:734
amrex_real Real
Floating Point Type for Fields.
Definition AMReX_REAL.H:80
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:1365
__host__ __device__ Dim3 length(Array4< T > const &a) noexcept
Return the spatial extents of an Array4 in Dim3 form.
Definition AMReX_Array4.H:1379
__host__ __device__ Dim3 lbound(Array4< T > const &a) noexcept
Return the inclusive lower bounds of an Array4 in Dim3 form.
Definition AMReX_Array4.H:1351
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
__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:186
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
int nComp(FabArrayBase const &fa)
Convenience wrapper that forwards to fa.nComp().
Definition AMReX_FabArrayBase.cpp:2860
BoxList boxDiff(const Box &b1in, const Box &b2)
Returns a BoxList defining the complement of b2 in b1in.
Definition AMReX_BoxList.cpp:599
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:236
void Abort(const std::string &msg)
Print a fatal-error message to stderr and abort execution.
Definition AMReX.cpp:242
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:289
__host__ __device__ constexpr std::size_t size() const noexcept
Total number of elements in the ArrayND's index region.
Definition AMReX_Array4.H:643
__host__ __device__ T * ptr(idx... i) const noexcept
Multi-index ptr() for accessing pointer to element.
Definition AMReX_Array4.H:561
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