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