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
153 [[nodiscard]] Real min (int comp,
154 int nghost = 0,
155 bool local = false) const;
162 [[nodiscard]] Real min (const Box& region,
163 int comp,
164 int nghost = 0,
165 bool local = false) const;
176 [[nodiscard]] Real max (int comp,
177 int nghost = 0,
178 bool local = false) const;
185 [[nodiscard]] Real max (const Box& region,
186 int comp,
187 int nghost = 0,
188 bool local = false) const;
195 [[nodiscard]] Real norm0 (int comp = 0, int nghost = 0, bool local = false, bool ignore_covered = false ) const;
197 [[nodiscard]] Real norminf (int comp = 0, int nghost = 0, bool local = false, bool ignore_covered = false ) const {
198 return norm0(comp,nghost,local,ignore_covered);
199 }
200
206 [[nodiscard]] Real norm0 (const iMultiFab& mask, int comp = 0, int nghost = 0, bool local = false) const;
208 [[nodiscard]] Real norminf (const iMultiFab& mask, int comp = 0, int nghost = 0, bool local = false) const {
209 return norm0(mask,comp,nghost,local);
210 }
211
217 [[nodiscard]] Real norm0 (int comp, int ncomp, IntVect const& nghost, bool local = false,
218 bool ignore_covered = false) const;
219
221
228 [[nodiscard]] Vector<Real> norm0 (const Vector<int>& comps, int nghost = 0, bool local = false, bool ignore_covered = false ) const;
230 [[nodiscard]] Vector<Real> norminf (const Vector<int>& comps, int nghost = 0, bool local = false, bool ignore_covered = false) const {
231 return norm0(comps,nghost,local,ignore_covered);
232 }
233
241 [[nodiscard]] Real norm1 (int comp, const Periodicity& period, bool ignore_covered = false) const;
249 [[nodiscard]] Real norm1 (int comp = 0, int ngrow = 0, bool local = false) const;
257 [[nodiscard]] Vector<Real> norm1 (const Vector<int>& comps, int ngrow = 0, bool local = false) const;
265 [[nodiscard]] Real norm2 (int comp = 0) const;
273 [[nodiscard]] Real norm2 (int comp, int numcomp) const;
281 [[nodiscard]] Real norm2 (int comp, const Periodicity& period) const;
289 [[nodiscard]] Vector<Real> norm2 (const Vector<int>& comps) const;
296 [[nodiscard]] Real sum (int comp = 0, bool local = false) const;
303 [[nodiscard]] Real sum (Box const& region, int comp = 0, bool local = false) const;
304
306
314 [[nodiscard]] Real sum_unique (int comp = 0,
315 bool local = false,
316 const Periodicity& period = Periodicity::NonPeriodic()) const;
325 [[nodiscard]] Real sum_unique (Box const& region, int comp = 0, bool local = false) const;
335 void plus (Real val,
336 int comp,
337 int num_comp,
338 int nghost = 0);
344 void plus (Real val,
345 const Box& region,
346 int comp,
347 int num_comp,
348 int nghost = 0);
355 void plus (Real val,
356 int nghost);
364 void plus (Real val,
365 const Box& region,
366 int nghost);
375 void mult (Real val,
376 int comp,
377 int num_comp,
378 int nghost = 0);
386 void mult (Real val,
387 const Box& region,
388 int comp,
389 int num_comp,
390 int nghost = 0);
397 void mult (Real val,
398 int nghost = 0);
406 void mult (Real val,
407 const Box& region,
408 int nghost = 0);
417 void invert (Real numerator,
418 int comp,
419 int num_comp,
420 int nghost = 0);
428 void invert (Real numerator,
429 const Box& region,
430 int comp,
431 int num_comp,
432 int nghost = 0);
439 void invert (Real numerator,
440 int nghost);
448 void invert (Real numerator,
449 const Box& region,
450 int nghost);
458 void negate (int comp,
459 int num_comp,
460 int nghost = 0);
466 void negate (const Box& region,
467 int comp,
468 int num_comp,
469 int nghost = 0);
475 void negate (int nghost = 0);
482 void negate (const Box& region,
483 int nghost = 0);
484
490 [[nodiscard]] IntVect minIndex (int comp,
491 int nghost = 0) const;
492
498 [[nodiscard]] IntVect maxIndex (int comp,
499 int nghost = 0) const;
509 void plus (const MultiFab& mf,
510 int strt_comp,
511 int num_comp,
512 int nghost);
522 void minus (const MultiFab& mf,
523 int strt_comp,
524 int num_comp,
525 int nghost);
536 void divide (const MultiFab& mf,
537 int strt_comp,
538 int num_comp,
539 int nghost);
546 static Real Dot (const MultiFab& x, int xcomp,
547 const MultiFab& y, int ycomp,
548 int numcomp, int nghost, bool local = false);
549
556 static Real Dot (const MultiFab& x, int xcomp,
557 int numcomp, int nghost, bool local = false);
558
562 static Real Dot (const iMultiFab& mask,
563 const MultiFab& x, int xcomp,
564 const MultiFab& y, int ycomp,
565 int numcomp, int nghost, bool local = false);
570 static void Add (MultiFab& dst,
571 const MultiFab& src,
572 int srccomp,
573 int dstcomp,
574 int numcomp,
575 int nghost);
576
577 static void Add (MultiFab& dst,
578 const MultiFab& src,
579 int srccomp,
580 int dstcomp,
581 int numcomp,
582 const IntVect& nghost);
583
589 [[nodiscard]] MultiFab deepCopy () const;
590
596 static void Copy (MultiFab& dst,
597 const MultiFab& src,
598 int srccomp,
599 int dstcomp,
600 int numcomp,
601 int nghost);
602
603 static void Copy (MultiFab& dst,
604 const MultiFab& src,
605 int srccomp,
606 int dstcomp,
607 int numcomp,
608 const IntVect& nghost);
609
615 static void Swap (MultiFab& dst,
616 MultiFab& src,
617 int srccomp,
618 int dstcomp,
619 int numcomp,
620 int nghost);
621
622 static void Swap (MultiFab& dst,
623 MultiFab& src,
624 int srccomp,
625 int dstcomp,
626 int numcomp,
627 const IntVect& nghost);
628
633 static void Subtract (MultiFab& dst,
634 const MultiFab& src,
635 int srccomp,
636 int dstcomp,
637 int numcomp,
638 int nghost);
639
640 static void Subtract (MultiFab& dst,
641 const MultiFab& src,
642 int srccomp,
643 int dstcomp,
644 int numcomp,
645 const IntVect& nghost);
650 static void Multiply (MultiFab& dst,
651 const MultiFab& src,
652 int srccomp,
653 int dstcomp,
654 int numcomp,
655 int nghost);
656
657 static void Multiply (MultiFab& dst,
658 const MultiFab& src,
659 int srccomp,
660 int dstcomp,
661 int numcomp,
662 const IntVect& nghost);
667 static void Divide (MultiFab& dst,
668 const MultiFab& src,
669 int srccomp,
670 int dstcomp,
671 int numcomp,
672 int nghost);
673
674 static void Divide (MultiFab& dst,
675 const MultiFab& src,
676 int srccomp,
677 int dstcomp,
678 int numcomp,
679 const IntVect& nghost);
683 static void Saxpy (MultiFab& dst,
684 Real a,
685 const MultiFab& src,
686 int srccomp,
687 int dstcomp,
688 int numcomp,
689 int nghost);
690
692
696 static void Xpay (MultiFab& dst,
697 Real a,
698 const MultiFab& src,
699 int srccomp,
700 int dstcomp,
701 int numcomp,
702 int nghost);
703
705
709 static void LinComb (MultiFab& dst,
710 Real a,
711 const MultiFab& x,
712 int xcomp,
713 Real b,
714 const MultiFab& y,
715 int ycomp,
716 int dstcomp,
717 int numcomp,
718 int nghost);
719
721
725 static void AddProduct (MultiFab& dst,
726 const MultiFab& src1,
727 int comp1,
728 const MultiFab& src2,
729 int comp2,
730 int dstcomp,
731 int numcomp,
732 int nghost);
733
734 static void AddProduct (MultiFab& dst,
735 const MultiFab& src1,
736 int comp1,
737 const MultiFab& src2,
738 int comp2,
739 int dstcomp,
740 int numcomp,
741 const IntVect& nghost);
742
752 [[nodiscard]] bool is_finite (bool local=false) const;
753 [[nodiscard]] bool is_finite (int scomp, int ncomp, int ngrow = 0, bool local=false) const;
754 [[nodiscard]] bool is_finite (int scomp, int ncomp, const IntVect& ngrow, bool local=false) const;
755
765 [[nodiscard]] bool contains_nan (bool local=false) const;
766 [[nodiscard]] bool contains_nan (int scomp, int ncomp, int ngrow = 0, bool local=false) const;
767 [[nodiscard]] bool contains_nan (int scomp, int ncomp, const IntVect& ngrow, bool local=false) const;
776 [[nodiscard]] bool contains_inf (bool local=false) const;
777 [[nodiscard]] bool contains_inf (int scomp, int ncomp, int ngrow = 0, bool local=false) const;
778 [[nodiscard]] bool contains_inf (int scomp, int ncomp, const IntVect& ngrow, bool local=false) const;
779
783 [[nodiscard]] std::unique_ptr<MultiFab> OverlapMask (const Periodicity& period = Periodicity::NonPeriodic()) const;
785 [[nodiscard]] std::unique_ptr<iMultiFab> OwnerMask (const Periodicity& period = Periodicity::NonPeriodic()) const;
786
788 void AverageSync (const Periodicity& period = Periodicity::NonPeriodic());
790 void WeightedSync (const MultiFab& wgt, const Periodicity& period = Periodicity::NonPeriodic());
808 void OverrideSync (const iMultiFab& msk, const Periodicity& period = Periodicity::NonPeriodic());
809
811
812 static void Initialize ();
813 static void Finalize ();
814
815private:
816 //
818 using CopyComTagsContainer = FabArrayBase::CopyComTagsContainer;
819 using MapOfCopyComTagContainers = FabArrayBase::MapOfCopyComTagContainers;
820
821 void initVal ();
822};
823
824#ifndef _MSC_VER
825inline void GccPlacaterMF ()
826{
827 std::allocator<MultiFab*> a1;
828 std::allocator<MultiFab const*> a2;
829 std::allocator<FabArray<FArrayBox>*> a3;
830 std::allocator<FabArray<FArrayBox> const*> a4;
831
836}
837#endif
838
839}
840
841#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:132
A collection of Boxes stored in an Array.
Definition AMReX_BoxArray.H:568
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:350
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
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
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
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:230
Real norminf(int comp=0, int nghost=0, bool local=false, bool ignore_covered=false) const
Definition AMReX_MultiFab.H:197
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:208
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:139
DefaultFabFactory< FArrayBox > FArrayBoxFactory
Definition AMReX_FArrayBox.H:460
void GccPlacaterMF()
Definition AMReX_MultiFab.H:825
FabArray memory allocation information.
Definition AMReX_FabArray.H:66