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