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

A FabSet is a group of FArrayBox's. The grouping is designed specifically to represent regions along the boundary of Box's, and are used to implement boundary conditions to discretized partial differential equations. More...

#include <AMReX_FabSet.H>

Public Types

using value_type = typename FabDataType< MF >::value_type
 
using FAB = typename FabDataType< MF >::fab_type
 

Public Member Functions

 FabSetT () noexcept=default
 The default constructor – you must later call define().
 
 FabSetT (const BoxArray &grids, const DistributionMapping &dmap, int ncomp)
 Construct a FabSetT<MF> covering grids with mapping dmap and ncomp components.
 
 ~FabSetT ()=default
 
 FabSetT (FabSetT< MF > &&rhs) noexcept=default
 
 FabSetT (const FabSetT< MF > &rhs)=delete
 
FabSetT< MF > & operator= (const FabSetT< MF > &rhs)=delete
 
FabSetT< MF > & operator= (FabSetT< MF > &&rhs)=delete
 
void define (const BoxArray &grids, const DistributionMapping &dmap, int ncomp)
 Define a FabSetT<MF> constructed via default constructor on grids/dmap with ncomp components.
 
FAB const & operator[] (const MFIter &mfi) const noexcept
 Access the FAB referenced by iterator mfi (const).
 
FABoperator[] (const MFIter &mfi) noexcept
 Access the FAB referenced by iterator mfi (mutable).
 
FAB const & operator[] (int i) const noexcept
 Access the FAB with index i (const).
 
FABoperator[] (int i) noexcept
 Access the FAB with index i (mutable).
 
Array4< value_type const > array (const MFIter &mfi) const noexcept
 Return an Array4 view for iterator mfi (const).
 
Array4< value_typearray (const MFIter &mfi) noexcept
 Return an Array4 view for iterator mfi (mutable).
 
Array4< value_type const > array (int i) const noexcept
 Return an Array4 view for FAB index i (const).
 
Array4< value_typearray (int i) noexcept
 Return an Array4 view for FAB index i (mutable).
 
Array4< value_type const > const_array (const MFIter &mfi) const noexcept
 Return a const Array4 view for iterator mfi.
 
Array4< value_type const > const_array (int i) const noexcept
 Return a const Array4 view for FAB index i.
 
MultiArray4< value_type const > arrays () const noexcept
 Return multi-fab const Array4 views to every FAB.
 
MultiArray4< value_typearrays () noexcept
 Return multi-fab mutable Array4 views to every FAB.
 
MultiArray4< value_type const > const_arrays () const noexcept
 Return multi-fab const Array4 views (alias of arrays()).
 
Box fabbox (int K) const noexcept
 Get the bounding Box of FAB index K.
 
int size () const noexcept
 Return number of FABs stored in this set.
 
const BoxArrayboxArray () const noexcept
 Return the BoxArray defining FAB locations.
 
const DistributionMappingDistributionMap () const noexcept
 Return the DistributionMapping describing FAB ownership.
 
MF & multiFab () noexcept
 Mutable access to the underlying MultiFab.
 
MF const & multiFab () const noexcept
 Const access to the underlying MultiFab.
 
int nComp () const noexcept
 Number of components carried by each FAB.
 
void clear ()
 Release storage and reset the internal MultiFab.
 
FabSetT< MF > & copyFrom (const FabSetT< MF > &src, int scomp, int dcomp, int ncomp)
 Copy components from another FabSet with identical layout.
 
FabSetT< MF > & copyFrom (const MF &src, int ngrow, int scomp, int dcomp, int ncomp, const Periodicity &period=Periodicity::NonPeriodic())
 Copy boundary data from a MultiFab intersecting the same boxes.
 
FabSetT< MF > & plusFrom (const FabSetT< MF > &src, int scomp, int dcomp, int ncomp)
 Accumulate data from another FabSet with identical layout.
 
FabSetT< MF > & plusFrom (const MF &src, int ngrow, int scomp, int dcomp, int ncomp, const Periodicity &period=Periodicity::NonPeriodic(), bool deterministic=false)
 Accumulate data from a MultiFab into this boundary set.
 
void copyTo (MF &dest, int ngrow, int scomp, int dcomp, int ncomp, const Periodicity &period=Periodicity::NonPeriodic()) const
 Copy data from this boundary set back into a MultiFab.
 
void plusTo (MF &dest, int ngrow, int scomp, int dcomp, int ncomp, const Periodicity &period=Periodicity::NonPeriodic()) const
 Add this boundary data back into a MultiFab.
 
void setVal (value_type val)
 Fill every component of every FAB with scalar value val.
 
void setVal (value_type val, int comp, int num_comp)
 Fill num_comp components starting at comp with scalar value val.
 
FabSetT< MF > & linComb (value_type a, value_type b, const FabSetT< MF > &src, int scomp, int dcomp, int ncomp)
 Linear combination: this := a*this + b*src.
 
FabSetT< MF > & linComb (value_type a, const MF &mfa, int a_comp, value_type b, const MF &mfb, int b_comp, int dcomp, int ncomp, int ngrow)
 Boundary linear combination from two MultiFabs: this := a*mfa + b*mfb.
 
void write (const std::string &name) const
 Write (used for writing to checkpoint).
 
void read (const std::string &name)
 Read (used for reading from checkpoint).
 

Static Public Member Functions

static void Copy (FabSetT< MF > &dst, const FabSetT< MF > &src)
 Local copy function.
 

Friends

class FabSetIter
 
class FluxRegister
 

Detailed Description

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

A FabSet is a group of FArrayBox's. The grouping is designed specifically to represent regions along the boundary of Box's, and are used to implement boundary conditions to discretized partial differential equations.

A FabSet is an array of pointers to FABs. The standard FAB operators, however, have been modified to be more useful for maintaining boundary conditions for partial differential equations discretized on boxes. Under normal circumstances, a FAB will be created for each face of a box. For a group of boxes, a FabSet will be the group of FABs at a particular orientation (ie. the lo-i side of each grid in a list).

Since a FabSet FAB will likely be used to bound a grid box, FArrayBox::resize() operations are disallowed. Also, to preserve flexibility in applicable boundary scenarios, intersecting FABs in the FabSet are not guaranteed to contain identical data–thus copy operations from a FabSet to any FAB-like structure may be order-dependent.

FabSets are used primarily as a data storage mechanism, and are manipulated by more sophisticated control classes.

Member Typedef Documentation

◆ FAB

template<typename MF >
using amrex::FabSetT< MF >::FAB = typename FabDataType<MF>::fab_type

◆ value_type

template<typename MF >
using amrex::FabSetT< MF >::value_type = typename FabDataType<MF>::value_type

Constructor & Destructor Documentation

◆ FabSetT() [1/4]

template<typename MF >
amrex::FabSetT< MF >::FabSetT ( )
defaultnoexcept

The default constructor – you must later call define().

◆ FabSetT() [2/4]

template<typename MF >
amrex::FabSetT< MF >::FabSetT ( const BoxArray grids,
const DistributionMapping dmap,
int  ncomp 
)

Construct a FabSetT<MF> covering grids with mapping dmap and ncomp components.

◆ ~FabSetT()

template<typename MF >
amrex::FabSetT< MF >::~FabSetT ( )
default

◆ FabSetT() [3/4]

template<typename MF >
amrex::FabSetT< MF >::FabSetT ( FabSetT< MF > &&  rhs)
defaultnoexcept

◆ FabSetT() [4/4]

template<typename MF >
amrex::FabSetT< MF >::FabSetT ( const FabSetT< MF > &  rhs)
delete

Member Function Documentation

◆ array() [1/4]

template<typename MF >
Array4< value_type const > amrex::FabSetT< MF >::array ( const MFIter mfi) const
inlinenoexcept

Return an Array4 view for iterator mfi (const).

◆ array() [2/4]

template<typename MF >
Array4< value_type > amrex::FabSetT< MF >::array ( const MFIter mfi)
inlinenoexcept

Return an Array4 view for iterator mfi (mutable).

◆ array() [3/4]

template<typename MF >
Array4< value_type const > amrex::FabSetT< MF >::array ( int  i) const
inlinenoexcept

Return an Array4 view for FAB index i (const).

◆ array() [4/4]

template<typename MF >
Array4< value_type > amrex::FabSetT< MF >::array ( int  i)
inlinenoexcept

Return an Array4 view for FAB index i (mutable).

◆ arrays() [1/2]

template<typename MF >
MultiArray4< value_type const > amrex::FabSetT< MF >::arrays ( ) const
inlinenoexcept

Return multi-fab const Array4 views to every FAB.

◆ arrays() [2/2]

template<typename MF >
MultiArray4< value_type > amrex::FabSetT< MF >::arrays ( )
inlinenoexcept

Return multi-fab mutable Array4 views to every FAB.

◆ boxArray()

template<typename MF >
const BoxArray & amrex::FabSetT< MF >::boxArray ( ) const
inlinenoexcept

Return the BoxArray defining FAB locations.

◆ clear()

template<typename MF >
void amrex::FabSetT< MF >::clear ( )
inline

Release storage and reset the internal MultiFab.

◆ const_array() [1/2]

template<typename MF >
Array4< value_type const > amrex::FabSetT< MF >::const_array ( const MFIter mfi) const
inlinenoexcept

Return a const Array4 view for iterator mfi.

◆ const_array() [2/2]

template<typename MF >
Array4< value_type const > amrex::FabSetT< MF >::const_array ( int  i) const
inlinenoexcept

Return a const Array4 view for FAB index i.

◆ const_arrays()

template<typename MF >
MultiArray4< value_type const > amrex::FabSetT< MF >::const_arrays ( ) const
inlinenoexcept

Return multi-fab const Array4 views (alias of arrays()).

◆ Copy()

template<typename MF >
void amrex::FabSetT< MF >::Copy ( FabSetT< MF > &  dst,
const FabSetT< MF > &  src 
)
static

Local copy function.

Parameters
dstDestination FabSet.
srcSource FabSet.

◆ copyFrom() [1/2]

template<typename MF >
FabSetT< MF > & amrex::FabSetT< MF >::copyFrom ( const FabSetT< MF > &  src,
int  scomp,
int  dcomp,
int  ncomp 
)

Copy components from another FabSet with identical layout.

Parameters
srcSource FabSet.
scompStarting component in src.
dcompDestination component index.
ncompNumber of components to copy.
Returns
Reference to this set.

◆ copyFrom() [2/2]

template<typename MF >
FabSetT< MF > & amrex::FabSetT< MF >::copyFrom ( const MF &  src,
int  ngrow,
int  scomp,
int  dcomp,
int  ncomp,
const Periodicity period = Periodicity::NonPeriodic() 
)

Copy boundary data from a MultiFab intersecting the same boxes.

Parameters
srcSource MultiFab.
ngrowNumber of ghost cells pulled from src.
scompSource component index.
dcompDestination component index.
ncompNumber of components.
periodPeriodicity description.

◆ copyTo()

template<typename MF >
void amrex::FabSetT< MF >::copyTo ( MF &  dest,
int  ngrow,
int  scomp,
int  dcomp,
int  ncomp,
const Periodicity period = Periodicity::NonPeriodic() 
) const

Copy data from this boundary set back into a MultiFab.

Parameters
destDestination MultiFab.
ngrowNumber of grow cells on dest to fill.
scompSource component index inside this set.
dcompDestination component index.
ncompNumber of components.
periodPeriodicity description.

◆ define()

template<typename MF >
void amrex::FabSetT< MF >::define ( const BoxArray grids,
const DistributionMapping dmap,
int  ncomp 
)

Define a FabSetT<MF> constructed via default constructor on grids/dmap with ncomp components.

◆ DistributionMap()

template<typename MF >
const DistributionMapping & amrex::FabSetT< MF >::DistributionMap ( ) const
inlinenoexcept

Return the DistributionMapping describing FAB ownership.

◆ fabbox()

template<typename MF >
Box amrex::FabSetT< MF >::fabbox ( int  K) const
inlinenoexcept

Get the bounding Box of FAB index K.

◆ linComb() [1/2]

template<typename MF >
FabSetT< MF > & amrex::FabSetT< MF >::linComb ( value_type  a,
const MF &  mfa,
int  a_comp,
value_type  b,
const MF &  mfb,
int  b_comp,
int  dcomp,
int  ncomp,
int  ngrow 
)

Boundary linear combination from two MultiFabs: this := a*mfa + b*mfb.

Parameters
aScale applied to mfa.
mfaFirst MultiFab input.
a_compComponent offset in mfa.
bScale applied to mfb.
mfbSecond MultiFab input.
b_compComponent offset in mfb.
dcompDestination component offset.
ncompNumber of components.
ngrowNumber of grow cells required in inputs.

◆ linComb() [2/2]

template<typename MF >
FabSetT< MF > & amrex::FabSetT< MF >::linComb ( value_type  a,
value_type  b,
const FabSetT< MF > &  src,
int  scomp,
int  dcomp,
int  ncomp 
)

Linear combination: this := a*this + b*src.

Parameters
aScaling applied to the current FabSet.
bScaling applied to src.
srcSource FabSet (must match layout).
scompSource component offset.
dcompDestination component offset.
ncompNumber of components.

◆ multiFab() [1/2]

template<typename MF >
MF const & amrex::FabSetT< MF >::multiFab ( ) const
inlinenoexcept

Const access to the underlying MultiFab.

◆ multiFab() [2/2]

template<typename MF >
MF & amrex::FabSetT< MF >::multiFab ( )
inlinenoexcept

Mutable access to the underlying MultiFab.

◆ nComp()

template<typename MF >
int amrex::FabSetT< MF >::nComp ( ) const
inlinenoexcept

Number of components carried by each FAB.

◆ operator=() [1/2]

template<typename MF >
FabSetT< MF > & amrex::FabSetT< MF >::operator= ( const FabSetT< MF > &  rhs)
delete

◆ operator=() [2/2]

template<typename MF >
FabSetT< MF > & amrex::FabSetT< MF >::operator= ( FabSetT< MF > &&  rhs)
delete

◆ operator[]() [1/4]

template<typename MF >
FAB const & amrex::FabSetT< MF >::operator[] ( const MFIter mfi) const
inlinenoexcept

Access the FAB referenced by iterator mfi (const).

◆ operator[]() [2/4]

template<typename MF >
FAB & amrex::FabSetT< MF >::operator[] ( const MFIter mfi)
inlinenoexcept

Access the FAB referenced by iterator mfi (mutable).

◆ operator[]() [3/4]

template<typename MF >
FAB const & amrex::FabSetT< MF >::operator[] ( int  i) const
inlinenoexcept

Access the FAB with index i (const).

◆ operator[]() [4/4]

template<typename MF >
FAB & amrex::FabSetT< MF >::operator[] ( int  i)
inlinenoexcept

Access the FAB with index i (mutable).

◆ plusFrom() [1/2]

template<typename MF >
FabSetT< MF > & amrex::FabSetT< MF >::plusFrom ( const FabSetT< MF > &  src,
int  scomp,
int  dcomp,
int  ncomp 
)

Accumulate data from another FabSet with identical layout.

Parameters
srcSource FabSet.
scompSource component index.
dcompDestination component index.
ncompNumber of components to accumulate.

◆ plusFrom() [2/2]

template<typename MF >
FabSetT< MF > & amrex::FabSetT< MF >::plusFrom ( const MF &  src,
int  ngrow,
int  scomp,
int  dcomp,
int  ncomp,
const Periodicity period = Periodicity::NonPeriodic(),
bool  deterministic = false 
)

Accumulate data from a MultiFab into this boundary set.

Parameters
srcSource MultiFab.
ngrowNumber of ghost cells read from src.
scompSource component index.
dcompDestination component index.
ncompNumber of components.
periodPeriodicity description.
deterministicWhether GPU operations must be deterministic.

◆ plusTo()

template<typename MF >
void amrex::FabSetT< MF >::plusTo ( MF &  dest,
int  ngrow,
int  scomp,
int  dcomp,
int  ncomp,
const Periodicity period = Periodicity::NonPeriodic() 
) const

Add this boundary data back into a MultiFab.

Parameters
destDestination MultiFab.
ngrowNumber of grow cells to cover on dest.
scompSource component index inside this set.
dcompDestination component index.
ncompNumber of components.
periodPeriodicity description.

◆ read()

template<typename MF >
void amrex::FabSetT< MF >::read ( const std::string &  name)

Read (used for reading from checkpoint).

Parameters
nameBase filename previously written by write().

◆ setVal() [1/2]

template<typename MF >
void amrex::FabSetT< MF >::setVal ( value_type  val)

Fill every component of every FAB with scalar value val.

◆ setVal() [2/2]

template<typename MF >
void amrex::FabSetT< MF >::setVal ( value_type  val,
int  comp,
int  num_comp 
)

Fill num_comp components starting at comp with scalar value val.

◆ size()

template<typename MF >
int amrex::FabSetT< MF >::size ( ) const
inlinenoexcept

Return number of FABs stored in this set.

◆ write()

template<typename MF >
void amrex::FabSetT< MF >::write ( const std::string &  name) const

Write (used for writing to checkpoint).

Parameters
nameBase filename to write.

Friends And Related Symbol Documentation

◆ FabSetIter

template<typename MF >
friend class FabSetIter
friend

◆ FluxRegister

template<typename MF >
friend class FluxRegister
friend

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