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
36 :
37 public FabArray<FArrayBox>
38{
39public:
40
47 MultiFab () noexcept;
48
57 explicit MultiFab (Arena* a) noexcept;
58
74 MultiFab (const BoxArray& bxs,
75 const DistributionMapping& dm,
76 int ncomp,
77 int ngrow,
78#ifdef AMREX_STRICT_MODE
79 const MFInfo& info,
80 const FabFactory<FArrayBox>& factory);
81#else
82 const MFInfo& info = MFInfo(),
83 const FabFactory<FArrayBox>& factory = FArrayBoxFactory());
84#endif
85
86 MultiFab (const BoxArray& bxs,
87 const DistributionMapping& dm,
88 int ncomp,
89 const IntVect& ngrow,
90#ifdef AMREX_STRICT_MODE
91 const MFInfo& info,
92 const FabFactory<FArrayBox>& factory);
93#else
94 const MFInfo& info = MFInfo(),
95 const FabFactory<FArrayBox>& factory = FArrayBoxFactory());
96#endif
97
105 MultiFab (const MultiFab& rhs, MakeType maketype, int scomp, int ncomp);
106
107 ~MultiFab ();
108
109 MultiFab (MultiFab&& rhs) noexcept;
110 MultiFab& operator= (MultiFab&& rhs) noexcept = default;
111
112 MultiFab (const MultiFab& rhs) = delete;
113 MultiFab& operator= (const MultiFab& rhs) = delete;
114
115 void define (const BoxArray& bxs,
116 const DistributionMapping& dm,
117 int nvar,
118 int ngrow,
119#ifdef AMREX_STRICT_MODE
120 const MFInfo& info,
121 const FabFactory<FArrayBox>& factory);
122#else
123 const MFInfo& info = MFInfo(),
124 const FabFactory<FArrayBox>& factory = FArrayBoxFactory());
125#endif
126
127 void define (const BoxArray& bxs,
128 const DistributionMapping& dm,
129 int nvar,
130 const IntVect& ngrow,
131#ifdef AMREX_STRICT_MODE
132 const MFInfo& info,
133 const FabFactory<FArrayBox>& factory);
134#else
135 const MFInfo& info = MFInfo(),
136 const FabFactory<FArrayBox>& factory = FArrayBoxFactory());
137#endif
138
139 MultiFab& operator= (Real r);
140
149 [[nodiscard]] Real min (int comp,
150 int nghost = 0,
151 bool local = false) const;
156 [[nodiscard]] Real min (const Box& region,
157 int comp,
158 int nghost = 0,
159 bool local = false) const;
168 [[nodiscard]] Real max (int comp,
169 int nghost = 0,
170 bool local = false) const;
175 [[nodiscard]] Real max (const Box& region,
176 int comp,
177 int nghost = 0,
178 bool local = false) const;
183 [[nodiscard]] Real norm0 (int comp = 0, int nghost = 0, bool local = false, bool ignore_covered = false ) const;
184 [[nodiscard]] Real norminf (int comp = 0, int nghost = 0, bool local = false, bool ignore_covered = false ) const {
185 return norm0(comp,nghost,local,ignore_covered);
186 }
187
188 [[nodiscard]] Real norm0 (const iMultiFab& mask, int comp = 0, int nghost = 0, bool local = false) const;
189 [[nodiscard]] Real norminf (const iMultiFab& mask, int comp = 0, int nghost = 0, bool local = false) const {
190 return norm0(mask,comp,nghost,local);
191 }
192
193 [[nodiscard]] Real norm0 (int comp, int ncomp, IntVect const& nghost, bool local = false,
194 bool ignore_covered = false) const;
195
197
202 [[nodiscard]] Vector<Real> norm0 (const Vector<int>& comps, int nghost = 0, bool local = false, bool ignore_covered = false ) const;
203 [[nodiscard]] Vector<Real> norminf (const Vector<int>& comps, int nghost = 0, bool local = false, bool ignore_covered = false) const {
204 return norm0(comps,nghost,local,ignore_covered);
205 }
206
212 [[nodiscard]] Real norm1 (int comp, const Periodicity& period, bool ignore_covered = false) const;
217 [[nodiscard]] Real norm1 (int comp = 0, int ngrow = 0, bool local = false) const;
222 [[nodiscard]] Vector<Real> norm1 (const Vector<int>& comps, int ngrow = 0, bool local = false) const;
227 [[nodiscard]] Real norm2 (int comp = 0) const;
232 [[nodiscard]] Real norm2 (int comp, int numcomp) const;
237 [[nodiscard]] Real norm2 (int comp, const Periodicity& period) const;
242 [[nodiscard]] Vector<Real> norm2 (const Vector<int>& comps) const;
246 [[nodiscard]] Real sum (int comp = 0, bool local = false) const;
250 [[nodiscard]] Real sum (Box const& region, int comp = 0, bool local = false) const;
251
253
258 [[nodiscard]] Real sum_unique (int comp = 0,
259 bool local = false,
260 const Periodicity& period = Periodicity::NonPeriodic()) const;
267 [[nodiscard]] Real sum_unique (Box const& region, int comp = 0, bool local = false) const;
277 void plus (Real val,
278 int comp,
279 int num_comp,
280 int nghost = 0);
286 void plus (Real val,
287 const Box& region,
288 int comp,
289 int num_comp,
290 int nghost = 0);
297 void plus (Real val,
298 int nghost);
306 void plus (Real val,
307 const Box& region,
308 int nghost);
317 void mult (Real val,
318 int comp,
319 int num_comp,
320 int nghost = 0);
328 void mult (Real val,
329 const Box& region,
330 int comp,
331 int num_comp,
332 int nghost = 0);
339 void mult (Real val,
340 int nghost = 0);
348 void mult (Real val,
349 const Box& region,
350 int nghost = 0);
359 void invert (Real numerator,
360 int comp,
361 int num_comp,
362 int nghost = 0);
370 void invert (Real numerator,
371 const Box& region,
372 int comp,
373 int num_comp,
374 int nghost = 0);
381 void invert (Real numerator,
382 int nghost);
390 void invert (Real numerator,
391 const Box& region,
392 int nghost);
400 void negate (int comp,
401 int num_comp,
402 int nghost = 0);
408 void negate (const Box& region,
409 int comp,
410 int num_comp,
411 int nghost = 0);
417 void negate (int nghost = 0);
424 void negate (const Box& region,
425 int nghost = 0);
426
427 [[nodiscard]] IntVect minIndex (int comp,
428 int nghost = 0) const;
429
430 [[nodiscard]] IntVect maxIndex (int comp,
431 int nghost = 0) const;
441 void plus (const MultiFab& mf,
442 int strt_comp,
443 int num_comp,
444 int nghost);
454 void minus (const MultiFab& mf,
455 int strt_comp,
456 int num_comp,
457 int nghost);
468 void divide (const MultiFab& mf,
469 int strt_comp,
470 int num_comp,
471 int nghost);
475 static Real Dot (const MultiFab& x, int xcomp,
476 const MultiFab& y, int ycomp,
477 int numcomp, int nghost, bool local = false);
478
482 static Real Dot (const MultiFab& x, int xcomp,
483 int numcomp, int nghost, bool local = false);
484
485 static Real Dot (const iMultiFab& mask,
486 const MultiFab& x, int xcomp,
487 const MultiFab& y, int ycomp,
488 int numcomp, int nghost, bool local = false);
493 static void Add (MultiFab& dst,
494 const MultiFab& src,
495 int srccomp,
496 int dstcomp,
497 int numcomp,
498 int nghost);
499
500 static void Add (MultiFab& dst,
501 const MultiFab& src,
502 int srccomp,
503 int dstcomp,
504 int numcomp,
505 const IntVect& nghost);
506
512 [[nodiscard]] MultiFab deepCopy () const;
513
519 static void Copy (MultiFab& dst,
520 const MultiFab& src,
521 int srccomp,
522 int dstcomp,
523 int numcomp,
524 int nghost);
525
526 static void Copy (MultiFab& dst,
527 const MultiFab& src,
528 int srccomp,
529 int dstcomp,
530 int numcomp,
531 const IntVect& nghost);
532
538 static void Swap (MultiFab& dst,
539 MultiFab& src,
540 int srccomp,
541 int dstcomp,
542 int numcomp,
543 int nghost);
544
545 static void Swap (MultiFab& dst,
546 MultiFab& src,
547 int srccomp,
548 int dstcomp,
549 int numcomp,
550 const IntVect& nghost);
551
556 static void Subtract (MultiFab& dst,
557 const MultiFab& src,
558 int srccomp,
559 int dstcomp,
560 int numcomp,
561 int nghost);
562
563 static void Subtract (MultiFab& dst,
564 const MultiFab& src,
565 int srccomp,
566 int dstcomp,
567 int numcomp,
568 const IntVect& nghost);
573 static void Multiply (MultiFab& dst,
574 const MultiFab& src,
575 int srccomp,
576 int dstcomp,
577 int numcomp,
578 int nghost);
579
580 static void Multiply (MultiFab& dst,
581 const MultiFab& src,
582 int srccomp,
583 int dstcomp,
584 int numcomp,
585 const IntVect& nghost);
590 static void Divide (MultiFab& dst,
591 const MultiFab& src,
592 int srccomp,
593 int dstcomp,
594 int numcomp,
595 int nghost);
596
597 static void Divide (MultiFab& dst,
598 const MultiFab& src,
599 int srccomp,
600 int dstcomp,
601 int numcomp,
602 const IntVect& nghost);
606 static void Saxpy (MultiFab& dst,
607 Real a,
608 const MultiFab& src,
609 int srccomp,
610 int dstcomp,
611 int numcomp,
612 int nghost);
613
615
619 static void Xpay (MultiFab& dst,
620 Real a,
621 const MultiFab& src,
622 int srccomp,
623 int dstcomp,
624 int numcomp,
625 int nghost);
626
628
632 static void LinComb (MultiFab& dst,
633 Real a,
634 const MultiFab& x,
635 int xcomp,
636 Real b,
637 const MultiFab& y,
638 int ycomp,
639 int dstcomp,
640 int numcomp,
641 int nghost);
642
644
648 static void AddProduct (MultiFab& dst,
649 const MultiFab& src1,
650 int comp1,
651 const MultiFab& src2,
652 int comp2,
653 int dstcomp,
654 int numcomp,
655 int nghost);
656
657 static void AddProduct (MultiFab& dst,
658 const MultiFab& src1,
659 int comp1,
660 const MultiFab& src2,
661 int comp2,
662 int dstcomp,
663 int numcomp,
664 const IntVect& nghost);
665
673 [[nodiscard]] bool is_finite (bool local=false) const;
674 [[nodiscard]] bool is_finite (int scomp, int ncomp, int ngrow = 0, bool local=false) const;
675 [[nodiscard]] bool is_finite (int scomp, int ncomp, const IntVect& ngrow, bool local=false) const;
676
684 [[nodiscard]] bool contains_nan (bool local=false) const;
685 [[nodiscard]] bool contains_nan (int scomp, int ncomp, int ngrow = 0, bool local=false) const;
686 [[nodiscard]] bool contains_nan (int scomp, int ncomp, const IntVect& ngrow, bool local=false) const;
693 [[nodiscard]] bool contains_inf (bool local=false) const;
694 [[nodiscard]] bool contains_inf (int scomp, int ncomp, int ngrow = 0, bool local=false) const;
695 [[nodiscard]] bool contains_inf (int scomp, int ncomp, const IntVect& ngrow, bool local=false) const;
696
700 [[nodiscard]] std::unique_ptr<MultiFab> OverlapMask (const Periodicity& period = Periodicity::NonPeriodic()) const;
702 [[nodiscard]] std::unique_ptr<iMultiFab> OwnerMask (const Periodicity& period = Periodicity::NonPeriodic()) const;
703
705 void AverageSync (const Periodicity& period = Periodicity::NonPeriodic());
707 void WeightedSync (const MultiFab& wgt, const Periodicity& period = Periodicity::NonPeriodic());
709 void OverrideSync (const iMultiFab& msk, const Periodicity& period = Periodicity::NonPeriodic());
710
712
713 static void Initialize ();
714 static void Finalize ();
715
716private:
717 //
721
722 void initVal ();
723};
724
725#ifndef _MSC_VER
726inline void GccPlacaterMF ()
727{
728 std::allocator<MultiFab*> a1;
729 std::allocator<MultiFab const*> a2;
730 std::allocator<FabArray<FArrayBox>*> a3;
731 std::allocator<FabArray<FArrayBox> const*> a4;
732
737}
738#endif
739
740}
741
742#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:551
Calculates the distribution of FABs to MPI processes.
Definition AMReX_DistributionMapping.H:41
A Fortran Array of REALs.
Definition AMReX_FArrayBox.H:229
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:345
Definition AMReX_FabFactory.H:50
A collection (stored as an array) of FArrayBox objects.
Definition AMReX_MultiFab.H:38
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:719
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:720
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
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 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
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
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
IntVect maxIndex(int comp, int nghost=0) const
Definition AMReX_MultiFab.cpp:1021
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
Vector< Real > norminf(const Vector< int > &comps, int nghost=0, bool local=false, bool ignore_covered=false) const
Definition AMReX_MultiFab.H:203
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
Real norminf(int comp=0, int nghost=0, bool local=false, bool ignore_covered=false) const
Definition AMReX_MultiFab.H:184
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
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
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
IntVect minIndex(int comp, int nghost=0) const
Definition AMReX_MultiFab.cpp:1013
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
Real norminf(const iMultiFab &mask, int comp=0, int nghost=0, bool local=false) const
Definition AMReX_MultiFab.H:189
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())
Sync up nodal data with owners overriding non-owners.
Definition AMReX_MultiFab.cpp:1578
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
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
Definition AMReX_iMultiFab.H:32
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:458
void GccPlacaterMF()
Definition AMReX_MultiFab.H:726
FabArray memory allocation information.
Definition AMReX_FabArray.H:66