Finite Volume Solver
prototype
A framework to build finite volume solvers for the AG Klein at the Freie Universität Berlin.
|
This class manages data and the numerical method to perform cut cell simulations with the AMReX library. More...
#include <IntegratorContext.hpp>
Classes | |
struct | LevelData |
Public Types | |
using | FeedbackFn = std::function< void(IntegratorContext &, int, Duration, std::pair< int, int >)> |
Public Member Functions | |
Constructors, Assignment Operators and Desctructor | |
IntegratorContext (std::shared_ptr< GriddingAlgorithm > gridding, HyperbolicMethod method) | |
Deeply copies a context and all its distributed data for all MPI ranks. More... | |
IntegratorContext (std::shared_ptr< GriddingAlgorithm > gridding, HyperbolicMethod method, int cell_gcw, int face_gcw) | |
Deeply copies a context and all its distributed data for all MPI ranks. More... | |
IntegratorContext (const IntegratorContext &) | |
Deeply copies a context and all its distributed data for all MPI ranks. More... | |
IntegratorContext & | operator= (const IntegratorContext &) |
Deeply copies a context and all its distributed data for all MPI ranks. More... | |
IntegratorContext (IntegratorContext &&) noexcept=default | |
Deeply copies a context and all its distributed data for all MPI ranks. More... | |
IntegratorContext & | operator= (IntegratorContext &&) noexcept=default |
Deeply copies a context and all its distributed data for all MPI ranks. More... | |
~IntegratorContext ()=default | |
Deeply copies a context and all its distributed data for all MPI ranks. More... | |
Member Accessors | |
const HyperbolicMethod & | GetHyperbolicMethod () const noexcept |
Returns the hyperbolic method member object. More... | |
const std::shared_ptr< GriddingAlgorithm > & | GetGriddingAlgorithm () const noexcept |
Returns a shared pointer to the underlying GriddingAlgorithm which owns the simulation. More... | |
const std::shared_ptr< CounterRegistry > & | GetCounterRegistry () const noexcept |
Returns a shared pointer to the counter registry. More... | |
const PatchHierarchy & | GetPatchHierarchy () const noexcept |
Returns a reference to const PatchHierarchy which is a member of the GriddingAlgorithm. More... | |
PatchHierarchy & | GetPatchHierarchy () noexcept |
Returns the hyperbolic method member object. More... | |
MPI_Comm | GetMpiCommunicator () const noexcept |
Returns the MPI communicator which is associated with this context. More... | |
void | SetFeedbackFunction (FeedbackFn feedback) |
Call this function upon every time integration. More... | |
Access Level-specific data | |
const AnyBoundaryCondition & | GetBoundaryCondition () const |
Returns the current boundary condition for the specified level. More... | |
AnyBoundaryCondition & | GetBoundaryCondition () |
Returns the current boundary condition for the specified level. More... | |
::amrex::EBFArrayBoxFactory & | GetEmbeddedBoundary (int level) |
Returns the current boundary condition for the specified level. More... | |
const ::amrex::EBFArrayBoxFactory & | GetEmbeddedBoundary (int level) const |
Returns the current boundary condition for the specified level. More... | |
::amrex::MultiFab & | GetData (int level) |
Returns the MultiFab associated with level data on the specifed level number. More... | |
::amrex::MultiFab & | GetReferenceStates (int level) |
Returns the MultiFab associated with flux data on the specifed level number and direction. More... | |
::amrex::MultiFab & | GetScratch (int level) |
Returns the MultiFab associated with level data with ghost cells on the specifed level number and direction. More... | |
::amrex::MultiCutFab & | GetBoundaryFluxes (int level) |
Returns the MultiFab associated with flux data on the specifed level number and direction. More... | |
::amrex::MultiFab & | GetFluxes (int level, Direction dir) |
Returns the MultiFab associated with flux data on the specifed level number and direction. More... | |
::amrex::MultiFab & | GetStabilizedFluxes (int level, Direction dir) |
Returns the MultiFab associated with flux data on the specifed level number and direction. More... | |
::amrex::MultiFab & | GetShieldedFromLeftFluxes (int level, Direction dir) |
Returns the MultiFab associated with flux data on the specifed level number and direction. More... | |
::amrex::MultiFab & | GetShieldedFromRightFluxes (int level, Direction dir) |
Returns the MultiFab associated with flux data on the specifed level number and direction. More... | |
::amrex::MultiFab & | GetDoublyShieldedFluxes (int level, Direction dir) |
Returns the MultiFab associated with flux data on the specifed level number and direction. More... | |
Duration | GetTimePoint (int level=0) const |
Returns the current time level for data at the specified refinement level and direction. More... | |
std::ptrdiff_t | GetCycles (int level=0) const |
Returns the current number of cycles for data at the specified refinement level and direction. More... | |
const ::amrex::Geometry & | GetGeometry (int level) const |
Returns the geometry object for the specified refinement level. More... | |
Observers | |
Returns true if the data exists for the specified level number. | |
bool | LevelExists (int level) const noexcept |
Returns the refinement ratio in the specified direction. More... | |
int | GetRatioToCoarserLevel (int level, Direction dir) const noexcept |
Returns the refinement ratio in the specified direction. More... | |
::amrex::IntVect | GetRatioToCoarserLevel (int level) const noexcept |
Returns the refinement ratio for all directions. More... | |
double | GetDx (int level, Direction dir) const noexcept |
Returns the refinement ratio in the specified direction. More... | |
::amrex::FabType | GetFabType (int level, const ::amrex::MFIter &mfi) const noexcept |
Returns the refinement ratio in the specified direction. More... | |
Modifiers | |
void | ResetHierarchyConfiguration (std::shared_ptr< GriddingAlgorithm > gridding) |
Replaces the underlying gridding algorithm with the specified one. More... | |
void | ResetHierarchyConfiguration (int level=0) |
Whenever the gridding algorithm changes the data hierarchy this function will regrid all distributed helper variables managed by the context. More... | |
void | SetCycles (std::ptrdiff_t cycle, int level) |
Sets the cycle count for a specific level number and direction. More... | |
void | SetTimePoint (Duration t, int level) |
Sets the time point for a specific level number and direction. More... | |
Member functions relevant for the level integrator algorithm. | |
void | PreAdvanceHierarchy () |
Compute and store reference states for cut-cells over all split cycles. More... | |
void | PostAdvanceHierarchy () |
Compute and store reference states for cut-cells over all split cycles. More... | |
int | PreAdvanceLevel (int level_num, Duration dt, std::pair< int, int > subcycle) |
On each first subcycle this will regrid the data if neccessary. More... | |
Result< void, TimeStepTooLarge > | PostAdvanceLevel (int level_num, Duration dt, std::pair< int, int > subcycle) |
Increases the internal time stamps and cycle counters for the specified level number and direction. More... | |
void | CopyDataToScratch (int level) |
Compute and store reference states for cut-cells over all split cycles. More... | |
void | CopyScratchToData (int level) |
Compute and store reference states for cut-cells over all split cycles. More... | |
void | ApplyBoundaryCondition (int level, Direction dir) |
Applies the boundary condition for the scratch space on level level in direcition dir . More... | |
void | ApplyBoundaryCondition (int level, Direction dir, AnyBoundaryCondition &bc) |
Compute and store reference states for cut-cells over all split cycles. More... | |
void | FillGhostLayerTwoLevels (int level, AnyBoundaryCondition &fbc, int coarse, AnyBoundaryCondition &cbc) |
Fills the ghost layer of the scratch data and interpolates in the coarse fine layer. More... | |
void | FillGhostLayerTwoLevels (int level, int coarse) |
Compute and store reference states for cut-cells over all split cycles. More... | |
void | FillGhostLayerSingleLevel (int level, AnyBoundaryCondition &bc) |
Fills the ghost layer of the scratch data and does nothing in the coarse fine layer. More... | |
void | FillGhostLayerSingleLevel (int level) |
Compute and store reference states for cut-cells over all split cycles. More... | |
Duration | ComputeStableDt (int level, Direction dir) |
Returns a estimate for a stable time step size which can be taken for specified level number in direction dir. More... | |
void | ComputeNumericFluxes (int level, Duration dt, Direction dir) |
Fill the flux MultiFab with numeric fluxes based on current states in scratch. More... | |
void | UpdateConservatively (int level, Duration dt, Direction dir) |
Apply a conservative time update for each conservative variable on the specified level number and direction. More... | |
void | CompleteFromCons (int level, Duration dt) |
Reconstruct complete state variables from conservative ones. More... | |
void | AccumulateCoarseFineFluxes (int level, double time_scale, Direction dir) |
Accumulate fluxes on the coarse fine interfaces for a specified fine level number. More... | |
void | ApplyFluxCorrection (int fine, int coarse, Duration dt) |
Replace the coarse fluxes by accumulated fine fluxes on the coarse fine interfaces. More... | |
void | ResetCoarseFineFluxes (int fine, int coarse) |
Resets all accumulates fluxes to zero. More... | |
void | CoarsenConservatively (int fine, int coarse) |
Coarsen scratch data from a fine level number to a coarse level number. More... | |
Public Attributes | |
std::shared_ptr< CounterRegistry > | registry_ |
bool | count_per_level |
Private Attributes | |
int | scratch_ghost_cell_width_ |
int | flux_ghost_cell_width_ |
std::shared_ptr< GriddingAlgorithm > | gridding_ |
std::vector< LevelData > | data_ |
HyperbolicMethod | method_ |
FeedbackFn | feedback_ |
This class manages data and the numerical method to perform cut cell simulations with the AMReX library.
using fub::amrex::cutcell::IntegratorContext::FeedbackFn = std::function<void(IntegratorContext&, int, Duration, std::pair<int, int>)> |
fub::amrex::cutcell::IntegratorContext::IntegratorContext | ( | std::shared_ptr< GriddingAlgorithm > | gridding, |
HyperbolicMethod | method | ||
) |
Deeply copies a context and all its distributed data for all MPI ranks.
fub::amrex::cutcell::IntegratorContext::IntegratorContext | ( | std::shared_ptr< GriddingAlgorithm > | gridding, |
HyperbolicMethod | method, | ||
int | cell_gcw, | ||
int | face_gcw | ||
) |
Deeply copies a context and all its distributed data for all MPI ranks.
fub::amrex::cutcell::IntegratorContext::IntegratorContext | ( | const IntegratorContext & | ) |
Deeply copies a context and all its distributed data for all MPI ranks.
|
defaultnoexcept |
Deeply copies a context and all its distributed data for all MPI ranks.
|
default |
Deeply copies a context and all its distributed data for all MPI ranks.
void fub::amrex::cutcell::IntegratorContext::AccumulateCoarseFineFluxes | ( | int | level, |
double | time_scale, | ||
Direction | dir | ||
) |
Accumulate fluxes on the coarse fine interfaces for a specified fine level number.
void fub::amrex::cutcell::IntegratorContext::ApplyBoundaryCondition | ( | int | level, |
Direction | dir | ||
) |
Applies the boundary condition for the scratch space on level level
in direcition dir
.
level | The refinement level on which the boundary condition shall be used. |
void fub::amrex::cutcell::IntegratorContext::ApplyBoundaryCondition | ( | int | level, |
Direction | dir, | ||
AnyBoundaryCondition & | bc | ||
) |
Compute and store reference states for cut-cells over all split cycles.
void fub::amrex::cutcell::IntegratorContext::ApplyFluxCorrection | ( | int | fine, |
int | coarse, | ||
Duration | dt | ||
) |
Replace the coarse fluxes by accumulated fine fluxes on the coarse fine interfaces.
void fub::amrex::cutcell::IntegratorContext::CoarsenConservatively | ( | int | fine, |
int | coarse | ||
) |
Coarsen scratch data from a fine level number to a coarse level number.
void fub::amrex::cutcell::IntegratorContext::CompleteFromCons | ( | int | level, |
Duration | dt | ||
) |
Reconstruct complete state variables from conservative ones.
void fub::amrex::cutcell::IntegratorContext::ComputeNumericFluxes | ( | int | level, |
Duration | dt, | ||
Direction | dir | ||
) |
Fill the flux MultiFab with numeric fluxes based on current states in scratch.
Returns a estimate for a stable time step size which can be taken for specified level number in direction dir.
void fub::amrex::cutcell::IntegratorContext::CopyDataToScratch | ( | int | level | ) |
Compute and store reference states for cut-cells over all split cycles.
void fub::amrex::cutcell::IntegratorContext::CopyScratchToData | ( | int | level | ) |
Compute and store reference states for cut-cells over all split cycles.
void fub::amrex::cutcell::IntegratorContext::FillGhostLayerSingleLevel | ( | int | level | ) |
Compute and store reference states for cut-cells over all split cycles.
void fub::amrex::cutcell::IntegratorContext::FillGhostLayerSingleLevel | ( | int | level, |
AnyBoundaryCondition & | bc | ||
) |
Fills the ghost layer of the scratch data and does nothing in the coarse fine layer.
void fub::amrex::cutcell::IntegratorContext::FillGhostLayerTwoLevels | ( | int | level, |
AnyBoundaryCondition & | fbc, | ||
int | coarse, | ||
AnyBoundaryCondition & | cbc | ||
) |
Fills the ghost layer of the scratch data and interpolates in the coarse fine layer.
void fub::amrex::cutcell::IntegratorContext::FillGhostLayerTwoLevels | ( | int | level, |
int | coarse | ||
) |
Compute and store reference states for cut-cells over all split cycles.
AnyBoundaryCondition& fub::amrex::cutcell::IntegratorContext::GetBoundaryCondition | ( | ) |
Returns the current boundary condition for the specified level.
const AnyBoundaryCondition& fub::amrex::cutcell::IntegratorContext::GetBoundaryCondition | ( | ) | const |
Returns the current boundary condition for the specified level.
::amrex::MultiCutFab& fub::amrex::cutcell::IntegratorContext::GetBoundaryFluxes | ( | int | level | ) |
Returns the MultiFab associated with flux data on the specifed level number and direction.
|
noexcept |
Returns a shared pointer to the counter registry.
std::ptrdiff_t fub::amrex::cutcell::IntegratorContext::GetCycles | ( | int | level = 0 | ) | const |
Returns the current number of cycles for data at the specified refinement level and direction.
::amrex::MultiFab& fub::amrex::cutcell::IntegratorContext::GetData | ( | int | level | ) |
Returns the MultiFab associated with level data on the specifed level number.
::amrex::MultiFab& fub::amrex::cutcell::IntegratorContext::GetDoublyShieldedFluxes | ( | int | level, |
Direction | dir | ||
) |
Returns the MultiFab associated with flux data on the specifed level number and direction.
|
noexcept |
Returns the refinement ratio in the specified direction.
::amrex::EBFArrayBoxFactory& fub::amrex::cutcell::IntegratorContext::GetEmbeddedBoundary | ( | int | level | ) |
Returns the current boundary condition for the specified level.
const ::amrex::EBFArrayBoxFactory& fub::amrex::cutcell::IntegratorContext::GetEmbeddedBoundary | ( | int | level | ) | const |
Returns the current boundary condition for the specified level.
|
noexcept |
Returns the refinement ratio in the specified direction.
::amrex::MultiFab& fub::amrex::cutcell::IntegratorContext::GetFluxes | ( | int | level, |
Direction | dir | ||
) |
Returns the MultiFab associated with flux data on the specifed level number and direction.
const ::amrex::Geometry& fub::amrex::cutcell::IntegratorContext::GetGeometry | ( | int | level | ) | const |
Returns the geometry object for the specified refinement level.
|
noexcept |
Returns a shared pointer to the underlying GriddingAlgorithm which owns the simulation.
|
noexcept |
Returns the hyperbolic method member object.
|
noexcept |
Returns the MPI communicator which is associated with this context.
|
noexcept |
Returns a reference to const PatchHierarchy which is a member of the GriddingAlgorithm.
|
noexcept |
Returns the hyperbolic method member object.
|
noexcept |
Returns the refinement ratio for all directions.
|
noexcept |
Returns the refinement ratio in the specified direction.
::amrex::MultiFab& fub::amrex::cutcell::IntegratorContext::GetReferenceStates | ( | int | level | ) |
Returns the MultiFab associated with flux data on the specifed level number and direction.
::amrex::MultiFab& fub::amrex::cutcell::IntegratorContext::GetScratch | ( | int | level | ) |
Returns the MultiFab associated with level data with ghost cells on the specifed level number and direction.
::amrex::MultiFab& fub::amrex::cutcell::IntegratorContext::GetShieldedFromLeftFluxes | ( | int | level, |
Direction | dir | ||
) |
Returns the MultiFab associated with flux data on the specifed level number and direction.
::amrex::MultiFab& fub::amrex::cutcell::IntegratorContext::GetShieldedFromRightFluxes | ( | int | level, |
Direction | dir | ||
) |
Returns the MultiFab associated with flux data on the specifed level number and direction.
::amrex::MultiFab& fub::amrex::cutcell::IntegratorContext::GetStabilizedFluxes | ( | int | level, |
Direction | dir | ||
) |
Returns the MultiFab associated with flux data on the specifed level number and direction.
Duration fub::amrex::cutcell::IntegratorContext::GetTimePoint | ( | int | level = 0 | ) | const |
Returns the current time level for data at the specified refinement level and direction.
|
noexcept |
Returns the refinement ratio in the specified direction.
IntegratorContext& fub::amrex::cutcell::IntegratorContext::operator= | ( | const IntegratorContext & | ) |
Deeply copies a context and all its distributed data for all MPI ranks.
|
defaultnoexcept |
Deeply copies a context and all its distributed data for all MPI ranks.
void fub::amrex::cutcell::IntegratorContext::PostAdvanceHierarchy | ( | ) |
Compute and store reference states for cut-cells over all split cycles.
Result<void, TimeStepTooLarge> fub::amrex::cutcell::IntegratorContext::PostAdvanceLevel | ( | int | level_num, |
Duration | dt, | ||
std::pair< int, int > | subcycle | ||
) |
Increases the internal time stamps and cycle counters for the specified level number and direction.
void fub::amrex::cutcell::IntegratorContext::PreAdvanceHierarchy | ( | ) |
Compute and store reference states for cut-cells over all split cycles.
int fub::amrex::cutcell::IntegratorContext::PreAdvanceLevel | ( | int | level_num, |
Duration | dt, | ||
std::pair< int, int > | subcycle | ||
) |
On each first subcycle this will regrid the data if neccessary.
void fub::amrex::cutcell::IntegratorContext::ResetCoarseFineFluxes | ( | int | fine, |
int | coarse | ||
) |
Resets all accumulates fluxes to zero.
void fub::amrex::cutcell::IntegratorContext::ResetHierarchyConfiguration | ( | int | level = 0 | ) |
Whenever the gridding algorithm changes the data hierarchy this function will regrid all distributed helper variables managed by the context.
[in] | level | The level number of the coarsest level which changed its shape. Regrid all levels finer than level. |
void fub::amrex::cutcell::IntegratorContext::ResetHierarchyConfiguration | ( | std::shared_ptr< GriddingAlgorithm > | gridding | ) |
Replaces the underlying gridding algorithm with the specified one.
void fub::amrex::cutcell::IntegratorContext::SetCycles | ( | std::ptrdiff_t | cycle, |
int | level | ||
) |
Sets the cycle count for a specific level number and direction.
|
inline |
Call this function upon every time integration.
void fub::amrex::cutcell::IntegratorContext::SetTimePoint | ( | Duration | t, |
int | level | ||
) |
Sets the time point for a specific level number and direction.
void fub::amrex::cutcell::IntegratorContext::UpdateConservatively | ( | int | level, |
Duration | dt, | ||
Direction | dir | ||
) |
Apply a conservative time update for each conservative variable on the specified level number and direction.
bool fub::amrex::cutcell::IntegratorContext::count_per_level |
|
private |
|
private |
|
private |
|
private |
|
private |
std::shared_ptr<CounterRegistry> fub::amrex::cutcell::IntegratorContext::registry_ |
|
private |