Block-Structured AMR Software Framework
 
Loading...
Searching...
No Matches
AMReX_MultiFab.H
Go to the documentation of this file.
1
2#ifndef BL_MULTIFAB_H
3#define BL_MULTIFAB_H
4#include <AMReX_Config.H>
5
6#include <AMReX_BLassert.H>
7#include <AMReX_BaseFab.H>
8#include <AMReX_FArrayBox.H>
9#include <AMReX_FabArray.H>
11#include <AMReX_Periodicity.H>
12#include <AMReX_NonLocalBC.H>
13
14#include <cstdint>
15
16namespace amrex
17{
18
20
21class iMultiFab;
22
38 :
39 public FabArray<FArrayBox>
40{
41public:
42
49 MultiFab () noexcept;
50
59 explicit MultiFab (Arena* a) noexcept;
60
76 MultiFab (const BoxArray& bxs,
77 const DistributionMapping& dm,
78 int ncomp,
79 int ngrow,
80#ifdef AMREX_STRICT_MODE
81 const MFInfo& info,
82 const FabFactory<FArrayBox>& factory);
83#else
84 const MFInfo& info = MFInfo(),
85 const FabFactory<FArrayBox>& factory = FArrayBoxFactory());
86#endif
87
88 MultiFab (const BoxArray& bxs,
89 const DistributionMapping& dm,
90 int ncomp,
91 const IntVect& ngrow,
92#ifdef AMREX_STRICT_MODE
93 const MFInfo& info,
94 const FabFactory<FArrayBox>& factory);
95#else
96 const MFInfo& info = MFInfo(),
97 const FabFactory<FArrayBox>& factory = FArrayBoxFactory());
98#endif
99
107 MultiFab (const MultiFab& rhs, MakeType maketype, int scomp, int ncomp);
108
109 ~MultiFab ();
110
111 MultiFab (MultiFab&& rhs) noexcept;
112 MultiFab& operator= (MultiFab&& rhs) noexcept = default;
113
114 MultiFab (const MultiFab& rhs) = delete;
115 MultiFab& operator= (const MultiFab& rhs) = delete;
116
117 void define (const BoxArray& bxs,
118 const DistributionMapping& dm,
119 int nvar,
120 int ngrow,
121#ifdef AMREX_STRICT_MODE
122 const MFInfo& info,
123 const FabFactory<FArrayBox>& factory);
124#else
125 const MFInfo& info = MFInfo(),
126 const FabFactory<FArrayBox>& factory = FArrayBoxFactory());
127#endif
128
129 void define (const BoxArray& bxs,
130 const DistributionMapping& dm,
131 int nvar,
132 const IntVect& ngrow,
133#ifdef AMREX_STRICT_MODE
134 const MFInfo& info,
135 const FabFactory<FArrayBox>& factory);
136#else
137 const MFInfo& info = MFInfo(),
138 const FabFactory<FArrayBox>& factory = FArrayBoxFactory());
139#endif
140
142
151 [[nodiscard]] Real min (int comp,
152 int nghost = 0,
153 bool local = false) const;
158 [[nodiscard]] Real min (const Box& region,
159 int comp,
160 int nghost = 0,
161 bool local = false) const;
170 [[nodiscard]] Real max (int comp,
171 int nghost = 0,
172 bool local = false) const;
177 [[nodiscard]] Real max (const Box& region,
178 int comp,
179 int nghost = 0,
180 bool local = false) const;
185 [[nodiscard]] Real norm0 (int comp = 0, int nghost = 0, bool local = false, bool ignore_covered = false ) const;
187 [[nodiscard]] Real norminf (int comp = 0, int nghost = 0, bool local = false, bool ignore_covered = false ) const {
188 return norm0(comp,nghost,local,ignore_covered);
189 }
190
191 [[nodiscard]] Real norm0 (const iMultiFab& mask, int comp = 0, int nghost = 0, bool local = false) const;
193 [[nodiscard]] Real norminf (const iMultiFab& mask, int comp = 0, int nghost = 0, bool local = false) const {
194 return norm0(mask,comp,nghost,local);
195 }
196
197 [[nodiscard]] Real norm0 (int comp, int ncomp, IntVect const& nghost, bool local = false,
198 bool ignore_covered = false) const;
199
201
206 [[nodiscard]] Vector<Real> norm0 (const Vector<int>& comps, int nghost = 0, bool local = false, bool ignore_covered = false ) const;
208 [[nodiscard]] Vector<Real> norminf (const Vector<int>& comps, int nghost = 0, bool local = false, bool ignore_covered = false) const {
209 return norm0(comps,nghost,local,ignore_covered);
210 }
211
217 [[nodiscard]] Real norm1 (int comp, const Periodicity& period, bool ignore_covered = false) const;
223 [[nodiscard]] Real norm1 (int comp = 0, int ngrow = 0, bool local = false) const;
229 [[nodiscard]] Vector<Real> norm1 (const Vector<int>& comps, int ngrow = 0, bool local = false) const;
235 [[nodiscard]] Real norm2 (int comp = 0) const;
241 [[nodiscard]] Real norm2 (int comp, int numcomp) const;
247 [[nodiscard]] Real norm2 (int comp, const Periodicity& period) const;
253 [[nodiscard]] Vector<Real> norm2 (const Vector<int>& comps) const;
258 [[nodiscard]] Real sum (int comp = 0, bool local = false) const;
263 [[nodiscard]] Real sum (Box const& region, int comp = 0, bool local = false) const;
264
266
272 [[nodiscard]] Real sum_unique (int comp = 0,
273 bool local = false,
274 const Periodicity& period = Periodicity::NonPeriodic()) const;
281 [[nodiscard]] Real sum_unique (Box const& region, int comp = 0, bool local = false) const;
291 void plus (Real val,
292 int comp,
293 int num_comp,
294 int nghost = 0);
300 void plus (Real val,
301 const Box& region,
302 int comp,
303 int num_comp,
304 int nghost = 0);
311 void plus (Real val,
312 int nghost);
320 void plus (Real val,
321 const Box& region,
322 int nghost);
331 void mult (Real val,
332 int comp,
333 int num_comp,
334 int nghost = 0);
342 void mult (Real val,
343 const Box& region,
344 int comp,
345 int num_comp,
346 int nghost = 0);
353 void mult (Real val,
354 int nghost = 0);
362 void mult (Real val,
363 const Box& region,
364 int nghost = 0);
373 void invert (Real numerator,
374 int comp,
375 int num_comp,
376 int nghost = 0);
384 void invert (Real numerator,
385 const Box& region,
386 int comp,
387 int num_comp,
388 int nghost = 0);
395 void invert (Real numerator,
396 int nghost);
404 void invert (Real numerator,
405 const Box& region,
406 int nghost);
414 void negate (int comp,
415 int num_comp,
416 int nghost = 0);
422 void negate (const Box& region,
423 int comp,
424 int num_comp,
425 int nghost = 0);
431 void negate (int nghost = 0);
438 void negate (const Box& region,
439 int nghost = 0);
440
442 [[nodiscard]] IntVect minIndex (int comp,
443 int nghost = 0) const;
444
446 [[nodiscard]] IntVect maxIndex (int comp,
447 int nghost = 0) const;
457 void plus (const MultiFab& mf,
458 int strt_comp,
459 int num_comp,
460 int nghost);
470 void minus (const MultiFab& mf,
471 int strt_comp,
472 int num_comp,
473 int nghost);
484 void divide (const MultiFab& mf,
485 int strt_comp,
486 int num_comp,
487 int nghost);
492 static Real Dot (const MultiFab& x, int xcomp,
493 const MultiFab& y, int ycomp,
494 int numcomp, int nghost, bool local = false);
495
500 static Real Dot (const MultiFab& x, int xcomp,
501 int numcomp, int nghost, bool local = false);
502
504 static Real Dot (const iMultiFab& mask,
505 const MultiFab& x, int xcomp,
506 const MultiFab& y, int ycomp,
507 int numcomp, int nghost, bool local = false);
512 static void Add (MultiFab& dst,
513 const MultiFab& src,
514 int srccomp,
515 int dstcomp,
516 int numcomp,
517 int nghost);
518
519 static void Add (MultiFab& dst,
520 const MultiFab& src,
521 int srccomp,
522 int dstcomp,
523 int numcomp,
524 const IntVect& nghost);
525
531 [[nodiscard]] MultiFab deepCopy () const;
532
538 static void Copy (MultiFab& dst,
539 const MultiFab& src,
540 int srccomp,
541 int dstcomp,
542 int numcomp,
543 int nghost);
544
545 static void Copy (MultiFab& dst,
546 const MultiFab& src,
547 int srccomp,
548 int dstcomp,
549 int numcomp,
550 const IntVect& nghost);
551
557 static void Swap (MultiFab& dst,
558 MultiFab& src,
559 int srccomp,
560 int dstcomp,
561 int numcomp,
562 int nghost);
563
564 static void Swap (MultiFab& dst,
565 MultiFab& src,
566 int srccomp,
567 int dstcomp,
568 int numcomp,
569 const IntVect& nghost);
570
575 static void Subtract (MultiFab& dst,
576 const MultiFab& src,
577 int srccomp,
578 int dstcomp,
579 int numcomp,
580 int nghost);
581
582 static void Subtract (MultiFab& dst,
583 const MultiFab& src,
584 int srccomp,
585 int dstcomp,
586 int numcomp,
587 const IntVect& nghost);
592 static void Multiply (MultiFab& dst,
593 const MultiFab& src,
594 int srccomp,
595 int dstcomp,
596 int numcomp,
597 int nghost);
598
599 static void Multiply (MultiFab& dst,
600 const MultiFab& src,
601 int srccomp,
602 int dstcomp,
603 int numcomp,
604 const IntVect& nghost);
609 static void Divide (MultiFab& dst,
610 const MultiFab& src,
611 int srccomp,
612 int dstcomp,
613 int numcomp,
614 int nghost);
615
616 static void Divide (MultiFab& dst,
617 const MultiFab& src,
618 int srccomp,
619 int dstcomp,
620 int numcomp,
621 const IntVect& nghost);
625 static void Saxpy (MultiFab& dst,
626 Real a,
627 const MultiFab& src,
628 int srccomp,
629 int dstcomp,
630 int numcomp,
631 int nghost);
632
634
638 static void Xpay (MultiFab& dst,
639 Real a,
640 const MultiFab& src,
641 int srccomp,
642 int dstcomp,
643 int numcomp,
644 int nghost);
645
647
651 static void LinComb (MultiFab& dst,
652 Real a,
653 const MultiFab& x,
654 int xcomp,
655 Real b,
656 const MultiFab& y,
657 int ycomp,
658 int dstcomp,
659 int numcomp,
660 int nghost);
661
663
667 static void AddProduct (MultiFab& dst,
668 const MultiFab& src1,
669 int comp1,
670 const MultiFab& src2,
671 int comp2,
672 int dstcomp,
673 int numcomp,
674 int nghost);
675
676 static void AddProduct (MultiFab& dst,
677 const MultiFab& src1,
678 int comp1,
679 const MultiFab& src2,
680 int comp2,
681 int dstcomp,
682 int numcomp,
683 const IntVect& nghost);
684
692 [[nodiscard]] bool is_finite (bool local=false) const;
693 [[nodiscard]] bool is_finite (int scomp, int ncomp, int ngrow = 0, bool local=false) const;
694 [[nodiscard]] bool is_finite (int scomp, int ncomp, const IntVect& ngrow, bool local=false) const;
695
703 [[nodiscard]] bool contains_nan (bool local=false) const;
704 [[nodiscard]] bool contains_nan (int scomp, int ncomp, int ngrow = 0, bool local=false) const;
705 [[nodiscard]] bool contains_nan (int scomp, int ncomp, const IntVect& ngrow, bool local=false) const;
712 [[nodiscard]] bool contains_inf (bool local=false) const;
713 [[nodiscard]] bool contains_inf (int scomp, int ncomp, int ngrow = 0, bool local=false) const;
714 [[nodiscard]] bool contains_inf (int scomp, int ncomp, const IntVect& ngrow, bool local=false) const;
715
719 [[nodiscard]] std::unique_ptr<MultiFab> OverlapMask (const Periodicity& period = Periodicity::NonPeriodic()) const;
721 [[nodiscard]] std::unique_ptr<iMultiFab> OwnerMask (const Periodicity& period = Periodicity::NonPeriodic()) const;
722
724 void AverageSync (const Periodicity& period = Periodicity::NonPeriodic());
726 void WeightedSync (const MultiFab& wgt, const Periodicity& period = Periodicity::NonPeriodic());
744 void OverrideSync (const iMultiFab& msk, const Periodicity& period = Periodicity::NonPeriodic());
745
747
748 static void Initialize ();
749 static void Finalize ();
750
751private:
752 //
756
757 void initVal ();
758};
759
760#ifndef _MSC_VER
761inline void GccPlacaterMF ()
762{
763 std::allocator<MultiFab*> a1;
764 std::allocator<MultiFab const*> a2;
765 std::allocator<FabArray<FArrayBox>*> a3;
766 std::allocator<FabArray<FArrayBox> const*> a4;
767
772}
773#endif
774
775}
776
777#endif /*BL_MULTIFAB_H*/
Array4< int const > mask
Definition AMReX_InterpFaceRegister.cpp:93
A virtual base class for objects that manage their own dynamic memory allocation.
Definition AMReX_Arena.H:105
A collection of Boxes stored in an Array.
Definition AMReX_BoxArray.H:567
Calculates the distribution of FABs to MPI processes.
Definition AMReX_DistributionMapping.H:43
A Fortran Array of REALs.
Definition AMReX_FArrayBox.H:231
CopyComTag::CopyComTagsContainer CopyComTagsContainer
Definition AMReX_FabArrayBase.H:220
CopyComTag::MapOfCopyComTagContainers MapOfCopyComTagContainers
Definition AMReX_FabArrayBase.H:221
An Array of FortranArrayBox(FAB)-like Objects.
Definition AMReX_FabArray.H:347
Definition AMReX_FabFactory.H:50
A collection (stored as an array) of FArrayBox objects.
Definition AMReX_MultiFab.H:40
static void LinComb(MultiFab &dst, Real a, const MultiFab &x, int xcomp, Real b, const MultiFab &y, int ycomp, int dstcomp, int numcomp, int nghost)
dst = a*x + b*y
Definition AMReX_MultiFab.cpp:306
void define(const BoxArray &bxs, const DistributionMapping &dm, int nvar, int ngrow, const MFInfo &info=MFInfo(), const FabFactory< FArrayBox > &factory=FArrayBoxFactory())
Definition AMReX_MultiFab.cpp:519
FabArrayBase::CopyComTagsContainer CopyComTagsContainer
Some useful typedefs.
Definition AMReX_MultiFab.H:754
void invert(Real numerator, int comp, int num_comp, int nghost=0)
Replaces the value of each cell in the specified subregion of the MultiFab with its reciprocal multip...
Definition AMReX_MultiFab.cpp:1437
void minus(const MultiFab &mf, int strt_comp, int num_comp, int nghost)
This function subtracts the values of the cells in mf from the corresponding cells of this MultiFab....
Definition AMReX_MultiFab.cpp:1378
FabArrayBase::MapOfCopyComTagContainers MapOfCopyComTagContainers
Definition AMReX_MultiFab.H:755
Real norm0(int comp=0, int nghost=0, bool local=false, bool ignore_covered=false) const
Returns the maximum absolute value contained in component comp of the MultiFab.
Definition AMReX_MultiFab.cpp:1035
static void Swap(MultiFab &dst, MultiFab &src, int srccomp, int dstcomp, int numcomp, int nghost)
Swap from src to dst including nghost ghost cells. The two MultiFabs MUST have the same underlying Bo...
Definition AMReX_MultiFab.cpp:213
std::unique_ptr< iMultiFab > OwnerMask(const Periodicity &period=Periodicity::NonPeriodic()) const
Owner is the grid with the lowest grid number containing the data.
Definition AMReX_MultiFab.cpp:1541
static void Add(MultiFab &dst, const MultiFab &src, int srccomp, int dstcomp, int numcomp, int nghost)
Add src to dst including nghost ghost cells. The two MultiFabs MUST have the same underlying BoxArray...
Definition AMReX_MultiFab.cpp:151
bool contains_nan(bool local=false) const
Are there any NaNs in the MF? This may return false, even if the MF contains NaNs,...
Definition AMReX_MultiFab.cpp:664
static void Saxpy(MultiFab &dst, Real a, const MultiFab &src, int srccomp, int dstcomp, int numcomp, int nghost)
dst += a*src
Definition AMReX_MultiFab.cpp:292
Real min(int comp, int nghost=0, bool local=false) const
Returns the minimum value contained in component comp of the MultiFab.
Definition AMReX_MultiFab.cpp:725
Real norm1(int comp, const Periodicity &period, bool ignore_covered=false) const
Returns the L1 norm of component comp over the MultiFab.
Definition AMReX_MultiFab.cpp:1143
static void Subtract(MultiFab &dst, const MultiFab &src, int srccomp, int dstcomp, int numcomp, int nghost)
Subtract src from dst including nghost ghost cells. The two MultiFabs MUST have the same underlying B...
Definition AMReX_MultiFab.cpp:232
static void Multiply(MultiFab &dst, const MultiFab &src, int srccomp, int dstcomp, int numcomp, int nghost)
Multiply dst by src including nghost ghost cells. The two MultiFabs MUST have the same underlying Box...
Definition AMReX_MultiFab.cpp:252
static void AddProduct(MultiFab &dst, const MultiFab &src1, int comp1, const MultiFab &src2, int comp2, int dstcomp, int numcomp, int nghost)
dst += src1*src2
Definition AMReX_MultiFab.cpp:315
bool is_finite(bool local=false) const
Are the numbers in the MF finite (i.e., neither nan nor inf)? This may return true,...
Definition AMReX_MultiFab.cpp:609
static void Initialize()
Definition AMReX_MultiFab.cpp:422
void AverageSync(const Periodicity &period=Periodicity::NonPeriodic())
Sync up nodal data via averaging.
Definition AMReX_MultiFab.cpp:1547
bool contains_inf(bool local=false) const
Are there any Infs in the MF? This may return false, even if the MF contains Infs,...
Definition AMReX_MultiFab.cpp:719
static void Xpay(MultiFab &dst, Real a, const MultiFab &src, int srccomp, int dstcomp, int numcomp, int nghost)
dst = src + a*dst
Definition AMReX_MultiFab.cpp:299
void initVal()
Definition AMReX_MultiFab.cpp:543
static void Divide(MultiFab &dst, const MultiFab &src, int srccomp, int dstcomp, int numcomp, int nghost)
Divide dst by src including nghost ghost cells. The two MultiFabs MUST have the same underlying BoxAr...
Definition AMReX_MultiFab.cpp:272
void mult(Real val, int comp, int num_comp, int nghost=0)
Scales the value of each cell in the specified subregion of the MultiFab by the scalar val (a[i] <- a...
Definition AMReX_MultiFab.cpp:1417
void plus(Real val, int comp, int num_comp, int nghost=0)
Adds the scalar value val to the value of each cell in the specified subregion of the MultiFab.
Definition AMReX_MultiFab.cpp:1390
static void Copy(MultiFab &dst, const MultiFab &src, int srccomp, int dstcomp, int numcomp, int nghost)
Copy from src to dst including nghost ghost cells. The two MultiFabs MUST have the same underlying Bo...
Definition AMReX_MultiFab.cpp:193
Real max(int comp, int nghost=0, bool local=false) const
Returns the maximum value contained in component comp of the MultiFab.
Definition AMReX_MultiFab.cpp:854
MultiFab() noexcept
Constructs an empty MultiFab.
Definition AMReX_MultiFab.cpp:443
void divide(const MultiFab &mf, int strt_comp, int num_comp, int nghost)
This function divides the values of the cells in mf from the corresponding cells of this MultiFab....
Definition AMReX_MultiFab.cpp:1384
void negate(int comp, int num_comp, int nghost=0)
Negates the value of each cell in the specified subregion of the MultiFab. The subregion consists of ...
Definition AMReX_MultiFab.cpp:1457
MultiFab & operator=(MultiFab &&rhs) noexcept=default
MultiFab deepCopy() const
Definition AMReX_MultiFab.cpp:171
void WeightedSync(const MultiFab &wgt, const Periodicity &period=Periodicity::NonPeriodic())
Sync up nodal data with weights.
Definition AMReX_MultiFab.cpp:1558
MultiFab(const MultiFab &rhs)=delete
~MultiFab()
Definition AMReX_MultiFab.cpp:504
static void Finalize()
Definition AMReX_MultiFab.cpp:438
std::unique_ptr< MultiFab > OverlapMask(const Periodicity &period=Periodicity::NonPeriodic()) const
Return a mask indicating how many duplicates are in each point.
Definition AMReX_MultiFab.cpp:1475
void OverrideSync(const iMultiFab &msk, const Periodicity &period=Periodicity::NonPeriodic())
Synchronize nodal data. The synchronization will override valid regions by the intersecting valid reg...
Definition AMReX_MultiFab.cpp:1578
This provides length of period for periodic domains. 0 means it is not periodic in that direction....
Definition AMReX_Periodicity.H:17
static const Periodicity & NonPeriodic() noexcept
Definition AMReX_Periodicity.cpp:52
This class is a thin wrapper around std::vector. Unlike vector, Vector::operator[] provides bound che...
Definition AMReX_Vector.H:28
A Collection of IArrayBoxes.
Definition AMReX_iMultiFab.H:34
amrex_real Real
Floating Point Type for Fields.
Definition AMReX_REAL.H:79
static Real Dot(const MultiFab &x, int xcomp, const MultiFab &y, int ycomp, int numcomp, int nghost, bool local=false)
Returns the dot product of two MultiFabs.
Definition AMReX_MultiFab.cpp:37
Real sum(int comp=0, bool local=false) const
Returns the sum of component "comp" over the MultiFab – no ghost cells are included.
Definition AMReX_MultiFab.cpp:1223
IntVect maxIndex(int comp, int nghost=0) const
Definition AMReX_MultiFab.cpp:1021
Vector< Real > norminf(const Vector< int > &comps, int nghost=0, bool local=false, bool ignore_covered=false) const
Definition AMReX_MultiFab.H:208
Real norminf(int comp=0, int nghost=0, bool local=false, bool ignore_covered=false) const
Definition AMReX_MultiFab.H:187
Real norm2(int comp=0) const
Returns the L2 norm of component comp over the MultiFab. No ghost cells are used.
Definition AMReX_MultiFab.cpp:1065
IntVect minIndex(int comp, int nghost=0) const
Definition AMReX_MultiFab.cpp:1013
Real norminf(const iMultiFab &mask, int comp=0, int nghost=0, bool local=false) const
Definition AMReX_MultiFab.H:193
Real sum_unique(int comp=0, bool local=false, const Periodicity &period=Periodicity::NonPeriodic()) const
Same as sum with local =false, but for non-cell-centered data, this only adds non-unique points that ...
Definition AMReX_MultiFab.cpp:1271
Definition AMReX_Amr.cpp:49
MakeType
Definition AMReX_MakeType.H:7
__host__ __device__ void ignore_unused(const Ts &...)
This shuts up the compiler about unused variables.
Definition AMReX.H:138
DefaultFabFactory< FArrayBox > FArrayBoxFactory
Definition AMReX_FArrayBox.H:460
void GccPlacaterMF()
Definition AMReX_MultiFab.H:761
FabArray memory allocation information.
Definition AMReX_FabArray.H:66