Block-Structured AMR Software Framework
Loading...
Searching...
No Matches
amrex::YAFluxRegisterT< MF > Class Template Reference

#include <AMReX_YAFluxRegister.H>

Inheritance diagram for amrex::YAFluxRegisterT< MF >:
amrex::EBFluxRegister

Public Types

enum  CellType : int { crse_cell = 0 , crse_fine_boundary_cell , fine_cell }
 
using T = typename MF::value_type
 
using FAB = typename MF::fab_type
 

Public Member Functions

 YAFluxRegisterT ()=default
 Construct an empty register; call define() before use.
 
 YAFluxRegisterT (const BoxArray &fba, const BoxArray &cba, const DistributionMapping &fdm, const DistributionMapping &cdm, const Geometry &fgeom, const Geometry &cgeom, const IntVect &ref_ratio, int fine_lev, int nvar)
 Fully construct the register for a fine level.
 
void define (const BoxArray &fba, const BoxArray &cba, const DistributionMapping &fdm, const DistributionMapping &cdm, const Geometry &fgeom, const Geometry &cgeom, const IntVect &ref_ratio, int fine_lev, int nvar)
 Define the register using the same arguments as the constructor.
 
void reset ()
 Reset stored fluxes and flags before a coarse advance.
 
void CrseAdd (const MFIter &mfi, const std::array< FAB const *, 3 > &flux, const Real *dx, Real dt, RunOn runon) noexcept
 Add coarse-level fluxes for the tile indicated by mfi.
 
void CrseAdd (const MFIter &mfi, const std::array< FAB const *, 3 > &flux, const Real *dx, Real dt, int srccomp, int destcomp, int numcomp, RunOn runon) noexcept
 Add coarse-level fluxes with component slicing.
 
void FineAdd (const MFIter &mfi, const std::array< FAB const *, 3 > &flux, const Real *dx, Real dt, RunOn runon) noexcept
 Add fine-level fluxes for the tile identified by mfi.
 
void FineAdd (const MFIter &mfi, const std::array< FAB const *, 3 > &a_flux, const Real *dx, Real dt, int srccomp, int destcomp, int numcomp, RunOn runon) noexcept
 Add fine-level fluxes with component slicing.
 
void Reflux (MF &state, int dc=0)
 Apply the accumulated flux divergence to the coarse MultiFab state.
 
void Reflux (MF &state, int srccomp, int destcomp, int numcomp)
 Component-aware version of Reflux writing components [destcomp,destcomp+numcomp).
 
bool CrseHasWork (const MFIter &mfi) const noexcept
 Return true if coarse tile mfi abuts a coarse/fine boundary.
 
bool FineHasWork (const MFIter &mfi) const noexcept
 Return true if fine tile mfi has flux registers to update.
 
MF & getFineData ()
 Mutable access to fine-side accumulation data.
 
MF & getCrseData ()
 Mutable access to coarse-side accumulation data.
 
void setCrseVolume (MF const *cvol)
 
void setDeterministic (bool flag)
 Enable deterministic mode for GPU operations via flag (uses deterministic reductions).
 
bool getDeterministic () const
 

Protected Attributes

MF m_crse_data
 
iMultiFab m_crse_flag
 
Vector< intm_crse_fab_flag
 
MF m_cfpatch
 This is built on crse/fine patches.
 
MF m_cfp_mask
 
Vector< Vector< FAB * > > m_cfp_fab
 The size of this is (# of local fine grids (# of crse/fine patches for that grid))
 
Vector< intm_cfp_localindex
 
Geometry m_fine_geom
 
Geometry m_crse_geom
 
IntVect m_ratio
 
int m_fine_level
 
int m_ncomp
 
MF const * m_cvol = nullptr
 
bool m_deterministic = false
 

Detailed Description

template<typename MF>
class amrex::YAFluxRegisterT< MF >

YAFluxRegister is yet another Flux Register for refluxing.

At the beginning of a coarse step, reset() is called. In MFIter for the coarse level advance, CrseAdd is called with coarse flux. The flux is not scaled. In MFIter for the fine level advance, FineAdd is called. After the fine level finished its time steps, Reflux is called to update the coarse cells next to the coarse/fine boundary.

Member Typedef Documentation

◆ FAB

template<typename MF >
using amrex::YAFluxRegisterT< MF >::FAB = typename MF::fab_type

◆ T

template<typename MF >
using amrex::YAFluxRegisterT< MF >::T = typename MF::value_type

Member Enumeration Documentation

◆ CellType

template<typename MF >
enum amrex::YAFluxRegisterT::CellType : int
Enumerator
crse_cell 
crse_fine_boundary_cell 
fine_cell 

Constructor & Destructor Documentation

◆ YAFluxRegisterT() [1/2]

template<typename MF >
amrex::YAFluxRegisterT< MF >::YAFluxRegisterT ( )
default

Construct an empty register; call define() before use.

◆ YAFluxRegisterT() [2/2]

template<typename MF >
amrex::YAFluxRegisterT< MF >::YAFluxRegisterT ( const BoxArray fba,
const BoxArray cba,
const DistributionMapping fdm,
const DistributionMapping cdm,
const Geometry fgeom,
const Geometry cgeom,
const IntVect ref_ratio,
int  fine_lev,
int  nvar 
)

Fully construct the register for a fine level.

Parameters
fbaFine BoxArray.
cbaCoarse BoxArray.
fdmDistribution mapping for fine grids.
cdmDistribution mapping for coarse grids.
fgeomFine Geometry.
cgeomCoarse Geometry.
ref_ratioRefinement ratio between fine and coarse.
fine_levFine AMR level index.
nvarNumber of components.

Member Function Documentation

◆ CrseAdd() [1/2]

template<typename MF >
void amrex::YAFluxRegisterT< MF >::CrseAdd ( const MFIter mfi,
const std::array< FAB const *, 3 > &  flux,
const Real dx,
Real  dt,
int  srccomp,
int  destcomp,
int  numcomp,
RunOn  runon 
)
noexcept

Add coarse-level fluxes with component slicing.

Parameters
mfiIterator describing the coarse tile.
fluxArray of per-direction flux FABs.
dxCell spacing array.
dtCoarse timestep.
srccompStarting component in the fluxes.
destcompDestination component in the register.
numcompNumber of components.
runonExecution context (host/device).

◆ CrseAdd() [2/2]

template<typename MF >
void amrex::YAFluxRegisterT< MF >::CrseAdd ( const MFIter mfi,
const std::array< FAB const *, 3 > &  flux,
const Real dx,
Real  dt,
RunOn  runon 
)
noexcept

Add coarse-level fluxes for the tile indicated by mfi.

Parameters
mfiIterator describing the coarse tile.
fluxArray of per-direction flux FABs.
dxCell spacing array.
dtCoarse timestep (with sign).
runonExecution context (host/device).

◆ CrseHasWork()

template<typename MF >
bool amrex::YAFluxRegisterT< MF >::CrseHasWork ( const MFIter mfi) const
inlinenoexcept

Return true if coarse tile mfi abuts a coarse/fine boundary.

◆ define()

template<typename MF >
void amrex::YAFluxRegisterT< MF >::define ( const BoxArray fba,
const BoxArray cba,
const DistributionMapping fdm,
const DistributionMapping cdm,
const Geometry fgeom,
const Geometry cgeom,
const IntVect ref_ratio,
int  fine_lev,
int  nvar 
)

Define the register using the same arguments as the constructor.

◆ FineAdd() [1/2]

template<typename MF >
void amrex::YAFluxRegisterT< MF >::FineAdd ( const MFIter mfi,
const std::array< FAB const *, 3 > &  a_flux,
const Real dx,
Real  dt,
int  srccomp,
int  destcomp,
int  numcomp,
RunOn  runon 
)
noexcept

Add fine-level fluxes with component slicing.

Parameters
mfiIterator describing the fine tile.
a_fluxArray of per-direction flux FABs.
dxFine cell spacing array.
dtFine timestep.
srccompStarting component in the fluxes.
destcompDestination component in the register.
numcompNumber of components.
runonExecution context (host/device).

◆ FineAdd() [2/2]

template<typename MF >
void amrex::YAFluxRegisterT< MF >::FineAdd ( const MFIter mfi,
const std::array< FAB const *, 3 > &  flux,
const Real dx,
Real  dt,
RunOn  runon 
)
noexcept

Add fine-level fluxes for the tile identified by mfi.

Parameters
mfiIterator describing the fine tile.
fluxArray of per-direction flux FABs.
dxFine cell spacing array.
dtFine timestep (per substep).
runonExecution context (host/device).

◆ FineHasWork()

template<typename MF >
bool amrex::YAFluxRegisterT< MF >::FineHasWork ( const MFIter mfi) const
inlinenoexcept

Return true if fine tile mfi has flux registers to update.

◆ getCrseData()

template<typename MF >
MF & amrex::YAFluxRegisterT< MF >::getCrseData ( )

Mutable access to coarse-side accumulation data.

◆ getDeterministic()

template<typename MF >
bool amrex::YAFluxRegisterT< MF >::getDeterministic ( ) const
inline

◆ getFineData()

template<typename MF >
MF & amrex::YAFluxRegisterT< MF >::getFineData ( )

Mutable access to fine-side accumulation data.

◆ Reflux() [1/2]

template<typename MF >
void amrex::YAFluxRegisterT< MF >::Reflux ( MF &  state,
int  dc = 0 
)

Apply the accumulated flux divergence to the coarse MultiFab state.

Parameters
stateCoarse solution MultiFab to update.
dcDestination component offset.

◆ Reflux() [2/2]

template<typename MF >
void amrex::YAFluxRegisterT< MF >::Reflux ( MF &  state,
int  srccomp,
int  destcomp,
int  numcomp 
)

Component-aware version of Reflux writing components [destcomp,destcomp+numcomp).

◆ reset()

template<typename MF >
void amrex::YAFluxRegisterT< MF >::reset ( )

Reset stored fluxes and flags before a coarse advance.

◆ setCrseVolume()

template<typename MF >
void amrex::YAFluxRegisterT< MF >::setCrseVolume ( MF const *  cvol)
inline

For curvilinear coordinates only, where the supplied flux has already been multiplied by area. Provide the coarse-volume MultiFab cvol and keep it alive for the register lifetime (YAFluxRegister does not copy it).

◆ setDeterministic()

template<typename MF >
void amrex::YAFluxRegisterT< MF >::setDeterministic ( bool  flag)
inline

Enable deterministic mode for GPU operations via flag (uses deterministic reductions).

Member Data Documentation

◆ m_cfp_fab

template<typename MF >
Vector<Vector<FAB*> > amrex::YAFluxRegisterT< MF >::m_cfp_fab
protected

The size of this is (# of local fine grids (# of crse/fine patches for that grid))

◆ m_cfp_localindex

template<typename MF >
Vector<int> amrex::YAFluxRegisterT< MF >::m_cfp_localindex
protected

◆ m_cfp_mask

template<typename MF >
MF amrex::YAFluxRegisterT< MF >::m_cfp_mask
protected

◆ m_cfpatch

template<typename MF >
MF amrex::YAFluxRegisterT< MF >::m_cfpatch
protected

This is built on crse/fine patches.

◆ m_crse_data

template<typename MF >
MF amrex::YAFluxRegisterT< MF >::m_crse_data
protected

◆ m_crse_fab_flag

template<typename MF >
Vector<int> amrex::YAFluxRegisterT< MF >::m_crse_fab_flag
protected

◆ m_crse_flag

template<typename MF >
iMultiFab amrex::YAFluxRegisterT< MF >::m_crse_flag
protected

◆ m_crse_geom

template<typename MF >
Geometry amrex::YAFluxRegisterT< MF >::m_crse_geom
protected

◆ m_cvol

template<typename MF >
MF const* amrex::YAFluxRegisterT< MF >::m_cvol = nullptr
protected

◆ m_deterministic

template<typename MF >
bool amrex::YAFluxRegisterT< MF >::m_deterministic = false
protected

◆ m_fine_geom

template<typename MF >
Geometry amrex::YAFluxRegisterT< MF >::m_fine_geom
protected

◆ m_fine_level

template<typename MF >
int amrex::YAFluxRegisterT< MF >::m_fine_level
protected

◆ m_ncomp

template<typename MF >
int amrex::YAFluxRegisterT< MF >::m_ncomp
protected

◆ m_ratio

template<typename MF >
IntVect amrex::YAFluxRegisterT< MF >::m_ratio
protected

The documentation for this class was generated from the following files: