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