inciter namespace

Inciter declarations and definitions.

Namespaces

namespace cmd
Inciter command line grammar definition.
namespace ctr
Inciter control facilitating user input to internal data transfer.
namespace deck
Inciter input deck facilitating user input for computing shock hydrodynamics.

Classes

class ALECG
ALECG Charm++ chare array used to advance PDEs in time with ALECG+RK.
struct AUSM
class CGPDE
Partial differential equation base for continuous Galerkin PDEs.
class CmdLineParser
Command-line parser for Inciter.
class CompFlowProblemGaussHump
CompFlow system of PDEs problem: GaussHump.
class CompFlowProblemNLEnergyGrowth
class CompFlowProblemRayleighTaylor
class CompFlowProblemRotatedSodShocktube
class CompFlowProblemSedovBlastwave
CompFlow system of PDEs problem: Sedov blast-wave.
class CompFlowProblemSodShocktube
class CompFlowProblemTaylorGreen
class CompFlowProblemUserDefined
CompFlow system of PDEs problem: user defined.
class CompFlowProblemVorticalFlow
class DG
DG Charm++ chare array used to advance PDEs in time with DG+RK.
class DGPDE
Partial differential equation base for discontinuous Galerkin PDEs.
class DiagCG
DiagCG Charm++ chare array used to advance PDEs in time with DiagCG+LW+FCT.
class Discretization
Discretization Charm++ chare array holding common functinoality to all discretization schemes.
class DistFCT
DistFCT Charm++ chare array used to advance PDEs in time with DiagCG+LW+FCT.
class ElemDiagnostics
ElemDiagnostics class used to compute diagnostics while integrating PDEs.
class FaceData
FaceData class holding face-connectivity data useful for DG discretization.
class FluxCorrector
struct HLLC
class InciterDriver
Inciter driver used polymorphically with tk::Driver.
class InciterPrint
InciterPrint : tk::Print.
class InputDeckParser
Control file parser for Inciter.
struct LaxFriedrichs
class MultiMatProblemInterfaceAdvection
class MultiMatProblemSodShocktube
class MultiMatProblemUserDefined
MultiMat system of PDEs problem: user defined.
class MultiMatProblemWaterAirShocktube
class NodeDiagnostics
NodeDiagnostics class used to compute diagnostics while integrating PDEs.
class Partitioner
class PDEStack
Partial differential equations stack.
class Refiner
Mesh refiner for interfacing the mesh refinement library.
template<template<class, class> class Eq>
struct registerCG
template<template<class, class> class Eq>
struct registerDG
template<template<class, class> class Eq, class Factory, class PDE>
struct registerPDE
Function object for registering a partial differential equation into the partial differential equation factory.
struct registerRiemannSolver
Functor to register a Riemann solver into the Riemann solver factory.
class RiemannSolver
Generic Riemann solver interface class for various Riemann solvers.
class Scheme
Base class for generic forwarding interface to discretization proxies.
class Sorter
Mesh sorter for global distributed mesh node reordering.
class Transporter
Transporter drives the time integration of transport equations.
class TransportProblemCylAdvect
Transport PDE problem: advection of cylinder.
class TransportProblemGaussHump
Transport PDE problem: advection of two-dimensional Gaussian hump.
class TransportProblemShearDiff
class TransportProblemSlotCyl
struct Upwind
Upwind Riemann solver.

Enums

enum Diag { L2SOL =0, L2ERR, LINFERR, ITER, TIME, DT }
Diagnostics labels.
enum ProgMesh { PART =0, DIST, REFINE, BND, COMM, MASK, REORD }
Indices for progress report on mesh preparation.
enum ProgWork { CREATE =0, BNDFACE, COMFAC, GHOST, ADJ }
Indices for progress report on workers preparation.

Typedefs

using GhostData = std::unordered_map<std::size_t, std::tuple<std::vector<std::size_t>, std::vector<tk::real>, std::array<tk::real, 3>, std::size_t, std::array<std::size_t, 4>>>
Data associated to a tetrahedron cell id used to communicate across faces.
using CompFlowProblems = brigand::list<CompFlowProblemUserDefined, CompFlowProblemVorticalFlow, CompFlowProblemNLEnergyGrowth, CompFlowProblemRayleighTaylor, CompFlowProblemTaylorGreen, CompFlowProblemSodShocktube, CompFlowProblemRotatedSodShocktube, CompFlowProblemSedovBlastwave, CompFlowProblemGaussHump>
List of all CompFlow Problem policies (defined in the includes above)
using RiemannFactory = std::map<ctr::FluxType, std::function<RiemannSolver()>>
using MultiMatProblems = brigand::list<MultiMatProblemUserDefined, MultiMatProblemSodShocktube, MultiMatProblemInterfaceAdvection, MultiMatProblemWaterAirShocktube>
List of all MultiMat Problem policies (defined in the includes above)
using CGFactory = std::map<ctr::PDEKey, std::function<CGPDE(const tk::ctr::ncomp_t&)>>
Factory for PDEs using continuous Galerkin discretization storing keys associated to their constructors.
using DGFactory = std::map<ctr::PDEKey, std::function<DGPDE(const tk::ctr::ncomp_t&)>>
Factory for PDEs using discontinuous Galerkin discretization storing keys associated to their constructors.
using TransportProblems = brigand::list<TransportProblemShearDiff, TransportProblemSlotCyl, TransportProblemGaussHump, TransportProblemCylAdvect>
List of all Transport Problem policies (defined in the includes above)

Functions

auto serialize(const std::vector<std::vector<tk::real>>& d) -> std::pair<int, std::unique_ptr<char[]>>
Serialize std::vector to raw memory stream.
auto mergeDiag(int nmsg, CkReductionMsg** msgs) -> CkReductionMsg*
Charm++ custom reducer for merging std::vectors during reduction across PEs.
template<typename T>
auto serialize(const std::vector<std::pair<int, T>>& m) -> std::pair<int, std::unique_ptr<char[]>>
template<typename T>
auto mergeFields(int nmsg, CkReductionMsg** msgs) -> CkReductionMsg*
Charm++ custom reducer for merging mesh node indices categorized by chares during reduction across PEs.
auto match(] tk::ctr::ncomp_t ncomp, tk::real t, tk::real dt, const tk::UnsMesh::Coords& coord, const std::vector<std::size_t>& gid, const std::unordered_map<std::size_t, std::size_t>& lid, const std::map<int, std::vector<std::size_t>>& sidenodes) -> std::unordered_map<std::size_t, std::vector<std::pair<bool, tk::real>>>
auto correctBC(const tk::Fields& a, const tk::Fields& dul, const std::unordered_map<std::size_t, std::vector<std::pair<bool, tk::real>>>& bc, const std::unordered_map<std::size_t, std::size_t>& lid) -> bool
Verify that the change in the solution at those nodes where Dirichlet boundary conditions are set is exactly the amount the BCs prescribe.
auto match(tk::ctr::ncomp_t ncomp, tk::real t, tk::real dt, const tk::UnsMesh::Coords& coord, const std::vector<std::size_t>& gid, const std::unordered_map<std::size_t, std::size_t>& lid, const std::map<int, std::vector<std::size_t>>& sidenodes) -> std::unordered_map<std::size_t, std::vector<std::pair<bool, tk::real>>>
Match user-specified boundary conditions at nodes for side sets.
void operator|(PUP::er& p, std::vector<CGPDE>& eqs)
Pack/Unpack selected partial differential equations using continuous Galerkin discretization.
void operator|(PUP::er& p, std::vector<DGPDE>& eqs)
Pack/Unpack selected partial differential equations using discontinuous Galerkin discretization.
void registerCompFlow(CGFactory& cf, DGFactory& df, std::set<ctr::PDEType>& cgt, std::set<ctr::PDEType>& dgt)
Register compressible flow PDEs into PDE factory.
auto infoCompFlow(std::map<ctr::PDEType, tk::ctr::ncomp_t>& cnt) -> std::vector<std::pair<std::string, std::string>>
Return information on the compressible flow PDE.
void registerMultiMat(DGFactory& df, std::set<ctr::PDEType>& dgt)
Register compressible flow PDEs into PDE factory.
auto infoMultiMat(std::map<ctr::PDEType, tk::ctr::ncomp_t>& cnt) -> std::vector<std::pair<std::string, std::string>>
Return information on the multi-material compressible flow PDE.
void registerTransport(CGFactory& cf, DGFactory& df, std::set<ctr::PDEType>& cgt, std::set<ctr::PDEType>& dgt)
Register transport PDEs into PDE factory.
auto infoTransport(std::map<ctr::PDEType, tk::ctr::ncomp_t>& cnt) -> std::vector<std::pair<std::string, std::string>>
Return information on the transport PDE.
template<class Eq>
auto eos_density(ncomp_t system, tk::real pr, tk::real temp, std::size_t imat = 0) -> tk::real
Calculate density from the material pressure and temperature using the stiffened-gas equation of state.
template<class Eq>
auto eos_pressure(ncomp_t system, tk::real arho, tk::real u, tk::real v, tk::real w, tk::real arhoE, tk::real alpha = 1.0, std::size_t imat = 0) -> tk::real
Calculate pressure from the material density, momentum and total energy using the stiffened-gas equation of state.
template<class Eq>
auto eos_soundspeed(ncomp_t system, tk::real arho, tk::real apr, tk::real alpha = 1.0, std::size_t imat = 0) -> tk::real
template<class Eq>
auto eos_totalenergy(ncomp_t system, tk::real rho, tk::real u, tk::real v, tk::real w, tk::real pr, std::size_t imat = 0) -> tk::real
Calculate material specific total energy from the material density, momentum and material pressure.
auto RiemannSolvers() -> RiemannFactory
Register available Riemann solvers into a factory.
void WENO_P1(const std::vector<int>& esuel, inciter::ncomp_t offset, tk::Fields& U)
Weighted Essentially Non-Oscillatory (WENO) limiter for DGP1.
void WENOMultiMat_P1(const std::vector<int>& esuel, inciter::ncomp_t offset, tk::Fields& U, tk::Fields& P, std::size_t nmat)
Weighted Essentially Non-Oscillatory (WENO) limiter for multi-material DGP1.
void Superbee_P1(const std::vector<int>& esuel, const std::vector<std::size_t>& inpoel, const std::vector<std::size_t>& ndofel, inciter::ncomp_t offset, const tk::UnsMesh::Coords& coord, tk::Fields& U)
Superbee limiter for DGP1.
void SuperbeeMultiMat_P1(const std::vector<int>& esuel, const std::vector<std::size_t>& inpoel, const std::vector<std::size_t>& ndofel, inciter::ncomp_t offset, const tk::UnsMesh::Coords& coord, tk::Fields& U, tk::Fields& P, std::size_t nmat)
Superbee limiter for multi-material DGP1.
void WENOFunction(const tk::Fields& U, const std::vector<int>& esuel, std::size_t e, inciter::ncomp_t c, std::size_t rdof, inciter::ncomp_t offset, tk::real cweight, std::array<std::vector<tk::real>, 3>& limU)
WENO limiter function calculation for P1 dofs.
auto SuperbeeFunction(const tk::Fields& U, const std::vector<int>& esuel, const std::vector<std::size_t>& inpoel, const tk::UnsMesh::Coords& coord, std::size_t e, std::size_t ndof, std::size_t rdof, std::size_t dof_el, inciter::ncomp_t offset, inciter::ncomp_t ncomp, tk::real beta_lim) -> std::vector<tk::real>
Superbee limiter function calculation for P1 dofs.
void consistentMultiMatLimiting_P1(std::size_t nmat, ncomp_t offset, std::size_t rdof, std::size_t e, tk::Fields& U, ] tk::Fields& P, std::vector<tk::real>& phic, ] std::vector<tk::real>& phip)
void consistentMultiMatLimiting_P1(std::size_t nmat, ncomp_t offset, std::size_t rdof, std::size_t e, tk::Fields& U, tk::Fields& P, std::vector<tk::real>& phic, std::vector<tk::real>& phip)
Consistent limiter modifications for P1 dofs.
auto volfracIdx(std::size_t, std::size_t kmat) -> std::size_t
auto densityIdx(std::size_t nmat, std::size_t kmat) -> std::size_t
auto momentumIdx(std::size_t nmat, std::size_t idir) -> std::size_t
auto energyIdx(std::size_t nmat, std::size_t kmat) -> std::size_t
auto velocityIdx(std::size_t nmat, std::size_t idir) -> std::size_t
auto pressureIdx(std::size_t, std::size_t kmat) -> std::size_t
auto volfracDofIdx(std::size_t nmat, std::size_t kmat, std::size_t ndof, std::size_t idof) -> std::size_t
Get the index of the required DOF of material volume fraction from the DG solution vector.
auto densityDofIdx(std::size_t nmat, std::size_t kmat, std::size_t ndof, std::size_t idof) -> std::size_t
Get the index of the required DOF of material continuity equation from the DG solution vector.
auto momentumDofIdx(std::size_t nmat, std::size_t idir, std::size_t ndof, std::size_t idof) -> std::size_t
Get the index of the required DOF of momentum equation component from the DG solution vector.
auto energyDofIdx(std::size_t nmat, std::size_t kmat, std::size_t ndof, std::size_t idof) -> std::size_t
Get the index of the required DOF of material total energy equation from the DG solution vector.
auto velocityDofIdx(std::size_t nmat, std::size_t idir, std::size_t ndof, std::size_t idof) -> std::size_t
Get the index of the required DOF of velocity component from the DG vector of primitives.
auto pressureDofIdx(std::size_t nmat, std::size_t kmat, std::size_t ndof, std::size_t idof) -> std::size_t
Get the index of the required DOF of material pressure from the DG vector of primitives.
template<typename V>
auto parameters(const V& v) -> std::string
Convert and return values from vector as string.
static void spectral_decay(std::size_t nunk, const std::vector<int>& esuel, const tk::Fields& unk, std::size_t ndof, std::size_t ndofmax, tk::real tolref, std::vector<std::size_t>& ndofel)
static void non_conformity(std::size_t nunk, std::size_t Nbfac, const std::vector<std::size_t>& inpoel, const tk::UnsMesh::Coords& coord, const std::vector<int>& esuel, const std::vector<int>& esuf, const std::vector<std::size_t>& inpofa, const tk::Fields& unk, std::size_t ndof, std::size_t ndofmax, std::vector<std::size_t>& ndofel)
void eval_ndof(std::size_t nunk, const tk::UnsMesh::Coords& coord, const std::vector<std::size_t>& inpoel, const inciter::FaceData& fd, const tk::Fields& unk, inciter::ctr::PrefIndicatorType indicator, std::size_t ndof, std::size_t ndofmax, tk::real tolref, std::vector<std::size_t>& ndofel)
Evaluate the adaptive indicator and mark the ndof for each element.

Variables

ctr::InputDeck g_inputdeck
ctr::InputDeck g_inputdeck_defaults
std::vector<CGPDE> g_cgpde
std::vector<DGPDE> g_dgpde
static const std::array<std::array<tk::real, 3>, 2> rkcoef
Runge-Kutta coefficients.
const std::size_t NUMDIAG
Number of entries in diagnostics vector (of vectors)
static const std::array<std::string, 7> ProgMeshPrefix
Prefixes for progress report on mesh preparation.
static const std::array<std::string, 5> ProgWorkPrefix
Prefixes for progress report on workers preparation.

Enum documentation

enum inciter::Diag

Diagnostics labels.

Enumerators
L2SOL

L2 norm of numerical solution.

L2ERR

L2 norm of numerical-analytic solution.

LINFERR

L_inf norm of numerical-analytic solution.

ITER

Iteration count.

TIME

Physical time.

DT

Typedef documentation

using inciter::RiemannFactory = std::map<ctr::FluxType, std::function<RiemannSolver()>>

Factory for Riemann solvers

This factory is used to store the constructors as a std::function of specific Riemann solvers that can be invoked at a later time compared to the point where the map is populated. The key is an enum, uniquely idenfitying a specific Riemann solver. The value is std::function storing a constructor to be invoked. The type of object stored in std::function is a generic (base) class constructor, which provides a polymorphyic interface (overridable functions) that specific (child) Riemann solvers override, yielding runtime polymorphism.

Function documentation

std::pair<int, std::unique_ptr<char[]>> inciter::serialize(const std::vector<std::vector<tk::real>>& d)

Serialize std::vector to raw memory stream.

Parameters
in Diagnostics vector of vectors (of eq components)
Returns Pair of the length and the raw stream containing the serialized vectors

CkReductionMsg* inciter::mergeDiag(int nmsg, CkReductionMsg** msgs)

Charm++ custom reducer for merging std::vectors during reduction across PEs.

Parameters
nmsg in Number of messages in msgs
msgs in Charm++ reduction message containing the serialized diagnostics
Returns Aggregated diagnostics built for further aggregation if needed

template<typename T>
std::pair<int, std::unique_ptr<char[]>> inciter::serialize(const std::vector<std::pair<int, T>>& m)

Parameters
in Chare mesh node indices to serialize
Returns Pair of the length and the raw stream containing the serialized data

Serialize mesh node indices categorized by chares to raw memory stream

template<typename T>
CkReductionMsg* inciter::mergeFields(int nmsg, CkReductionMsg** msgs)

Charm++ custom reducer for merging mesh node indices categorized by chares during reduction across PEs.

Parameters
nmsg in Number of messages in msgs
msgs in Charm++ reduction message containing the serialized mesh node indices categorized by chares
Returns Aggregated mesh node indices categorized by chares built for further aggregation if needed

During aggregation the outer vector is simply concatenated.

std::unordered_map<std::size_t, std::vector<std::pair<bool, tk::real>>> inciter::match(] tk::ctr::ncomp_t ncomp, tk::real t, tk::real dt, const tk::UnsMesh::Coords& coord, const std::vector<std::size_t>& gid, const std::unordered_map<std::size_t, std::size_t>& lid, const std::map<int, std::vector<std::size_t>>& sidenodes)

Parameters
ncomp in Number of scalar components in PDE system
in Physical time at which to query boundary conditions
dt in Time step size (for querying BC increments in time)
coord in Mesh node coordinates
gid in Global node IDs
lid in Local node IDs associated to global node IDs
sidenodes in Map storing global mesh node IDs mapped to side set ids
Returns Vector of pairs of bool and boundary condition value associated to mesh node IDs at which the user has set Dirichlet boundary conditions for all systems of PDEs integrated. The bool indicates whether the BC is set at the node for that component: if true, the real value is the increment (from t to dt) in the BC specified for a component.

Boundary conditions (BC), mathematically speaking, are applied on finite surfaces. These finite surfaces are given by element sets (i.e., a list of elements). This function queries Dirichlet boundary condition values from all PDEs in the system of systems of PDEs integrated at the node lists associated to side set IDs, given by sidenodes. Each PDE system returns a BC data structure. Note that the BC mesh nodes that this function results in, stored in dirbc, only contains those nodes that are supplied via sidenodes, i.e., in parallel only a part of the mesh is worked on.

bool inciter::correctBC(const tk::Fields& a, const tk::Fields& dul, const std::unordered_map<std::size_t, std::vector<std::pair<bool, tk::real>>>& bc, const std::unordered_map<std::size_t, std::size_t>& lid)

Verify that the change in the solution at those nodes where Dirichlet boundary conditions are set is exactly the amount the BCs prescribe.

Parameters
in Limited antidiffusive element contributions (from FCT)
dul in Low order solution increment
bc in Vector of boundary conditions (true/false + BC value) for all scalar components integrated associated of all systems to global node ID
lid in Local node IDs associated to global node IDs
Returns True if solution is correct at Dirichlet boundary condition nodes

We loop through the map that associates a vector of of boundary conditions (true/false, indicating whether the BC is set + BC value if true) for all scalar components integrated associated of all systems to global node IDs. Then for all scalar components of all systems of systems of PDEs integrated if a BC is to be set for a given component, we compute the low order solution increment + the anti-diffusive element contributions (in FCT), which is the current solution increment (to be used to update the solution at time n in FCT) at that node. This solution increment must equal the BC prescribed at the given node as we solve for solution increments. If not, the BCs are not set correctly, which is an error.

void inciter::operator|(PUP::er& p, std::vector<CGPDE>& eqs)

Pack/Unpack selected partial differential equations using continuous Galerkin discretization.

This Pack/Unpack method (re-)creates the PDE factory since it needs to (re-)bind function pointers on different processing elements. Therefore we circumvent Charm's usual pack/unpack for this type, and thus sizing does not make sense: sizing is a no-op. We could initialize the factory in InciterDriver's constructor and let this function re-create the stack only when unpacking, but that leads to repeating the same code twice: once in InciterDriver's constructor, once here. Another option is to use this pack/unpack routine to both initially create (when packing) and to re-create (when unpacking) the factory, which eliminates the need for pre-creating the object in InciterDriver's constructor and therefore eliminates the repeated code. This explains the guard for sizing: the code below is called for packing only (in serial) and packing and unpacking (in parallel).

void inciter::operator|(PUP::er& p, std::vector<DGPDE>& eqs)

Pack/Unpack selected partial differential equations using discontinuous Galerkin discretization.

This Pack/Unpack method (re-)creates the PDE factory since it needs to (re-)bind function pointers on different processing elements. Therefore we circumvent Charm's usual pack/unpack for this type, and thus sizing does not make sense: sizing is a no-op. We could initialize the factory in InciterDriver's constructor and let this function re-create the stack only when unpacking, but that leads to repeating the same code twice: once in InciterDriver's constructor, once here. Another option is to use this pack/unpack routine to both initially create (when packing) and to re-create (when unpacking) the factory, which eliminates the need for pre-creating the object in InciterDriver's constructor and therefore eliminates the repeated code. This explains the guard for sizing: the code below is called for packing only (in serial) and packing and unpacking (in parallel).

void inciter::registerCompFlow(CGFactory& cf, DGFactory& df, std::set<ctr::PDEType>& cgt, std::set<ctr::PDEType>& dgt)

Register compressible flow PDEs into PDE factory.

Parameters
cf in/out Continuous Galerkin PDE factory to register to
df in/out Discontinuous Galerkin PDE factory to register to
cgt in/out Counters for equation types registered into CG factory
dgt in/out Counters for equation types registered into DG factory

std::vector<std::pair<std::string, std::string>> inciter::infoCompFlow(std::map<ctr::PDEType, tk::ctr::ncomp_t>& cnt)

Return information on the compressible flow PDE.

Parameters
cnt in/out std::map of counters for all PDE types
Returns vector of string pairs describing the PDE configuration

void inciter::registerMultiMat(DGFactory& df, std::set<ctr::PDEType>& dgt)

Register compressible flow PDEs into PDE factory.

Parameters
df in/out Discontinuous Galerkin PDE factory to register to
dgt in/out Counters for equation types registered into DG factory

std::vector<std::pair<std::string, std::string>> inciter::infoMultiMat(std::map<ctr::PDEType, tk::ctr::ncomp_t>& cnt)

Return information on the multi-material compressible flow PDE.

Parameters
cnt in/out std::map of counters for all PDE types
Returns vector of string pairs describing the PDE configuration

void inciter::registerTransport(CGFactory& cf, DGFactory& df, std::set<ctr::PDEType>& cgt, std::set<ctr::PDEType>& dgt)

Register transport PDEs into PDE factory.

Parameters
cf in/out Continuous Galerkin PDE factory to register to
df in/out Discontinuous Galerkin PDE factory to register to
cgt in/out Counters for equation types registered into CG factory
dgt in/out Counters for equation types registered into DG factory

std::vector<std::pair<std::string, std::string>> inciter::infoTransport(std::map<ctr::PDEType, tk::ctr::ncomp_t>& cnt)

Return information on the transport PDE.

Parameters
cnt in/out std::map of counters for all PDE types
Returns vector of string pairs describing the PDE configuration

template<class Eq>
tk::real inciter::eos_density(ncomp_t system, tk::real pr, tk::real temp, std::size_t imat = 0)

Calculate density from the material pressure and temperature using the stiffened-gas equation of state.

Template parameters
Eq Equation type to operate on, e.g., tag::compflow, tag::multimat
Parameters
system in Equation system index
pr in Material pressure
temp in Material temperature
imat in Material-id who's EoS is required. Default is 0, so that for the single-material system, this argument can be left unspecified by the calling code
Returns Material density calculated using the stiffened-gas EoS

template<class Eq>
tk::real inciter::eos_pressure(ncomp_t system, tk::real arho, tk::real u, tk::real v, tk::real w, tk::real arhoE, tk::real alpha = 1.0, std::size_t imat = 0)

Calculate pressure from the material density, momentum and total energy using the stiffened-gas equation of state.

Template parameters
Eq Equation type to operate on, e.g., tag::compflow, tag::multimat
Parameters
system in Equation system index
arho in Material partial density (alpha_k * rho_k)
in X-velocity
in Y-velocity
in Z-velocity
arhoE in Material total energy (alpha_k * rho_k * E_k)
alpha in Material volume fraction. Default is 1.0, so that for the single-material system, this argument can be left unspecified by the calling code
imat in Material-id who's EoS is required. Default is 0, so that for the single-material system, this argument can be left unspecified by the calling code
Returns Material partial pressure (alpha_k * p_k) calculated using the stiffened-gas EoS

template<class Eq>
tk::real inciter::eos_soundspeed(ncomp_t system, tk::real arho, tk::real apr, tk::real alpha = 1.0, std::size_t imat = 0)

Template parameters
Eq Equation type to operate on, e.g., tag::compflow, tag::multimat
Parameters
system in Equation system index
arho in Material partial density (alpha_k * rho_k)
apr in Material partial pressure (alpha_k * p_k)
alpha in Material volume fraction. Default is 1.0, so that for the single-material system, this argument can be left unspecified by the calling code
imat in Material-id who's EoS is required. Default is 0, so that for the single-material system, this argument can be left unspecified by the calling code
Returns Material speed of sound using the stiffened-gas EoS

Calculate speed of sound from the material density and material pressure

template<class Eq>
tk::real inciter::eos_totalenergy(ncomp_t system, tk::real rho, tk::real u, tk::real v, tk::real w, tk::real pr, std::size_t imat = 0)

Calculate material specific total energy from the material density, momentum and material pressure.

Template parameters
Eq Equation type to operate on, e.g., tag::compflow, tag::multimat
Parameters
system in Equation system index
rho in Material density
in X-velocity
in Y-velocity
in Z-velocity
pr in Material pressure
imat in Material-id who's EoS is required. Default is 0, so that for the single-material system, this argument can be left unspecified by the calling code
Returns Material specific total energy using the stiffened-gas EoS

RiemannFactory inciter::RiemannSolvers()

Register available Riemann solvers into a factory.

Returns Riemann solver factory

void inciter::WENO_P1(const std::vector<int>& esuel, inciter::ncomp_t offset, tk::Fields& U)

Weighted Essentially Non-Oscillatory (WENO) limiter for DGP1.

Parameters
esuel in Elements surrounding elements
offset in Index for equation systems
in/out High-order solution vector which gets limited

This WENO function should be called for transport and compflow

void inciter::WENOMultiMat_P1(const std::vector<int>& esuel, inciter::ncomp_t offset, tk::Fields& U, tk::Fields& P, std::size_t nmat)

Weighted Essentially Non-Oscillatory (WENO) limiter for multi-material DGP1.

Parameters
esuel in Elements surrounding elements
offset in Index for equation systems
in/out High-order solution vector which gets limited
in/out High-order vector of primitives which gets limited
nmat in Number of materials in this PDE system

This WENO function should be called for multimat

void inciter::Superbee_P1(const std::vector<int>& esuel, const std::vector<std::size_t>& inpoel, const std::vector<std::size_t>& ndofel, inciter::ncomp_t offset, const tk::UnsMesh::Coords& coord, tk::Fields& U)

Superbee limiter for DGP1.

Parameters
esuel in Elements surrounding elements
inpoel in Element connectivity
ndofel in Vector of local number of degrees of freedom
offset in Index for equation systems
coord in Array of nodal coordinates
in/out High-order solution vector which gets limited

This Superbee function should be called for transport and compflow

void inciter::SuperbeeMultiMat_P1(const std::vector<int>& esuel, const std::vector<std::size_t>& inpoel, const std::vector<std::size_t>& ndofel, inciter::ncomp_t offset, const tk::UnsMesh::Coords& coord, tk::Fields& U, tk::Fields& P, std::size_t nmat)

Superbee limiter for multi-material DGP1.

Parameters
esuel in Elements surrounding elements
inpoel in Element connectivity
ndofel in Vector of local number of degrees of freedom
offset in Index for equation systems
coord in Array of nodal coordinates
in/out High-order solution vector which gets limited
in/out High-order vector of primitives which gets limited
nmat in Number of materials in this PDE system

This Superbee function should be called for multimat

void inciter::WENOFunction(const tk::Fields& U, const std::vector<int>& esuel, std::size_t e, inciter::ncomp_t c, std::size_t rdof, inciter::ncomp_t offset, tk::real cweight, std::array<std::vector<tk::real>, 3>& limU)

WENO limiter function calculation for P1 dofs.

Parameters
in High-order solution vector which is to be limited
esuel in Elements surrounding elements
in Id of element whose solution is to be limited
in Index of component which is to be limited
rdof in Maximum number of reconstructed degrees of freedom
offset in Index for equation systems
cweight in Weight of the central stencil
limU in/out Limited gradients of component c

std::vector<tk::real> inciter::SuperbeeFunction(const tk::Fields& U, const std::vector<int>& esuel, const std::vector<std::size_t>& inpoel, const tk::UnsMesh::Coords& coord, std::size_t e, std::size_t ndof, std::size_t rdof, std::size_t dof_el, inciter::ncomp_t offset, inciter::ncomp_t ncomp, tk::real beta_lim)

Superbee limiter function calculation for P1 dofs.

Parameters
in High-order solution vector which is to be limited
esuel in Elements surrounding elements
inpoel in Element connectivity
coord in Array of nodal coordinates
in Id of element whose solution is to be limited
ndof in Maximum number of degrees of freedom
rdof in Maximum number of reconstructed degrees of freedom
dof_el in Local number of degrees of freedom
offset in Index for equation systems
ncomp in Number of scalar components in this PDE system
beta_lim in Parameter which is equal to 2 for Superbee and 1 for minmod limiter
Returns phi Limiter function for solution in element e

void inciter::consistentMultiMatLimiting_P1(std::size_t nmat, ncomp_t offset, std::size_t rdof, std::size_t e, tk::Fields& U, ] tk::Fields& P, std::vector<tk::real>& phic, ] std::vector<tk::real>& phip)

Parameters
nmat in Number of materials in this PDE system
offset in Index for equation system
rdof in Total number of reconstructed dofs
in Element being checked for consistency
in/out Second-order solution vector which gets modified near material interfaces for consistency
in/out Second-order vector of primitive quantities which gets modified near material interfaces for consistency
phic in/out Vector of limiter functions for the conserved quantities
phip in/out Vector of limiter functions for the primitive quantities

std::size_t inciter::volfracIdx(std::size_t, std::size_t kmat)

Parameters
kmat in Index of required material
Returns Index of the required material volume fraction

Get the index of the required material volume fraction

std::size_t inciter::densityIdx(std::size_t nmat, std::size_t kmat)

Parameters
nmat in Number of materials
kmat in Index of required material
Returns Index of the required material continuity equation

Get the index of the required material continuity equation

std::size_t inciter::momentumIdx(std::size_t nmat, std::size_t idir)

Parameters
nmat in Number of materials
idir in Required component direction; 0: X-component, 1: Y-component, 2: Z-component.
Returns Index of the required momentum equation component

Get the index of the required momentum equation component

std::size_t inciter::energyIdx(std::size_t nmat, std::size_t kmat)

Parameters
nmat in Number of materials
kmat in Index of required material
Returns Index of the required material total energy equation

Get the index of the required material total energy equation

std::size_t inciter::velocityIdx(std::size_t nmat, std::size_t idir)

Parameters
nmat in Number of materials
idir in Required component direction; 0: X-component, 1: Y-component, 2: Z-component.
Returns Index of the required velocity component from vector of primitives

Get the index of the required velocity component from vector of primitives

std::size_t inciter::pressureIdx(std::size_t, std::size_t kmat)

Parameters
kmat in Index of required material
Returns Index of the required material pressure from vector of primitives

Get the index of the required material pressure from vector of primitives

std::size_t inciter::volfracDofIdx(std::size_t nmat, std::size_t kmat, std::size_t ndof, std::size_t idof)

Get the index of the required DOF of material volume fraction from the DG solution vector.

Parameters
nmat in Number of materials
kmat in Index of required material
ndof in Number of solution DOFs stored in DG solution vector
idof in Index of required solution DOF from DG solution vector
Returns Index of the required material volume fraction

This function is used to get the index of the required DOF in the solution vector, which is of type tk::Fields.

std::size_t inciter::densityDofIdx(std::size_t nmat, std::size_t kmat, std::size_t ndof, std::size_t idof)

Get the index of the required DOF of material continuity equation from the DG solution vector.

Parameters
nmat in Number of materials
kmat in Index of required material
ndof in Number of solution DOFs stored in DG solution vector
idof in Index of required solution DOF from DG solution vector
Returns Index of the required material continuity equation

This function is used to get the index of the required DOF in the solution vector, which is of type tk::Fields.

std::size_t inciter::momentumDofIdx(std::size_t nmat, std::size_t idir, std::size_t ndof, std::size_t idof)

Get the index of the required DOF of momentum equation component from the DG solution vector.

Parameters
nmat in Number of materials
idir in Required component direction; 0: X-component, 1: Y-component, 2: Z-component.
ndof in Number of solution DOFs stored in DG solution vector
idof in Index of required solution DOF from DG solution vector
Returns Index of the required momentum equation component

This function is used to get the index of the required DOF in the solution vector, which is of type tk::Fields.

std::size_t inciter::energyDofIdx(std::size_t nmat, std::size_t kmat, std::size_t ndof, std::size_t idof)

Get the index of the required DOF of material total energy equation from the DG solution vector.

Parameters
nmat in Number of materials
kmat in Index of required material
ndof in Number of solution DOFs stored in DG solution vector
idof in Index of required solution DOF from DG solution vector
Returns Index of the required material total energy equation

This function is used to get the index of the required DOF in the solution vector, which is of type tk::Fields.

std::size_t inciter::velocityDofIdx(std::size_t nmat, std::size_t idir, std::size_t ndof, std::size_t idof)

Get the index of the required DOF of velocity component from the DG vector of primitives.

Parameters
nmat in Number of materials
idir in Required component direction; 0: X-component, 1: Y-component, 2: Z-component.
ndof in Number of solution DOFs stored in DG solution vector
idof in Index of required solution DOF from DG solution vector
Returns Index of the required velocity component from vector of primitives

This function is used to get the index of the required DOF in the solution vector, which is of type tk::Fields.

std::size_t inciter::pressureDofIdx(std::size_t nmat, std::size_t kmat, std::size_t ndof, std::size_t idof)

Get the index of the required DOF of material pressure from the DG vector of primitives.

Parameters
nmat in Number of materials
kmat in Index of required material
ndof in Number of solution DOFs stored in DG solution vector
idof in Index of required solution DOF from DG solution vector
Returns Index of the required material pressure from vector of primitives

This function is used to get the index of the required DOF in the solution vector, which is of type tk::Fields.

template<typename V>
std::string inciter::parameters(const V& v)

Convert and return values from vector as string.

Parameters
in Vector whose components to return as a string
Returns Concatenated string of values read from a vector

static void inciter::spectral_decay(std::size_t nunk, const std::vector<int>& esuel, const tk::Fields& unk, std::size_t ndof, std::size_t ndofmax, tk::real tolref, std::vector<std::size_t>& ndofel)

Parameters
nunk in Number of unknowns
esuel in Elements surrounding elements
unk in Array of unknowns
ndof in Number of degrees of freedom in the solution
ndofmax in Max number of degrees of freedom for p-refinement
tolref in Tolerance for p-refinement
ndofel in/out Vector of local number of degrees of freedome

Evaluate the spectral-decay indicator and mark the ndof for each element The spectral decay indicator, implemented in this functiopn, calculates the difference between the projections of the numerical solutions on finite element space of order p and p-1.

static void inciter::non_conformity(std::size_t nunk, std::size_t Nbfac, const std::vector<std::size_t>& inpoel, const tk::UnsMesh::Coords& coord, const std::vector<int>& esuel, const std::vector<int>& esuf, const std::vector<std::size_t>& inpofa, const tk::Fields& unk, std::size_t ndof, std::size_t ndofmax, std::vector<std::size_t>& ndofel)

Parameters
nunk in Number of unknowns
Nbfac in Number of internal faces
inpoel in Element-node connectivity
coord in Array of nodal coordinates
esuel in Elements surrounding elements
esuf
inpofa in Face-node connectivity
unk in Array of unknowns
ndof in Number of degrees of freedom in the solution
ndofmax in Max number of degrees of freedom for p-refinement
ndofel in/out Vector of local number of degrees of freedome

Evaluate the non-conformity indicator and mark the ndof for each element The non-conformity indicator, this function implements, evaluates the jump in the numerical solutions as a measure of the numerical error.

void inciter::eval_ndof(std::size_t nunk, const tk::UnsMesh::Coords& coord, const std::vector<std::size_t>& inpoel, const inciter::FaceData& fd, const tk::Fields& unk, inciter::ctr::PrefIndicatorType indicator, std::size_t ndof, std::size_t ndofmax, tk::real tolref, std::vector<std::size_t>& ndofel)

Evaluate the adaptive indicator and mark the ndof for each element.

Parameters
nunk in Number of unknowns
coord in Array of nodal coordinates
inpoel in Element-node connectivity
fd in Face connectivity and boundary conditions object
unk in Array of unknowns
indicator in p-refinement indicator type
ndof in Number of degrees of freedom in the solution
ndofmax in Max number of degrees of freedom for p-refinement
tolref in Tolerance for p-refinement
ndofel in/out Vector of local number of degrees of freedome

Evaluate the adaptive indicator and mark the ndof for each element

Variable documentation

ctr::InputDeck inciter::g_inputdeck

Input deck filled by parser, containing all input data

This object is in global scope, it contains all of user input, and thus it is made available to all PEs for convenience reasons. The runtime system distributes it to all PEs during initialization. Once distributed, the object does not change.

ctr::InputDeck inciter::g_inputdeck_defaults

Global-scope data. Initialized by the main chare and distibuted to all PEs by the Charm++ runtime system. Though semantically not const, all these global data should be considered read-only. See also http://charm.cs.illinois.edu/manuals/html/charm++/manual.html. The data below is global-scope because they must be available to all PEs which could be on different machines. Defaults of input deck, facilitates detection what is set by user

This object is in global scope, it contains the default of all possible user input, and thus it is made available to all PEs for convenience reasons. The runtime system distributes it to all PEs during initialization. Once distributed, the object does not change.

std::vector<CGPDE> inciter::g_cgpde

Partial differential equations using continuous Galerkin selected by user

This vector is in global scope, because it holds polymorphic objects, and thus must be distributed to all PEs during initialization. Once distributed by the runtime system, the objects do not change.

std::vector<DGPDE> inciter::g_dgpde

Partial differential equations using discontinuous Galerkin selected by user

This vector is in global scope, because it holds polymorphic objects, and thus must be distributed to all PEs during initialization. Once distributed by the runtime system, the objects do not change.