Block-Structured AMR Software Framework
amrex::NonLocalBC Namespace Reference

Namespaces

 detail
 

Classes

struct  HasInverseMemFn
 Type trait that tests if T has an InverseImage class member function. More...
 
struct  IsIndexMapping
 Tests if a given type IndexMap is usable as an index mapping between two index based coordinate systems. More...
 
struct  MultiBlockIndexMapping
 This struct describes an affine index transformation for two coordinate systems. More...
 
struct  MultiBlockCommMetaData
 This is the index mapping based on the DTOS MultiBlockDestToSrc. More...
 
struct  IsFabProjection
 This type trait tests if a type P is a projection for FAB. More...
 
struct  Identity
 This class acts as a default no-op operator. More...
 
struct  MapComponents
 This class takes a projection and a component map and combines them to form a new projection. More...
 
struct  SwapComponents
 This is a permutation where only two components are swapped. More...
 
struct  SwapComponents< I, -1 >
 
struct  SwapComponents<-1, J >
 
struct  SwapComponents<-1, -1 >
 
struct  CommData
 This class holds data buffers for either immediate MPI send or recv calls. More...
 
struct  CommHandler
 This class stores both recv and send buffers with an associated MPI tag. More...
 
struct  IsDataPacking
 This type trait tests if a given type DP satisfies the DataPacking concept for type FAB. More...
 
struct  PackComponents
 Contains information about which components take part of the data transaction. More...
 
struct  ApplyDtosAndProjectionOnReciever
 This class specializes behaviour on local copies and unpacking receive buffers. More...
 
struct  NoLocalCopy
 
struct  DoLocalCopy
 

Typedefs

template<typename T , typename... Args>
using Inverse_t = decltype(std::declval< T >().Inverse(std::declval< Args >()...))
 Return type of an InverseImage class member function. More...
 
using DynamicSwapComponents = SwapComponents<-1, -1 >
 
template<typename... Args>
using PrepareSendBuffers_t = decltype(PrepareSendBuffers(std::declval< Args >()...))
 
template<typename... Args>
using PrepareRecvBuffers_t = decltype(PrepareRecvBuffers(std::declval< Args >()...))
 
template<typename... Args>
using PackSendBuffers_t = decltype(PackSendBuffers(std::declval< Args >()...))
 
template<typename... Args>
using UnpackRecvBuffers_t = decltype(UnpackRecvBuffers(std::declval< Args >()...))
 
template<typename... Args>
using LocalCopy_t = decltype(LocalCopy(std::declval< Args >()...))
 

Functions

void PrepareCommBuffers (CommData &comm, const FabArrayBase::MapOfCopyComTagContainers &cctc, int n_components, std::size_t object_size, std::size_t align)
 Fill all class member variables of comm but the request and the stats vector. More...
 
void PostRecvs (CommData &recv, int mpi_tag)
 Initiate all recvieves with MPI_Irecv calls associated with tag mpi_tag. More...
 
void PostSends (CommData &send, int mpi_tag)
 Initiate all sends with MPI_Isend calls associated with tag mpi_tag. More...
 
template MultiBlockCommMetaData ParallelCopy (FabArray< FArrayBox > &dest, const Box &destbox, const FabArray< FArrayBox > &src, int destcomp, int srccomp, int numcomp, const IntVect &ngrow, MultiBlockIndexMapping const &, Identity const &)
 
template<typename DTOS >
std::enable_if_t< IsCallableR< Dim3, DTOS, Dim3 >::value, IntVectApply (DTOS const &dtos, const IntVect &iv)
 Applies the Dim3 to Dim3 mapping onto IntVects. More...
 
template<typename DTOS >
std::enable_if_t< IsCallableR< Dim3, DTOS, Dim3 >::value &&!IsCallableR< IndexType, DTOS, IndexType >::value, BoxImage (DTOS const &dtos, const Box &box)
 Applies the Dim3 to Dim3 mapping onto Boxes but does not change the index type. More...
 
template<typename DTOS >
std::enable_if_t< IsCallableR< Dim3, DTOS, Dim3 >::value &&IsCallableR< IndexType, DTOS, IndexType >::value, BoxImage (DTOS const &dtos, const Box &box)
 Applies the Dim3 to Dim3 mapping onto Boxes and maps the index type. More...
 
template<typename DTOS >
std::enable_if_t< HasInverseMemFn< DTOS >::value, IntVectApplyInverse (DTOS const &dtos, const IntVect &iv)
 Applies the Dim3 to Dim3 invserse mapping onto IntVects. More...
 
template<typename DTOS >
std::enable_if_t< HasInverseMemFn< DTOS >::value &&!IsCallableR< IndexType, DTOS, IndexType >::value, BoxInverseImage (DTOS const &dtos, const Box &box)
 Applies the inverse Dim3 to Dim3 mapping onto Boxes without changing the index type. More...
 
template<typename DTOS >
std::enable_if_t< HasInverseMemFn< DTOS >::value &&IsCallableR< IndexType, DTOS, IndexType >::value, BoxInverseImage (DTOS const &dtos, const Box &box)
 Applies the inverse Dim3 to Dim3 mapping onto Boxes. More...
 
template<class FAB , class DTOS = Identity, class Proj = Identity>
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >)> local_copy_cpu (FabArray< FAB > &dest, const FabArray< FAB > &src, int dcomp, int scomp, int ncomp, FabArrayBase::CopyComTagsContainer const &local_tags, DTOS const &dtos=DTOS{}, Proj const &proj=Proj{}) noexcept
 
template<class FAB , class DTOS = Identity, class Proj = Identity>
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >)> unpack_recv_buffer_cpu (FabArray< FAB > &mf, int dcomp, int ncomp, Vector< char * > const &recv_data, Vector< std::size_t > const &recv_size, Vector< FabArrayBase::CopyComTagsContainer const * > const &recv_cctc, DTOS const &dtos=DTOS{}, Proj const &proj=Proj{}) noexcept
 
template<class FAB , class DTOS = Identity, class Proj = Identity>
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >)> local_copy_gpu (FabArray< FAB > &dest, const FabArray< FAB > &src, int dcomp, int scomp, int ncomp, FabArrayBase::CopyComTagsContainer const &local_tags, DTOS const &dtos=DTOS{}, Proj const &proj=Proj{}) noexcept
 
template<class FAB , class DTOS = Identity, class Proj = Identity>
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >)> unpack_recv_buffer_gpu (FabArray< FAB > &mf, int scomp, int ncomp, Vector< char * > const &recv_data, Vector< std::size_t > const &recv_size, Vector< FabArrayBase::CopyComTagsContainer const * > const &recv_cctc, DTOS const &dtos=DTOS{}, Proj const &proj=Proj{})
 
template<typename FAB >
std::enable_if_t< IsBaseFab< FAB >::value > LocalCopy (const PackComponents &components, FabArray< FAB > &dest, const FabArray< FAB > &src, const FabArrayBase::CopyComTagsContainer &local_tags)
 Dispatch local copies to the default behaviour that knows no DTOS nor projection. More...
 
template<typename FAB >
std::enable_if_t< IsBaseFab< FAB >::value > PrepareSendBuffers (const PackComponents &components, FabArray< FAB > &dest, const FabArray< FAB > &src, CommData &comm, const FabArrayBase::MapOfCopyComTagContainers &cctc)
 Calls PrepareComBuffers. More...
 
template<typename FAB >
std::enable_if_t< IsBaseFab< FAB >::value > PrepareRecvBuffers (const PackComponents &components, FabArray< FAB > &dest, const FabArray< FAB > &src, CommData &comm, const FabArrayBase::MapOfCopyComTagContainers &cctc)
 Calls PrepareComBuffers. More...
 
template<typename FAB >
std::enable_if_t< IsBaseFab< FAB >::value > PackSendBuffers (const PackComponents &components, const FabArray< FAB > &src, CommData &send)
 Serializes FAB data without any knowledge of a DTOS nor a projection. More...
 
template<typename FAB >
std::enable_if_t< IsBaseFab< FAB >::value > UnpackRecvBuffers (const PackComponents &components, FabArray< FAB > &dest, const CommData &recv)
 De-serializes FAB data without any knowledge of a DTOS nor a projection. More...
 
template<typename FAB , typename DTOS , typename FabProj >
std::enable_if_t< IsBaseFab< FAB >::value > LocalCopy (const ApplyDtosAndProjectionOnReciever< DTOS, FabProj > &packing, FabArray< FAB > &dest, const FabArray< FAB > &src, const FabArrayBase::CopyComTagsContainer &local_tags)
 Do local copies of FABs using DTOS and projection. More...
 
template<typename FAB , typename DTOS , typename FabProj >
std::enable_if_t< IsBaseFab< FAB >::value > UnpackRecvBuffers (const ApplyDtosAndProjectionOnReciever< DTOS, FabProj > &packing, FabArray< FAB > &dest, const CommData &recv)
 Copy from received data in the buffer to destination FABs using DTOS and projection. More...
 
template<typename FAB , typename DataPacking , typename = std::enable_if_t<IsBaseFab<FAB>::value>, typename = std::enable_if_t<IsDataPacking<DataPacking, FAB>::value>>
AMREX_NODISCARD CommHandler ParallelCopy_nowait (NoLocalCopy, FabArray< FAB > &dest, const FabArray< FAB > &src, const FabArrayBase::CommMetaData &cmd, const DataPacking &data_packing)
 
template<typename FAB , typename DataPacking , typename = std::enable_if_t<IsBaseFab<FAB>::value>, typename = std::enable_if_t<IsDataPacking<DataPacking, FAB>::value>>
AMREX_NODISCARD CommHandler ParallelCopy_nowait (FabArray< FAB > &dest, const FabArray< FAB > &src, const FabArrayBase::CommMetaData &cmd, const DataPacking &data_packing)
 
template<typename FAB , typename DataPacking >
std::enable_if_t< IsBaseFab< FAB >) &&IsDataPacking< DataPacking, FAB >)> ParallelCopy_finish (FabArray< FAB > &dest, CommHandler handler, const FabArrayBase::CommMetaData &cmd, const DataPacking &data_packing)
 
template<typename FAB , typename DataPacking >
std::enable_if_t< IsBaseFab< FAB >) &&IsDataPacking< DataPacking, FAB >)> ParallelCopy_finish (DoLocalCopy, FabArray< FAB > &dest, const FabArray< FAB > &src, CommHandler handler, const FabArrayBase::CommMetaData &cmd, const DataPacking &data_packing)
 
template<typename FAB , typename DTOS = Identity, typename Proj = Identity>
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >)> ParallelCopy (FabArray< FAB > &dest, const FabArray< FAB > &src, const FabArrayBase::CommMetaData &cmd, SrcComp srccomp, DestComp destcomp, NumComps numcomp, DTOS const &dtos=DTOS{}, Proj const &proj=Proj{})
 Call ParallelCopy_nowait followed by ParallelCopy_finish, strong typed version. More...
 
template<typename FAB , typename DTOS = Identity, typename Proj = Identity>
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >)> ParallelCopy (FabArray< FAB > &dest, const FabArray< FAB > &src, const FabArrayBase::CommMetaData &cmd, int srccomp, int destcomp, int numcomp, DTOS const &dtos=DTOS{}, Proj const &proj=Proj{})
 Call ParallelCopy_nowait followed by ParallelCopy_finish. More...
 
template<typename FAB , typename DTOS = Identity, typename Proj = Identity>
std::enable_if_t< IsBaseFab< FAB >) &&IsIndexMapping< DTOS >) &&IsFabProjection< Proj, FAB >), MultiBlockCommMetaDataParallelCopy (FabArray< FAB > &dest, const Box &destbox, const FabArray< FAB > &src, SrcComp srccomp, DestComp destcomp, NumComps numcomp, const IntVect &ngrow, DTOS const &dtos=DTOS{}, Proj const &proj=Proj{})
 Call ParallelCopy_nowait followed by ParallelCopy_finish, strong typed version. More...
 
template<typename FAB , typename DTOS = Identity, typename Proj = Identity>
std::enable_if_t< IsBaseFab< FAB >) &&IsIndexMapping< DTOS >) &&IsFabProjection< Proj, FAB >), MultiBlockCommMetaDataParallelCopy (FabArray< FAB > &dest, const Box &destbox, const FabArray< FAB > &src, int srccomp, int destcomp, int numcomp, const IntVect &ngrow, DTOS const &dtos=DTOS{}, Proj const &proj=Proj{})
 Call ParallelCopy_nowait followed by ParallelCopy_finish. More...
 
template<class FAB >
std::enable_if_t< IsBaseFab< FAB >::value > Rotate90 (FabArray< FAB > &mf, int scomp, int ncomp, IntVect const &nghost, Box const &domain)
 
template<class FAB >
std::enable_if_t< IsBaseFab< FAB >::value > Rotate90 (FabArray< FAB > &mf, Box const &domain)
 
template<class FAB >
std::enable_if_t< IsBaseFab< FAB >::value > Rotate180 (FabArray< FAB > &mf, int scomp, int ncomp, IntVect const &nghost, Box const &domain)
 
template<class FAB >
std::enable_if_t< IsBaseFab< FAB >::value > Rotate180 (FabArray< FAB > &mf, Box const &domain)
 
template<class FAB >
std::enable_if_t< IsBaseFab< FAB >::value > FillPolar (FabArray< FAB > &mf, int scomp, int ncomp, IntVect const &nghost, Box const &domain)
 
template<class FAB >
std::enable_if_t< IsBaseFab< FAB >::value > FillPolar (FabArray< FAB > &mf, Box const &domain)
 
template<typename FAB , typename DTOS , typename Proj = Identity>
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >), CommHandlerFillBoundary_nowait (FabArray< FAB > &mf, const FabArrayBase::CommMetaData &cmd, int scomp, int ncomp, DTOS const &dtos, Proj const &proj=Proj{})
 Start communication to fill boundary. More...
 
template<typename FAB , typename DTOS , typename Proj = Identity>
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >)> FillBoundary_finish (CommHandler handler, FabArray< FAB > &mf, const FabArrayBase::CommMetaData &cmd, int scomp, int ncomp, DTOS const &dtos, Proj const &proj=Proj{})
 Finish communication started by FillBoundary_nowait. More...
 
template<typename FAB , typename DTOS , typename Proj = Identity>
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >)> FillBoundary (FabArray< FAB > &mf, const FabArrayBase::CommMetaData &cmd, int scomp, int ncomp, DTOS const &dtos, Proj const &proj=Proj{})
 Fill ghost cells for FabArray/MultiFab. More...
 
template<typename FAB , typename DTOS >
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >), FabArrayBase::CommMetaDatamakeFillBoundaryMetaData (FabArray< FAB > &mf, IntVect const &nghost, Geometry const &geom, DTOS const &dtos)
 Make metadata for FillBoundary. More...
 

Variables

static constexpr Identity identity {}
 
template<int I, int J>
static constexpr SwapComponents< I, J > swap_indices {}
 
static constexpr struct amrex::NonLocalBC::NoLocalCopy no_local_copy
 
static constexpr struct amrex::NonLocalBC::DoLocalCopy do_local_copy
 

Typedef Documentation

◆ DynamicSwapComponents

◆ Inverse_t

template<typename T , typename... Args>
using amrex::NonLocalBC::Inverse_t = typedef decltype(std::declval<T>().Inverse(std::declval<Args>()...))

Return type of an InverseImage class member function.

◆ LocalCopy_t

template<typename... Args>
using amrex::NonLocalBC::LocalCopy_t = typedef decltype(LocalCopy(std::declval<Args>()...))

◆ PackSendBuffers_t

template<typename... Args>
using amrex::NonLocalBC::PackSendBuffers_t = typedef decltype(PackSendBuffers(std::declval<Args>()...))

◆ PrepareRecvBuffers_t

template<typename... Args>
using amrex::NonLocalBC::PrepareRecvBuffers_t = typedef decltype(PrepareRecvBuffers(std::declval<Args>()...))

◆ PrepareSendBuffers_t

template<typename... Args>
using amrex::NonLocalBC::PrepareSendBuffers_t = typedef decltype(PrepareSendBuffers(std::declval<Args>()...))

◆ UnpackRecvBuffers_t

template<typename... Args>
using amrex::NonLocalBC::UnpackRecvBuffers_t = typedef decltype(UnpackRecvBuffers(std::declval<Args>()...))

Function Documentation

◆ Apply()

template<typename DTOS >
std::enable_if_t<IsCallableR<Dim3, DTOS, Dim3>::value, IntVect> amrex::NonLocalBC::Apply ( DTOS const &  dtos,
const IntVect iv 
)

Applies the Dim3 to Dim3 mapping onto IntVects.

This is used to map indices from the dest index space into the source index space.

IntVect is being embedded in Dim3 by trailing zeros. Dim3 is being projected to IntVect by dopping z for AMREX_SPACEDIM = 2 or z and y components for AMREX_SPACEDIM = 1.

Parameters
[in]ivThe IntVect that lives in the destination index space.
Returns
Returns IntVect{dtos(Dim3{iv})}

◆ ApplyInverse()

template<typename DTOS >
std::enable_if_t<HasInverseMemFn<DTOS>::value, IntVect> amrex::NonLocalBC::ApplyInverse ( DTOS const &  dtos,
const IntVect iv 
)

Applies the Dim3 to Dim3 invserse mapping onto IntVects.

This is used to map indices from the src index space into the dest index space.

IntVect is being embedded in Dim3 by trailing zeros. Dim3 is being projected to IntVect by dopping z for AMREX_SPACEDIM = 2 or z and y components for AMREX_SPACEDIM = 1.

Parameters
[in]ivThe IntVect that lives in the src index space.
Returns
Returns IntVect{dtos.Inverse(Dim3{iv})}

◆ FillBoundary()

template<typename FAB , typename DTOS , typename Proj = Identity>
std::enable_if_t<IsBaseFab<FAB>) && IsCallableR<Dim3,DTOS,Dim3>) && IsFabProjection<Proj,FAB>)> amrex::NonLocalBC::FillBoundary ( FabArray< FAB > &  mf,
const FabArrayBase::CommMetaData cmd,
int  scomp,
int  ncomp,
DTOS const &  dtos,
Proj const &  proj = Proj{} 
)

Fill ghost cells for FabArray/MultiFab.

This fills ghost cells of a FabArray. This function is supposed to be used together with makeFillBoundaryMetaData as follows.

    auto cmd = makeFillBoundaryMetaData(mf, mf.nGrowVect, geom, dtos);
    // The metadata cmd can be cached and reused on a MultiFab/FabArray with
    // the same BoxArray and DistributionMapping.
    FillBoundary(mf, cmd, scomp, ncomp, dtos, proj);

The FillBoundary capability here is more flexible than FabArray's FillBoundary member functions, which only fill ghost cells inside the domain and ghost cells at periodic boundaries. The FillBoundary here can be used to fill non-local domain boundaries (e.g., in spherical and cylindrical coordinates) given appropriate index and component mappings.

Template Parameters
FABMultiFab/FabArray type
DTOSIndex mapping from destination from source. See SphThetaPhiRIndexMapping for an example.
ProjComponent mapping from source to destination. See SphThetaPhiRComponentMapping for an example.
Parameters
mfFabArray/MultiFab whose ghost cells need to be filled.
cmdcommunication metadata.
scompstarting component.
ncompnumber of components.
dtosindex mapping.
projcomponent mapping.

◆ FillBoundary_finish()

template<typename FAB , typename DTOS , typename Proj = Identity>
std::enable_if_t<IsBaseFab<FAB>) && IsCallableR<Dim3,DTOS,Dim3>) && IsFabProjection<Proj,FAB>)> amrex::NonLocalBC::FillBoundary_finish ( CommHandler  handler,
FabArray< FAB > &  mf,
const FabArrayBase::CommMetaData cmd,
int  scomp,
int  ncomp,
DTOS const &  dtos,
Proj const &  proj = Proj{} 
)

Finish communication started by FillBoundary_nowait.

This finishes the communication to fill ghost cells of a FabArray. This function is supposed to be used together with FillBoundary_nowait and makeFillBoundaryMetaData as follows.

    auto cmd = makeFillBoundaryMetaData(mf, mf.nGrowVect, geom, dtos);
    // The metadata cmd can be cached and reused on a MultiFab/FabArray with
    // the same BoxArray and DistributionMapping.
    auto handler = FillBoundary_nowait(mf, cmd, scomp, ncomp, dtos, proj);
    // Independent computation can be performed.
    FillBoundary_finish(std::move(handler), mf, cmd, scomp, ncomp, dtos, proj);

The FillBoundary capability here is more flexible than FabArray's FillBoundary member functions, which only fill ghost cells inside the domain and ghost cells at periodic boundaries. The FillBoundary here can be used to fill non-local domain boundaries (e.g., in spherical and cylindrical coordinates) given appropriate index and component mappings.

Template Parameters
FABMultiFab/FabArray type
DTOSIndex mapping from destination from source. See SphThetaPhiRIndexMapping for an example.
ProjComponent mapping from source to destination. See SphThetaPhiRComponentMapping for an example.
Parameters
handlerCommHandler returned by FillBoundary_nowait.
mfFabArray/MultiFab whose ghost cells need to be filled.
cmdcommunication metadata.
scompstarting component.
ncompnumber of components.
dtosindex mapping
projcomponent mapping

◆ FillBoundary_nowait()

template<typename FAB , typename DTOS , typename Proj = Identity>
std::enable_if_t<IsBaseFab<FAB>) && IsCallableR<Dim3,DTOS,Dim3>) && IsFabProjection<Proj,FAB>), CommHandler> amrex::NonLocalBC::FillBoundary_nowait ( FabArray< FAB > &  mf,
const FabArrayBase::CommMetaData cmd,
int  scomp,
int  ncomp,
DTOS const &  dtos,
Proj const &  proj = Proj{} 
)

Start communication to fill boundary.

This starts communication to fill ghost cells of a FabArray. This function is supposed to be used together with FillBoundary_finish and makeFillBoundaryMetaData as follows.

    auto cmd = makeFillBoundaryMetaData(mf, mf.nGrowVect, geom, dtos);
    // The metadata cmd can be cached and reused on a MultiFab/FabArray with
    // the same BoxArray and DistributionMapping.
    auto handler = FillBoundary_nowait(mf, cmd, scomp, ncomp, dtos, proj);
    // Independent computation can be performed.
    FillBoundary_finish(std::move(handler), mf, cmd, scomp, ncomp, dtos, proj);

The FillBoundary capability here is more flexible than FabArray's FillBoundary member functions, which only fill ghost cells inside the domain and ghost cells at periodic boundaries. The FillBoundary here can be used to fill non-local domain boundaries (e.g., in spherical and cylindrical coordinates) given appropriate index and component mappings.

Template Parameters
FABMultiFab/FabArray type
DTOSIndex mapping from destination from source. See SphThetaPhiRIndexMapping for an example.
ProjComponent mapping from source to destination. See SphThetaPhiRComponentMapping for an example.
Parameters
mfFabArray/MultiFab whose ghost cells need to be filled.
cmdcommunication metadata.
scompstarting component.
ncompnumber of components.
dtosindex mapping.
projcomponent mapping.
Returns
a CommHandler object needed for calling FillBoundary_finish.

◆ FillPolar() [1/2]

template<class FAB >
std::enable_if_t<IsBaseFab<FAB>::value> amrex::NonLocalBC::FillPolar ( FabArray< FAB > &  mf,
Box const &  domain 
)

◆ FillPolar() [2/2]

template<class FAB >
std::enable_if_t<IsBaseFab<FAB>::value> amrex::NonLocalBC::FillPolar ( FabArray< FAB > &  mf,
int  scomp,
int  ncomp,
IntVect const &  nghost,
Box const &  domain 
)

◆ Image() [1/2]

template<typename DTOS >
std::enable_if_t<IsCallableR<Dim3, DTOS, Dim3>::value && !IsCallableR<IndexType, DTOS, IndexType>::value, Box> amrex::NonLocalBC::Image ( DTOS const &  dtos,
const Box box 
)

Applies the Dim3 to Dim3 mapping onto Boxes but does not change the index type.

This function assumes monotonicity of dtos in each component.

Parameters
[in]boxThe box that lives the destination index space
Returns
Returns the smallest Box in the source index space that contains images of Apply(dtos, box.smallEnd()) and Apply(dtos, box.bigEnd()).

◆ Image() [2/2]

template<typename DTOS >
std::enable_if_t<IsCallableR<Dim3, DTOS, Dim3>::value && IsCallableR<IndexType, DTOS, IndexType>::value, Box> amrex::NonLocalBC::Image ( DTOS const &  dtos,
const Box box 
)

Applies the Dim3 to Dim3 mapping onto Boxes and maps the index type.

This function assumes monotonicity of dtos in each component.

Parameters
[in]boxThe box that lives the destination index space
Returns
Returns the smallest Box in the source index space that contains images of Apply(dtos, box.smallEnd()) and Apply(dtos, box.bigEnd()).

◆ InverseImage() [1/2]

template<typename DTOS >
std::enable_if_t<HasInverseMemFn<DTOS>::value && !IsCallableR<IndexType, DTOS, IndexType>::value, Box> amrex::NonLocalBC::InverseImage ( DTOS const &  dtos,
const Box box 
)

Applies the inverse Dim3 to Dim3 mapping onto Boxes without changing the index type.

This function assumes monotonicity in each component.

Parameters
[in]boxThe box that lives the source index space
Returns
Returns the smallest Box in the destination index space that contains images of ApplyInverse(box.smallEnd()) and ApplyInverse(box.bigEnd()).

◆ InverseImage() [2/2]

template<typename DTOS >
std::enable_if_t<HasInverseMemFn<DTOS>::value && IsCallableR<IndexType, DTOS, IndexType>::value, Box> amrex::NonLocalBC::InverseImage ( DTOS const &  dtos,
const Box box 
)

Applies the inverse Dim3 to Dim3 mapping onto Boxes.

This function assumes monotonicity in each component.

Parameters
[in]boxThe box that lives the source index space
Returns
Returns the smallest Box in the destination index space that contains images of ApplyInverse(box.smallEnd()) and ApplyInverse(box.bigEnd()).

◆ local_copy_cpu()

template<class FAB , class DTOS = Identity, class Proj = Identity>
std::enable_if_t<IsBaseFab<FAB>) && IsCallableR<Dim3, DTOS, Dim3>) && IsFabProjection<Proj, FAB>)> amrex::NonLocalBC::local_copy_cpu ( FabArray< FAB > &  dest,
const FabArray< FAB > &  src,
int  dcomp,
int  scomp,
int  ncomp,
FabArrayBase::CopyComTagsContainer const &  local_tags,
DTOS const &  dtos = DTOS{},
Proj const &  proj = Proj{} 
)
noexcept

◆ local_copy_gpu()

template<class FAB , class DTOS = Identity, class Proj = Identity>
std::enable_if_t<IsBaseFab<FAB>) && IsCallableR<Dim3, DTOS, Dim3>) && IsFabProjection<Proj, FAB>)> amrex::NonLocalBC::local_copy_gpu ( FabArray< FAB > &  dest,
const FabArray< FAB > &  src,
int  dcomp,
int  scomp,
int  ncomp,
FabArrayBase::CopyComTagsContainer const &  local_tags,
DTOS const &  dtos = DTOS{},
Proj const &  proj = Proj{} 
)
noexcept

◆ LocalCopy() [1/2]

template<typename FAB , typename DTOS , typename FabProj >
std::enable_if_t<IsBaseFab<FAB>::value> amrex::NonLocalBC::LocalCopy ( const ApplyDtosAndProjectionOnReciever< DTOS, FabProj > &  packing,
FabArray< FAB > &  dest,
const FabArray< FAB > &  src,
const FabArrayBase::CopyComTagsContainer local_tags 
)

Do local copies of FABs using DTOS and projection.

◆ LocalCopy() [2/2]

template<typename FAB >
std::enable_if_t<IsBaseFab<FAB>::value> amrex::NonLocalBC::LocalCopy ( const PackComponents components,
FabArray< FAB > &  dest,
const FabArray< FAB > &  src,
const FabArrayBase::CopyComTagsContainer local_tags 
)

Dispatch local copies to the default behaviour that knows no DTOS nor projection.

◆ makeFillBoundaryMetaData()

template<typename FAB , typename DTOS >
std::enable_if_t<IsBaseFab<FAB>) && IsCallableR<Dim3,DTOS,Dim3>), FabArrayBase::CommMetaData> amrex::NonLocalBC::makeFillBoundaryMetaData ( FabArray< FAB > &  mf,
IntVect const &  nghost,
Geometry const &  geom,
DTOS const &  dtos 
)

Make metadata for FillBoundary.

Template Parameters
FABMultiFab/FabArray type
DTOSIndex mapping from destination from source. See SphThetaPhiRIndexMapping for an example.
Parameters
mfFabArray/MultiFab whose ghost cells need to be filled.
nghostnumber of ghost cells to be filled.
geoma Geometry object that contains the domain information.
dtosindex mapping.
Returns
communication metadata

◆ PackSendBuffers()

template<typename FAB >
std::enable_if_t<IsBaseFab<FAB>::value> amrex::NonLocalBC::PackSendBuffers ( const PackComponents components,
const FabArray< FAB > &  src,
CommData send 
)

Serializes FAB data without any knowledge of a DTOS nor a projection.

◆ ParallelCopy() [1/5]

template<typename FAB , typename DTOS = Identity, typename Proj = Identity>
std::enable_if_t<IsBaseFab<FAB>) && IsIndexMapping<DTOS>) && IsFabProjection<Proj, FAB>),MultiBlockCommMetaData> amrex::NonLocalBC::ParallelCopy ( FabArray< FAB > &  dest,
const Box destbox,
const FabArray< FAB > &  src,
int  srccomp,
int  destcomp,
int  numcomp,
const IntVect ngrow,
DTOS const &  dtos = DTOS{},
Proj const &  proj = Proj{} 
)

Call ParallelCopy_nowait followed by ParallelCopy_finish.

This function constructs a new MultiCommMetaData from the given DTOS, destbox and ngrow.

Parameters
[out]destThe Multifab that is going to be filled with received data.
[in]destboxThe index box in the destination space that will be filled by data from src. The source box that describes the dependencies will be computed by the specified DTOS.
[in]srcThe Multifab that is used to fill the send buffers.
[in]srccompThe first component in src that will be copied to dest.
[in]destcompThe first component in dest that will get written by src.
[in]ncompThe number of successive components that will be copied.
[in]ngrowThe amount of ghost cells that will be taking into consideration. Note, even if destbox contains indices outside the domain we need to specify an appropriate ngrow that covers the amount of ghost cells that we want to copy.
[in]dtosAn index mapping that maps indices from destination space to source space and from source space to destination space.
[in]projA transformation function that might change the data when it is being copied.
Returns
Returns the CommMetaData object that can be cached for future calls to ParallelCopy.

◆ ParallelCopy() [2/5]

template<typename FAB , typename DTOS = Identity, typename Proj = Identity>
std::enable_if_t<IsBaseFab<FAB>) && IsIndexMapping<DTOS>) && IsFabProjection<Proj, FAB>),MultiBlockCommMetaData> amrex::NonLocalBC::ParallelCopy ( FabArray< FAB > &  dest,
const Box destbox,
const FabArray< FAB > &  src,
SrcComp  srccomp,
DestComp  destcomp,
NumComps  numcomp,
const IntVect ngrow,
DTOS const &  dtos = DTOS{},
Proj const &  proj = Proj{} 
)

Call ParallelCopy_nowait followed by ParallelCopy_finish, strong typed version.

This function constructs a new MultiCommMetaData from the given DTOS, destbox and ngrow.

Parameters
[out]destThe Multifab that is going to be filled with received data.
[in]destboxThe index box in the destination space that will be filled by data from src. The source box that describes the dependencies will be computed by the specified DTOS.
[in]srcThe Multifab that is used to fill the send buffers.
[in]srccompThe first component in src that will be copied to dest.
[in]destcompThe first component in dest that will get written by src.
[in]ncompThe number of successive components that will be copied.
[in]ngrowThe amount of ghost cells that will be taking into consideration. Note, even if destbox contains indices outside the domain we need to specify an appropriate ngrow that covers the amount of ghost cells that we want to copy.
[in]dtosAn index mapping that maps indices from destination space to source space and from source space to destination space.
[in]projA transformation function that might change the data when it is being copied.
Returns
Returns the CommMetaData object that can be cached for future calls to ParallelCopy.

◆ ParallelCopy() [3/5]

template<typename FAB , typename DTOS = Identity, typename Proj = Identity>
std::enable_if_t<IsBaseFab<FAB>) && IsCallableR<Dim3, DTOS, Dim3>) && IsFabProjection<Proj, FAB>)> amrex::NonLocalBC::ParallelCopy ( FabArray< FAB > &  dest,
const FabArray< FAB > &  src,
const FabArrayBase::CommMetaData cmd,
int  srccomp,
int  destcomp,
int  numcomp,
DTOS const &  dtos = DTOS{},
Proj const &  proj = Proj{} 
)

Call ParallelCopy_nowait followed by ParallelCopy_finish.

This function overload uses an already cached CommMetaData. This CommMetaData needs to be compatible with the specified DTOS and projection, otherwise undefined behaviour occurs.

Parameters
[out]destThe Multifab that is going to be filled with received data.
[in]srcThe Multifab that is used to fill the send buffers.
[in]cmdThe communication meta data object holds spatial information about FAB boxes that need to be filled and copied from.
[in]srccompThe first component in src that will be copied to dest.
[in]destcompThe first component in dest that will get written by src.
[in]ncompThe number of successive components that will be copied.
[in]dtosAn index mapping that maps indices from destination space to source space and from source space to destination space.
[in]projA transformation function that might change the data when it is being copied.
Returns
Nothing.

◆ ParallelCopy() [4/5]

template<typename FAB , typename DTOS = Identity, typename Proj = Identity>
std::enable_if_t<IsBaseFab<FAB>) && IsCallableR<Dim3, DTOS, Dim3>) && IsFabProjection<Proj, FAB>)> amrex::NonLocalBC::ParallelCopy ( FabArray< FAB > &  dest,
const FabArray< FAB > &  src,
const FabArrayBase::CommMetaData cmd,
SrcComp  srccomp,
DestComp  destcomp,
NumComps  numcomp,
DTOS const &  dtos = DTOS{},
Proj const &  proj = Proj{} 
)

Call ParallelCopy_nowait followed by ParallelCopy_finish, strong typed version.

This function overload uses an already cached CommMetaData. This CommMetaData needs to be compatible with the specified DTOS and projection, otherwise undefined behaviour occurs.

Parameters
[out]destThe Multifab that is going to be filled with received data.
[in]srcThe Multifab that is used to fill the send buffers.
[in]cmdThe communication meta data object holds spatial information about FAB boxes that need to be filled and copied from.
[in]srccompThe first component in src that will be copied to dest.
[in]destcompThe first component in dest that will get written by src.
[in]ncompThe number of successive components that will be copied.
[in]dtosAn index mapping that maps indices from destination space to source space and from source space to destination space.
[in]projA transformation function that might change the data when it is being copied.
Returns
Nothing.

◆ ParallelCopy() [5/5]

template MultiBlockCommMetaData amrex::NonLocalBC::ParallelCopy ( FabArray< FArrayBox > &  dest,
const Box destbox,
const FabArray< FArrayBox > &  src,
int  destcomp,
int  srccomp,
int  numcomp,
const IntVect ngrow,
MultiBlockIndexMapping const &  ,
Identity const &   
)

◆ ParallelCopy_finish() [1/2]

template<typename FAB , typename DataPacking >
std::enable_if_t<IsBaseFab<FAB>) && IsDataPacking<DataPacking, FAB>)> amrex::NonLocalBC::ParallelCopy_finish ( DoLocalCopy  ,
FabArray< FAB > &  dest,
const FabArray< FAB > &  src,
CommHandler  handler,
const FabArrayBase::CommMetaData cmd,
const DataPacking &  data_packing 
)

Blockingly wait for all communication to be done and fill the local FABs with received data.

This function overload performs local copies, i.e. from this MPI process to itself. It will block the current thread until all MPI recv and send requests are done and calls the DataPacking object to unpack the received buffers.

Parameters
[out]destThe Multifab that is going to be filled with received data.
[in]srcThe Multifab that is used to get the data for the local copies from.
[in]handlerThis object holds all data buffers that need to be kept alive as long as the data transaction is not done.
[in]cmdThe communication meta data object holds spatial information about FAB boxes that need to be filled.
[in]data_packingA CPO that controls behaviour of unpacking the received data.
Returns
Nothing.

◆ ParallelCopy_finish() [2/2]

template<typename FAB , typename DataPacking >
std::enable_if_t<IsBaseFab<FAB>) && IsDataPacking<DataPacking, FAB>)> amrex::NonLocalBC::ParallelCopy_finish ( FabArray< FAB > &  dest,
CommHandler  handler,
const FabArrayBase::CommMetaData cmd,
const DataPacking &  data_packing 
)

Blockingly wait for all communication to be done and fill the local FABs with received data.

This function overload performs no local copies, i.e. from this MPI process to itself. It will block the current thread until all MPI recv and send requests are done and calls the DataPacking object to unpack the received buffers.

Parameters
[out]destThe Multifab that is going to be filled with received data.
[in]handlerThis object holds all data buffers that need to be kept alive as long as the data transaction is not done.
[in]cmdThe communication meta data object holds spatial information about FAB boxes that need to be filled.
[in]data_packingA CPO that controls behaviour of unpacking the received buffer to the destination FabArray.
Returns
Nothing.

◆ ParallelCopy_nowait() [1/2]

template<typename FAB , typename DataPacking , typename = std::enable_if_t<IsBaseFab<FAB>::value>, typename = std::enable_if_t<IsDataPacking<DataPacking, FAB>::value>>
AMREX_NODISCARD CommHandler amrex::NonLocalBC::ParallelCopy_nowait ( FabArray< FAB > &  dest,
const FabArray< FAB > &  src,
const FabArrayBase::CommMetaData cmd,
const DataPacking &  data_packing 
)

Initiate recv and send calls for MPI and return after doing local work.

DataPacking is a customization point object to control the behaviour of packing and unpacking send or recv data buffers. It is used to perform interpolation or data transformations on either sender or receiver side.

This function performs a data packing on sender side and we expect a call to Parallel_finish that performs data unpacking on the receiver side.

Parameters
[out]destThe Multifab that is going to be filled with received data.
[in]srcThe Multifab that is used to fill the send buffers.
[in]cmdThe communication meta data object holds spatial information about FAB boxes that need to be filled and copied from.
[in]data_packingA CPO that controls behaviour of preparing buffers and packing the source data into the send buffers.
Returns
Returns a CommHandler object that owns context and memory buffers for the whole life time of the MPI transaction.

◆ ParallelCopy_nowait() [2/2]

template<typename FAB , typename DataPacking , typename = std::enable_if_t<IsBaseFab<FAB>::value>, typename = std::enable_if_t<IsDataPacking<DataPacking, FAB>::value>>
AMREX_NODISCARD CommHandler amrex::NonLocalBC::ParallelCopy_nowait ( NoLocalCopy  ,
FabArray< FAB > &  dest,
const FabArray< FAB > &  src,
const FabArrayBase::CommMetaData cmd,
const DataPacking &  data_packing 
)

Initiate recv and send calls for MPI and immediately return without doing any work.

DataPacking is a customization point object to control the behaviour of packing and unpacking send or recv data buffers. It is used to perform interpolation or data transformations on either sender or receiver side.

This function performs a data packing on sender side and we expect a call to Parallel_finish that performs data unpacking on the receiver side.

Parameters
[out]destThe Multifab that is going to be filled with received data.
[in]srcThe Multifab that is used to fill the send buffers.
[in]cmdThe communication meta data object holds spatial information about FAB boxes that need to be filled and copied from.
[in]data_packingA CPO that controls behaviour of preparing buffers and packing the source data into the send buffers.
Returns
Returns a CommHandler object that owns context and memory buffers for the whole life time of the MPI transaction.

◆ PostRecvs()

void amrex::NonLocalBC::PostRecvs ( CommData recv,
int  mpi_tag 
)

Initiate all recvieves with MPI_Irecv calls associated with tag mpi_tag.

◆ PostSends()

void amrex::NonLocalBC::PostSends ( CommData send,
int  mpi_tag 
)

Initiate all sends with MPI_Isend calls associated with tag mpi_tag.

◆ PrepareCommBuffers()

void amrex::NonLocalBC::PrepareCommBuffers ( CommData comm,
const FabArrayBase::MapOfCopyComTagContainers cctc,
int  n_components,
std::size_t  object_size,
std::size_t  align 
)

Fill all class member variables of comm but the request and the stats vector.

◆ PrepareRecvBuffers()

template<typename FAB >
std::enable_if_t<IsBaseFab<FAB>::value> amrex::NonLocalBC::PrepareRecvBuffers ( const PackComponents components,
FabArray< FAB > &  dest,
const FabArray< FAB > &  src,
CommData comm,
const FabArrayBase::MapOfCopyComTagContainers cctc 
)

Calls PrepareComBuffers.

◆ PrepareSendBuffers()

template<typename FAB >
std::enable_if_t<IsBaseFab<FAB>::value> amrex::NonLocalBC::PrepareSendBuffers ( const PackComponents components,
FabArray< FAB > &  dest,
const FabArray< FAB > &  src,
CommData comm,
const FabArrayBase::MapOfCopyComTagContainers cctc 
)

Calls PrepareComBuffers.

◆ Rotate180() [1/2]

template<class FAB >
std::enable_if_t<IsBaseFab<FAB>::value> amrex::NonLocalBC::Rotate180 ( FabArray< FAB > &  mf,
Box const &  domain 
)

◆ Rotate180() [2/2]

template<class FAB >
std::enable_if_t<IsBaseFab<FAB>::value> amrex::NonLocalBC::Rotate180 ( FabArray< FAB > &  mf,
int  scomp,
int  ncomp,
IntVect const &  nghost,
Box const &  domain 
)

◆ Rotate90() [1/2]

template<class FAB >
std::enable_if_t<IsBaseFab<FAB>::value> amrex::NonLocalBC::Rotate90 ( FabArray< FAB > &  mf,
Box const &  domain 
)

◆ Rotate90() [2/2]

template<class FAB >
std::enable_if_t<IsBaseFab<FAB>::value> amrex::NonLocalBC::Rotate90 ( FabArray< FAB > &  mf,
int  scomp,
int  ncomp,
IntVect const &  nghost,
Box const &  domain 
)

◆ unpack_recv_buffer_cpu()

template<class FAB , class DTOS = Identity, class Proj = Identity>
std::enable_if_t<IsBaseFab<FAB>) && IsCallableR<Dim3, DTOS, Dim3>) && IsFabProjection<Proj, FAB>)> amrex::NonLocalBC::unpack_recv_buffer_cpu ( FabArray< FAB > &  mf,
int  dcomp,
int  ncomp,
Vector< char * > const &  recv_data,
Vector< std::size_t > const &  recv_size,
Vector< FabArrayBase::CopyComTagsContainer const * > const &  recv_cctc,
DTOS const &  dtos = DTOS{},
Proj const &  proj = Proj{} 
)
noexcept

◆ unpack_recv_buffer_gpu()

template<class FAB , class DTOS = Identity, class Proj = Identity>
std::enable_if_t<IsBaseFab<FAB>) && IsCallableR<Dim3, DTOS, Dim3>) && IsFabProjection<Proj, FAB>)> amrex::NonLocalBC::unpack_recv_buffer_gpu ( FabArray< FAB > &  mf,
int  scomp,
int  ncomp,
Vector< char * > const &  recv_data,
Vector< std::size_t > const &  recv_size,
Vector< FabArrayBase::CopyComTagsContainer const * > const &  recv_cctc,
DTOS const &  dtos = DTOS{},
Proj const &  proj = Proj{} 
)

◆ UnpackRecvBuffers() [1/2]

template<typename FAB , typename DTOS , typename FabProj >
std::enable_if_t<IsBaseFab<FAB>::value> amrex::NonLocalBC::UnpackRecvBuffers ( const ApplyDtosAndProjectionOnReciever< DTOS, FabProj > &  packing,
FabArray< FAB > &  dest,
const CommData recv 
)

Copy from received data in the buffer to destination FABs using DTOS and projection.

◆ UnpackRecvBuffers() [2/2]

template<typename FAB >
std::enable_if_t<IsBaseFab<FAB>::value> amrex::NonLocalBC::UnpackRecvBuffers ( const PackComponents components,
FabArray< FAB > &  dest,
const CommData recv 
)

De-serializes FAB data without any knowledge of a DTOS nor a projection.

Variable Documentation

◆ do_local_copy

constexpr struct amrex::NonLocalBC::DoLocalCopy amrex::NonLocalBC::do_local_copy
static

◆ identity

constexpr Identity amrex::NonLocalBC::identity {}
staticconstexpr

◆ no_local_copy

constexpr struct amrex::NonLocalBC::NoLocalCopy amrex::NonLocalBC::no_local_copy
static

◆ swap_indices

template<int I, int J>
constexpr SwapComponents<I, J> amrex::NonLocalBC::swap_indices {}
staticconstexpr