Block-Structured AMR Software Framework
AMReX_FluxRegister.H
Go to the documentation of this file.
1 #ifndef AMREX_FLUXREGISTER_H_
2 #define AMREX_FLUXREGISTER_H_
3 #include <AMReX_Config.H>
4 
5 #include <AMReX_BndryRegister.H>
6 #include <AMReX_Geometry.H>
7 #include <AMReX_Array.H>
8 
9 namespace amrex {
10 
11 
18  :
19  public BndryRegister
20 {
21 public:
22 
26  FluxRegister();
27 
37  FluxRegister (const BoxArray& fine_boxes,
38  const DistributionMapping& dm,
39  const IntVect& ref_ratio,
40  int fine_lev,
41  int nvar);
42 
46  ~FluxRegister () = default;
47 
48  FluxRegister (FluxRegister&& rhs) noexcept = default;
49 
50  FluxRegister (const FluxRegister& rhs) = delete;
51  FluxRegister& operator= (const FluxRegister& rhs) = delete;
53 
55  enum FrOp {COPY = 0, ADD = 1};
56 
67  void define (const BoxArray& fine_boxes,
68  const DistributionMapping& dm,
69  const IntVect& ref_ratio,
70  int fine_lev,
71  int nvar);
72 
73  void clear ();
74 
75 
79  const IntVect& refRatio () const noexcept;
80 
84  int fineLevel () const noexcept;
85 
89  int crseLevel () const noexcept;
90 
94  int nComp () const noexcept;
95 
99  const BoxArray& coarsenedBoxes () const noexcept;
100 
106  Real SumReg (int comp) const;
107 
120  void CrseInit (const MultiFab& mflx,
121  const MultiFab& area,
122  int dir,
123  int srccomp,
124  int destcomp,
125  int numcomp,
126  Real mult = -1.0,
127  FrOp op = FluxRegister::COPY);
128 
140  void CrseInit (const MultiFab& mflx,
141  int dir,
142  int srccomp,
143  int destcomp,
144  int numcomp,
145  Real mult = -1.0,
146  FrOp op = FluxRegister::COPY);
147 
163  void CrseAdd (const MultiFab& mflx,
164  const MultiFab& area,
165  int dir,
166  int srccomp,
167  int destcomp,
168  int numcomp,
169  Real mult,
170  const Geometry& geom);
171 
183  void CrseAdd (const MultiFab& mflx,
184  int dir,
185  int srccomp,
186  int destcomp,
187  int numcomp,
188  Real mult,
189  const Geometry& geom);
190 
203  void FineAdd (const MultiFab& mflx,
204  int dir,
205  int srccomp,
206  int destcomp,
207  int numcomp,
208  Real mult);
209 
221  void FineAdd (const MultiFab& mflx,
222  const MultiFab& area,
223  int dir,
224  int srccomp,
225  int destcomp,
226  int numcomp,
227  Real mult);
228 
242  void FineAdd (const FArrayBox& flux,
243  int dir,
244  int boxno,
245  int srccomp,
246  int destcomp,
247  int numcomp,
248  Real mult,
249  RunOn runon) noexcept;
250 
263  void FineAdd (const FArrayBox& flux,
264  const FArrayBox& area,
265  int dir,
266  int boxno,
267  int srccomp,
268  int destcomp,
269  int numcomp,
270  Real mult,
271  RunOn runon) noexcept;
272 
283  void FineSetVal (int dir,
284  int boxno,
285  int destcomp,
286  int numcomp,
287  Real val,
288  RunOn runon) noexcept;
289 
301  void Reflux (MultiFab& mf,
302  const MultiFab& volume,
303  Real scale,
304  int scomp,
305  int dcomp,
306  int nc,
307  const Geometry& crse_geom);
308 
309  void Reflux (MultiFab& mf,
310  const MultiFab& volume,
311  int dir,
312  Real scale,
313  int scomp,
314  int dcomp,
315  int nc,
316  const Geometry& crse_geom);
317 
328  void Reflux (MultiFab& mf,
329  Real scale,
330  int scomp,
331  int dcomp,
332  int nc,
333  const Geometry& crse_geom);
334 
335  void Reflux (MultiFab& mf,
336  int dir,
337  Real scale,
338  int scomp,
339  int dcomp,
340  int nc,
341  const Geometry& crse_geom);
342 
353  void OverwriteFlux (Array<MultiFab*,AMREX_SPACEDIM> const& crse_fluxes,
354  Real scale, int srccomp, int destcomp, int numcomp,
355  const Geometry& crse_geom);
356 
357 
363  void ClearInternalBorders (const Geometry& crse_geom);
364 
371  void write (const std::string& name, std::ostream& os) const;
372 
379  void read (const std::string& name, std::istream& is);
380 
381 // public for cuda
382 
383  void Reflux (MultiFab& mf, const MultiFab& volume, Orientation face,
384  Real scale, int scomp, int dcomp, int nc, const Geometry& geom);
385 
386 private:
387 
390 
393 
395  int ncomp;
396 };
397 
398 }
399 
400 #endif /*_FLUXREGISTER_H_*/
A BndryRegister organizes FabSets bounding each grid in a BoxArray. A FabSet is maintained for each b...
Definition: AMReX_BndryRegister.H:41
A collection of Boxes stored in an Array.
Definition: AMReX_BoxArray.H:530
Calculates the distribution of FABs to MPI processes.
Definition: AMReX_DistributionMapping.H:41
A Fortran Array of REALs.
Definition: AMReX_FArrayBox.H:229
Flux Register.
Definition: AMReX_FluxRegister.H:20
~FluxRegister()=default
The destructor.
FrOp
An enum that says whether to add or copy src data to members.
Definition: AMReX_FluxRegister.H:55
@ ADD
Definition: AMReX_FluxRegister.H:55
@ COPY
Definition: AMReX_FluxRegister.H:55
void write(const std::string &name, std::ostream &os) const
Write (used for writing to checkpoint)
Definition: AMReX_FluxRegister.cpp:881
void read(const std::string &name, std::istream &is)
Read (used for reading from checkpoint)
Definition: AMReX_FluxRegister.cpp:897
void Reflux(MultiFab &mf, const MultiFab &volume, Real scale, int scomp, int dcomp, int nc, const Geometry &crse_geom)
Apply flux correction. Note that this takes the coarse Geometry.
Definition: AMReX_FluxRegister.cpp:523
int ncomp
Number of state components.
Definition: AMReX_FluxRegister.H:395
void OverwriteFlux(Array< MultiFab *, AMREX_SPACEDIM > const &crse_fluxes, Real scale, int srccomp, int destcomp, int numcomp, const Geometry &crse_geom)
Overwrite the coarse flux at the coarse/fine interface (and the interface only) with the fine flux st...
Definition: AMReX_FluxRegister.cpp:732
int fine_level
Current level + 1.
Definition: AMReX_FluxRegister.H:392
FluxRegister(FluxRegister &&rhs) noexcept=default
void clear()
Definition: AMReX_FluxRegister.cpp:89
void CrseAdd(const MultiFab &mflx, const MultiFab &area, int dir, int srccomp, int destcomp, int numcomp, Real mult, const Geometry &geom)
Add coarse fluxes to the flux register. This is different from CrseInit with FluxRegister::ADD....
Definition: AMReX_FluxRegister.cpp:283
FluxRegister(const FluxRegister &rhs)=delete
void FineSetVal(int dir, int boxno, int destcomp, int numcomp, Real val, RunOn runon) noexcept
Set flux correction data for a fine box (given by boxno) to a given value. This routine used by FLASH...
Definition: AMReX_FluxRegister.cpp:496
int fineLevel() const noexcept
Returns the level number of the fine level.
Definition: AMReX_FluxRegister.cpp:34
IntVect ratio
Refinement ratio.
Definition: AMReX_FluxRegister.H:389
Real SumReg(int comp) const
Returns the sum of the registers.
Definition: AMReX_FluxRegister.cpp:95
FluxRegister()
The default constructor.
Definition: AMReX_FluxRegister.cpp:11
const BoxArray & coarsenedBoxes() const noexcept
The coarsened boxes.
Definition: AMReX_FluxRegister.cpp:52
void ClearInternalBorders(const Geometry &crse_geom)
Set internal borders to zero.
Definition: AMReX_FluxRegister.cpp:644
void FineAdd(const MultiFab &mflx, int dir, int srccomp, int destcomp, int numcomp, Real mult)
Increment flux correction with fine data.
Definition: AMReX_FluxRegister.cpp:359
void CrseInit(const MultiFab &mflx, const MultiFab &area, int dir, int srccomp, int destcomp, int numcomp, Real mult=-1.0, FrOp op=FluxRegister::COPY)
Initialize flux correction with coarse data.
Definition: AMReX_FluxRegister.cpp:160
const IntVect & refRatio() const noexcept
Returns the refinement ratio.
Definition: AMReX_FluxRegister.cpp:28
FluxRegister & operator=(const FluxRegister &rhs)=delete
void define(const BoxArray &fine_boxes, const DistributionMapping &dm, const IntVect &ref_ratio, int fine_lev, int nvar)
Initialize after using default constructor. This version allows setting the DistributionMapping.
Definition: AMReX_FluxRegister.cpp:58
int crseLevel() const noexcept
Returns the level number of the coarse level (fineLevel()-1).
Definition: AMReX_FluxRegister.cpp:40
int nComp() const noexcept
The number of components.
Definition: AMReX_FluxRegister.cpp:46
Rectangular problem domain geometry.
Definition: AMReX_Geometry.H:73
A collection (stored as an array) of FArrayBox objects.
Definition: AMReX_MultiFab.H:38
Encapsulation of the Orientation of the Faces of a Box.
Definition: AMReX_Orientation.H:29
Definition: AMReX_Amr.cpp:49
RunOn
Definition: AMReX_GpuControl.H:69
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE IntVectND< dim > scale(const IntVectND< dim > &p, int s) noexcept
Returns a IntVectND obtained by multiplying each of the components of this IntVectND by s.
Definition: AMReX_IntVect.H:1006
std::array< T, N > Array
Definition: AMReX_Array.H:23