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