Peano 4
Loading...
Searching...
No Matches
exahype2::aderdg Namespace Reference

Namespaces

namespace  tests
 

Functions

tarch::la::Vector< 2, double > getQuadraturePoint (const tarch::la::Vector< 2, double > &cellCentre, const tarch::la::Vector< 2, double > &cellSize, const tarch::la::Vector< 2, int > &index, int polynomialOrder, const double *__restrict__ quadraturePoints)
 In ExaHyPE's Finite Volume setup, a cell hosts a patch of Finite Volumes.
 
tarch::la::Vector< 3, double > getQuadraturePoint (const tarch::la::Vector< 3, double > &cellCentre, const tarch::la::Vector< 3, double > &cellSize, const tarch::la::Vector< 3, int > &index, int polynomialOrder, const double *__restrict__ quadraturePoints)
 
template<typename T >
void computeCellErrorIntegral (std::function< void(const tarch::la::Vector< Dimensions, double > &position, const double t, const double dt, double *sol) > exactSolution, const tarch::la::Vector< Dimensions, double > &cellCentre, const tarch::la::Vector< Dimensions, double > &cellSize, double t, double dt, int Order, const double *__restrict__ quadratureWeights, const double *__restrict__ quadraturePoints, const int unknowns, const int auxiliary_variables, T *__restrict__ Q, double errors[3])
 
double corrector_adjustSolution_computeMaxEigenvalue_body_AoS (std::function< void(double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t) > adjustSolution, std::function< double(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, const int normal) > maxAbsoluteEigenvalue, double *__restrict__ UOut, const double *const __restrict__ nodes, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const int nodesPerAxis, const int unknowns, const int strideQ, const int callMaxEigenvalue, const int scalarIndex)
 Allow user to modify solution at the given coordinates and time and compute the max eigenvalue afterwards.
 
void corrector_addFluxContributions_body_AoS (std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux, double *__restrict__ UOut, const double *__restrict__ QIn, double *__restrict__ FAux, const double *__restrict__ nodes, const double *__restrict__ weights, const double *__restrict__ Kxi, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const double dt, const int nodesPerAxis, const int unknowns, const int strideQ, const int scalarIndex)
 Add space-time volume flux contributions to the solution.
 
void corrector_addSourceContributions_body_AoS (std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, double *__restrict__ S) > algebraicSource, double *__restrict__ UOut, double *__restrict__ SAux, const double *__restrict__ QIn, const double *__restrict__ nodes, const double *__restrict__ weights, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const double dt, const int nodesPerAxis, const int unknowns, const int strideQ, const int scalarIndex)
 Add source contributions to the solution.
 
void corrector_addNcpContributions_body_AoS (std::function< void(const double *const __restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct, double *__restrict__ UOut, double *__restrict__ gradQAux, double *__restrict__ SAux, const double *__restrict__ QIn, const double *__restrict__ nodes, const double *__restrict__ weights, const double *__restrict__ dudx, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const double dt, const int nodesPerAxis, const int unknowns, const int strideQ, const int scalarIndex)
 Add nonconservative product contributions to the solution.
 
GPUCallableMethod void corrector_addRiemannContributions_body_AoS (double *__restrict__ UOut, const double *const __restrict__ riemannResultIn, const double *const __restrict__ weights, const double *const __restrict__ FLCoeff, const double dx, const double dt, const int nodesPerAxis, const int unknowns, const int strideQ, const int strideF, const int scalarIndexCell)
 Add Riemann flux contributions to the solution.
 
void corrector_addCellContributions_loop_AoS (std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux, std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, double *__restrict__ S) > algebraicSource, std::function< void(const double *const __restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct, double *__restrict__ UOut, const double *const __restrict__ QIn, const double *const __restrict__ nodes, const double *const __restrict__ weights, const double *const __restrict__ Kxi, const double *const __restrict__ dudx, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const double dt, const int order, const int unknowns, const int auxiliaryVariables, const bool callFlux, const bool callSource, const bool callNonconservativeProduct)
 Add cell-local contributions to new solution or update vector.
 
double corrector_addRiemannContributions_loop_AoS (std::function< void(double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t) > adjustSolution, std::function< double(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, const int normal) > maxAbsoluteEigenvalue, double *__restrict__ UOut, const double *const __restrict__ riemannResultIn, const double *const __restrict__ nodes, const double *const __restrict__ weights, const double *const __restrict__ FLCoeff, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const double dt, const int order, const int unknowns, const int auxiliaryVariables, const bool callMaxEigenvalue)
 Add cell-local contributions to new solution or update vector.
 
std::string plotNode (double Q[], int unknowns)
 Render a DG nodes components as a string.
 
void validatePatch (const double *__restrict__ UIn, int unknowns, int auxiliaryVariables, int order, const std::string &location="")
 Just runs over the patch and ensures that no entry is non or infinite.
 
void validateSpacetimePatch (const double *__restrict__ QIn, int unknowns, int auxiliaryVariables, int order, const std::string &location="")
 Just runs over the space-time patch and ensures that no entry is non or infinite.
 
GPUCallableMethod int getNodesPerCell (int nodesPerAxis)
 
GPUCallableMethod int getSpaceTimeNodesPerCell (int nodesPerAxis)
 
GPUCallableMethod tarch::la::Vector< Dimensions+1, intgetStrides (int nodesPerAxis, bool strides4SpaceTimeQuantity=true)
 If the stride is used for a space-time quantity, then this function generates:
 
GPUCallableMethod int lineariseIndex (const tarch::la::Vector< Dimensions+1, int > &index, const tarch::la::Vector< Dimensions+1, int > &strides)
 
GPUCallableMethod tarch::la::Vector< Dimensions+1, intdelineariseIndex (int scalarIndex, const tarch::la::Vector< Dimensions+1, int > &strides)
 
GPUCallableMethod tarch::la::Vector< Dimensions+1, double > getCoordinates (const tarch::la::Vector< Dimensions+1, int > &index, const tarch::la::Vector< Dimensions, double > &centre, const double dx, const double t, const double dt, const double *__restrict__ nodes)
 Compute coordinates from cell geometry and quadrature index.
 
GPUCallableMethod tarch::la::Vector< Dimensions+1, double > getCoordinatesOnFace (const tarch::la::Vector< Dimensions+1, int > &indexOnFace, const tarch::la::Vector< Dimensions, double > &faceCentre, const int direction, const double dx, const double t, const double dt, const double *__restrict__ nodes)
 Compute coordinates from cell geometry and quadrature index.
 
GPUCallableMethod int mapCellIndexToScalarFaceIndex (const tarch::la::Vector< Dimensions+1, int > &indexCell, const int direction, const int nodesPerAxis)
 Map cell index to an index for a single face.
 
GPUCallableMethod int mapCellIndexToScalarHullIndex (const tarch::la::Vector< Dimensions+1, int > &indexCell, const int direction, const int orientation, const int nodesPerAxis)
 Map cell index to an index for the whole hull, i.e.
 
GPUCallableMethod int mapSpaceTimeFaceIndexToScalarCellIndex (const tarch::la::Vector< Dimensions+1, int > &indexFace, const int direction, const int id, const int nodesPerAxis)
 
GPUCallableMethod void gradient_AoS (const double *__restrict__ const QIn, const double *__restrict__ const dudx, const double invDx, const int nodesPerAxis, const int strideQ, const int scalarIndex, double *__restrict__ gradQAux)
 Computes the gradient at a given (space-time) quadrature point.
 
GPUCallableMethod tarch::la::Vector< Dimensions, double > mapToReferenceCoordinates (const tarch::la::Vector< Dimensions, double > &x, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx)
 Map point in physical domain to corresponding coordinates in reference domain.
 
GPUCallableMethod void interpolate_AoS (const double *__restrict__ const UIn, const double *__restrict__ const nodes, const double *__restrict__ const barycentricWeights, const tarch::la::Vector< Dimensions, double > &referenceCoodinates, const int nodesPerAxis, const int strideQ, double *__restrict__ pointwiseQOut)
 Evalutes DG polynomial at arbitrary reference space coordinate in reference domain [0,1]^d.
 
void riemann_setBoundaryState_body_AoS (std::function< void(double *__restrict__ QInside, double *__restrict__ OOutside, const tarch::la::Vector< Dimensions, double > &x, double t, int normal) > boundaryState, double *__restrict__ QOut, double *__restrict__ QIn, const double *const __restrict__ nodes, const tarch::la::Vector< Dimensions, double > &faceCentre, const double dx, const double t, const double dt, const int nodesPerAxis, const int unknowns, const int strideQ, const int direction, const int scalarIndexFace)
 Set boundary state at given coordinates and time.
 
double riemann_maxAbsoluteEigenvalue_body_AoS (std::function< double(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, const int normal) > maxAbsoluteEigenvalue, const double *const __restrict__ QLR, const double *const __restrict__ nodes, const tarch::la::Vector< Dimensions, double > &faceCentre, const double dx, const double t, const double dt, const int nodesPerAxis, const int unknowns, const int strideQ, const int direction, const int scalarIndexFace)
 Determines the maximum absolute value among the eigenvalues computed at a space-time quadrature point on a face.
 
void riemann_setBoundaryState_loop_AoS (std::function< void(const double *const __restrict__ QInside, double *__restrict__ OOutside, const tarch::la::Vector< Dimensions, double > &x, double t, int normal) > boundaryState, double *QHullOut[Dimensions *2], const double *const __restrict__ nodes, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const double dt, const int order, const int unknowns, const int auxiliaryVariables, const bool atBoundary[Dimensions *2])
 Prescribe boundary values (QOut) and overwrite interior values (QInside).
 
void riemann_maxAbsoluteEigenvalue_loop_AoS (std::function< double(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, const int direction) > maxAbsoluteEigenvalue, double maxEigenvaluePerFaceOut[Dimensions *2], double *const QHullIn[Dimensions *2], const double *const __restrict__ nodes, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const double dt, const int order, const int unknowns, const int auxiliaryVariables)
 Per face, determine the eigenvalue with the largest absolute value among the boundary-extrapolated predictor values of the current cell and its face neighbour(s).
 
void rusanovNonlinear_body_AoS (std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux, std::function< void(double *__restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct, double *__restrict__ riemannResultOut, double *__restrict__ FLAux, double *__restrict__ FRAux, double *__restrict__ QAvgAux, double *__restrict__ dQAux, double *__restrict__ SAux, const double *const __restrict__ QLIn, const double *const __restrict__ QRIn, const double smax, const double *const __restrict__ nodes, const double *const __restrict__ weights, const tarch::la::Vector< Dimensions, double > &faceCentre, const double dx, const double t, const double dt, const int nodesPerAxis, const int unknowns, const int strideQ, const int strideF, const int direction, const bool orientationToCell, const bool callFlux, const bool callNonconservativeProduct, const int scalarIndexFace)
 
void rusanovNonlinear_loop_AoS (std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux, std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > boundaryFlux, std::function< void(double *__restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct, std::function< void(double *__restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > boundaryNonconservativeProduct, double *__restrict__ riemannResultOut, const double *const __restrict__ QHullIn[Dimensions *2], const double maxEigenvaluePerFace[Dimensions *2], const double *const __restrict__ nodes, const double *const __restrict__ weights, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const double dt, const int order, const int unknowns, const int auxiliaryVariables, const bool atBoundary[Dimensions *2], const bool callFlux, const bool callNonconservativeProduct)
 Apply a Rusanov Riemann solver to the pair of extrapolated predictor values from two neighbouring cells.
 
STP Picard loop kernel bodies
void spaceTimePredictor_initialGuess_body_AoS (double *__restrict__ Qout, const double *__restrict__ UIn, const int nodesPerAxis, const int strideQ, const int scalarIndex)
 Initial guess for Q before STP Picard loop.
 
void spaceTimePredictor_PicardLoop_initialiseRhs_AoS (double *__restrict__ rhsOut, const double *const __restrict__ UIn, const double *const __restrict__ FLCoeff, const int nodesPerAxis, const int strideQ, const int strideRhs, const int scalarIndex)
 Initialise RHS of Picard iterations.
 
void spaceTimePredictor_PicardLoop_addContributions_body_AoS (std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux, std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, double *__restrict__ S) > algebraicSource, std::function< void(const double *const __restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct, double *__restrict__ rhsOut, double *__restrict__ FAux, double *__restrict__ gradQAux, double *__restrict__ SAux, const double *__restrict__ QIn, const double *__restrict__ nodes, const double *__restrict__ weights, const double *__restrict__ Kxi, const double *__restrict__ dudx, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const double dt, const int nodesPerAxis, const int unknowns, const int strideQ, const int strideRhs, const bool callFlux, const bool callSource, const bool callNonconservativeProduct, const int scalarIndex)
 Add flux, NCP, and source contributions to RHS during STP Picard loop.
 
void spaceTimePredictor_PicardLoop_addFluxContributionsToRhs_body_AoS (std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux, double *__restrict__ rhsOut, double *__restrict__ FAux, const double *__restrict__ QIn, const double *__restrict__ nodes, const double *__restrict__ weights, const double *__restrict__ Kxi, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const double dt, const int nodesPerAxis, const int unknowns, const int strideQ, const int strideRhs, const int scalarIndex)
 Add flux contributions to RHS during STP Picard loop.
 
void spaceTimePredictor_PicardLoop_addSourceContributionToRhs_body_AoS (std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, double *__restrict__ S) > algebraicSource, double *__restrict__ rhsOut, double *__restrict__ SAux, const double *__restrict__ QIn, const double *__restrict__ nodes, const double *__restrict__ weights, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const double dt, const int nodesPerAxis, const int unknowns, const int strideQ, const int strideRhs, const int scalarIndex)
 Add source contributions to RHS during STP Picard loop.
 
void spaceTimePredictor_PicardLoop_addNcpContributionToRhs_body_AoS (std::function< void(const double *const __restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct, double *__restrict__ rhsOut, double *__restrict__ gradQAux, double *__restrict__ SAux, const double *const __restrict__ QIn, const double *const __restrict__ nodes, const double *const __restrict__ weights, const double *const __restrict__ dudx, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const double dt, const int nodesPerAxis, const int unknowns, const int strideQ, const int strideRhs, const int scalarIndex)
 
void spaceTimePredictor_PicardLoop_invert_body_AoS (double *__restrict__ QOut, double &squaredResiduumOut, const double *const __restrict__ rhsIn, const double *const __restrict__ iK1, const int nodesPerAxis, const int unknowns, const int strideQ, const int strideRhs, const int scalarIndex)
 Invert the Picard space-time system to compute the next solution.
 
void spaceTimePredictor_extrapolate_body_AoS (double *__restrict__ QHullOut[Dimensions *2], const double *const __restrict__ QIn, const double *const __restrict__ FLRCoeff[2], const int nodesPerAxis, const int strideQLR, const int strideQ, const int scalarIndexHull)
 Extrapolate the predictor onto DoF associated with the faces of a cell (cell hull).
 
void spaceTimePredictor_extrapolate_Lobatto_body_AoS (double *__restrict__ QHullOut[Dimensions *2], const double *const __restrict__ QIn, const int nodesPerAxis, const int strideQLR, const int strideQ, const int scalarIndexHull)
 Extrapolate the predictor onto DoF associated with the faces of a cell (cell hull).
 
void spaceTimePredictor_extrapolateInTime_body_AoS (double *__restrict__ UOut, const double *const __restrict__ QIn, const double *const __restrict__ FRCoeff, const double nodesPerAxis, const int strideQ, const int scalarIndexCell)
 Extrapolates the predictor to time t+dt.
 
void spaceTimePredictor_extrapolateInTime_Lobatto_body_AoS (double *__restrict__ UOut, const double *const __restrict__ QIn, const double *const __restrict__ FRCoeff, const double nodesPerAxis, const int strideQ, const int scalarIndexCell)
 Extrapolates the predictor to time t+dt.
 
void spaceTimePredictor_PicardLoop_loop_AoS (std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux, std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, double *__restrict__ S) > algebraicSource, std::function< void(const double *const __restrict__ Q, const double *const __restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct, double *__restrict__ QOut, const double *const __restrict__ UIn, const double *const __restrict__ nodes, const double *const __restrict__ weights, const double *const __restrict__ Kxi, const double *const __restrict__ iK1, const double *const __restrict__ FLCoeff, const double *const __restrict__ dudx, const tarch::la::Vector< Dimensions, double > &cellCentre, const double dx, const double t, const double dt, const int order, const int unknowns, const int auxiliaryVariables, const double atol, const bool callFlux, const bool callSource, const bool callNonconservativeProduct)
 Compute the space-time predictor (Qout) from the current solution (UIn).
 
void spaceTimePredictor_extrapolate_loop_AoS (double *__restrict__ QHullOut[Dimensions *2], const double *const __restrict__ QIn, const double *const __restrict__ FLCoeff, const double *const __restrict__ FRCoeff, const int order, const int unknowns, const int auxiliaryVariables)
 Extrapolate all of cell's predictor DOF to the hull of the element.
 
void spaceTimePredictor_extrapolate_Lobatto_loop_AoS (double *__restrict__ QHullOut[Dimensions *2], const double *const __restrict__ QIn, const double *const __restrict__ FLCoeff, const double *const __restrict__ FRCoeff, const int order, const int unknowns, const int auxiliaryVariables)
 Extrapolate all of cell's predictor DOF to the hull of the element.
 
void spaceTimePredictor_extrapolateInTime_loop_AoS (double *__restrict__ UOut, const double *const __restrict__ QIn, const double *const __restrict__ FRCoeff, const int order, const int unknowns, const int auxiliaryVariables)
 Extrapolates the predictor to t+dt.
 
void spaceTimePredictor_extrapolateInTime_Lobatto_loop_AoS (double *__restrict__ UOut, const double *const __restrict__ QIn, const double *const __restrict__ FRCoeff, const int order, const int unknowns, const int auxiliaryVariables)
 Extrapolates the predictor to t+dt.
 

Function Documentation

◆ computeCellErrorIntegral()

template<typename T >
void exahype2::aderdg::computeCellErrorIntegral ( std::function< void(const tarch::la::Vector< Dimensions, double > &position, const double t, const double dt, double *sol) > exactSolution,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const tarch::la::Vector< Dimensions, double > & cellSize,
double t,
double dt,
int Order,
const double *__restrict__ quadratureWeights,
const double *__restrict__ quadraturePoints,
const int unknowns,
const int auxiliary_variables,
T *__restrict__ Q,
double errors[3] )

Definition at line 2 of file AderUtils.cpph.

References dfor, and getQuadraturePoint().

Here is the call graph for this function:

◆ corrector_addCellContributions_loop_AoS()

void exahype2::aderdg::corrector_addCellContributions_loop_AoS ( std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux,
std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, double *__restrict__ S) > algebraicSource,
std::function< void(const double *const __restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct,
double *__restrict__ UOut,
const double *const __restrict__ QIn,
const double *const __restrict__ nodes,
const double *const __restrict__ weights,
const double *const __restrict__ Kxi,
const double *const __restrict__ dudx,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const double dt,
const int order,
const int unknowns,
const int auxiliaryVariables,
const bool callFlux,
const bool callSource,
const bool callNonconservativeProduct )

Add cell-local contributions to new solution or update vector.

Parameters
[in]flux
[in]algebraicSource
[in]nonconservativeProduct
[in,out]UOut
[in]QIn
[in]weightsquadrature weights; size: (order+1)
[in]nodesquadrature nodes; size: (order+1)
[in]Kxistiffness matrix; size: (order+1)*(order+1)
[in]dudxderivative operator; size: (order+1)*(order+1)
[in]cellCentre
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]ttime stamp
[in]dttime step size
[in]orderthe DG approximation order, which corresponds to order+1 DG nodes/Lagrange basis functions per coordinate axis
[in]unknownsthe number of PDE unknowns that we evolve
[in]auxiliaryVariablesother quantities such as material parameters that we do not evolve
[in]callFlux
[in]callSource
[in]callNonconservativeProduct

Definition at line 266 of file CorrectorAoS.cpp.

References corrector_addFluxContributions_body_AoS(), corrector_addNcpContributions_body_AoS(), corrector_addSourceContributions_body_AoS(), dx, and getNodesPerCell().

Referenced by exahype2::aderdg::tests::ADERDGTest::runADERDGStep().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ corrector_addFluxContributions_body_AoS()

void exahype2::aderdg::corrector_addFluxContributions_body_AoS ( std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux,
double *__restrict__ UOut,
const double *__restrict__ QIn,
double *__restrict__ FAux,
const double *__restrict__ nodes,
const double *__restrict__ weights,
const double *__restrict__ Kxi,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const double dt,
const int nodesPerAxis,
const int unknowns,
const int strideQ,
const int scalarIndex )

Add space-time volume flux contributions to the solution.

Writes to a single spatial degree of freedom.

Note
Directly run after Picard iterations in order to not store predictor.
Parameters
[in]flux
[in,out]UOut
[in]QIn
[in,out]FAux
[in]nodesquadrature nodes; size: (order+1)
[in]weightsquadrature weights; size: (order+1)
[in]Kxistiffness matrix; size: (order+1)*(order+1)
[in]cellCentre
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]ttime stamp
[in]dttime step size
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]unknownsthe number of PDE unknowns that we evolve
[in]strideQ
[in]scalarIndex

Definition at line 64 of file CorrectorAoS.cpp.

References delineariseIndex(), dx, getCoordinates(), getStrides(), and it.

Referenced by corrector_addCellContributions_loop_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ corrector_addNcpContributions_body_AoS()

void exahype2::aderdg::corrector_addNcpContributions_body_AoS ( std::function< void(const double *const __restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct,
double *__restrict__ UOut,
double *__restrict__ gradQAux,
double *__restrict__ SAux,
const double *__restrict__ QIn,
const double *__restrict__ nodes,
const double *__restrict__ weights,
const double *__restrict__ dudx,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const double dt,
const int nodesPerAxis,
const int unknowns,
const int strideQ,
const int scalarIndex )

Add nonconservative product contributions to the solution.

Writes to a single spatial degree of freedom.

Note
Directly run after Picard iterations in order to not store predictor.
Parameters
[in]nonconservativeProduct
[in,out]UOut
[in]QIn
[in,out]gradQAux
[in,out]SAux
[in]nodesquadrature nodes; size: (order+1)
[in]weightsquadrature weights; size: (order+1)
[in]dudxderivative operator; size: (order+1)*(order+1)
[in]cellCentre
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]ttime stamp
[in]dttime step size
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]unknownsthe number of PDE unknowns that we evolve
[in]strideQ
[in]scalarIndex

Definition at line 166 of file CorrectorAoS.cpp.

References delineariseIndex(), dx, getCoordinates(), getStrides(), gradient_AoS(), and it.

Referenced by corrector_addCellContributions_loop_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ corrector_addRiemannContributions_body_AoS()

GPUCallableMethod void exahype2::aderdg::corrector_addRiemannContributions_body_AoS ( double *__restrict__ UOut,
const double *const __restrict__ riemannResultIn,
const double *const __restrict__ weights,
const double *const __restrict__ FLCoeff,
const double dx,
const double dt,
const int nodesPerAxis,
const int unknowns,
const int strideQ,
const int strideF,
const int scalarIndexCell )

Add Riemann flux contributions to the solution.

Writes to a single spatial degree of freedom.

Parameters
[in,out]UOut
[in]riemannResultIn
[in]weightsquadrature weights; size: (order+1)
[in]FCoeff
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]unknownsthe number of PDE unknowns that we evolve
[in]strideQ
[in]scalarIndex

Definition at line 223 of file CorrectorAoS.cpp.

References delineariseIndex(), dx, FLCoeff, getStrides(), and mapCellIndexToScalarHullIndex().

Referenced by corrector_addRiemannContributions_loop_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ corrector_addRiemannContributions_loop_AoS()

double exahype2::aderdg::corrector_addRiemannContributions_loop_AoS ( std::function< void(double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t) > adjustSolution,
std::function< double(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, const int normal) > maxAbsoluteEigenvalue,
double *__restrict__ UOut,
const double *const __restrict__ riemannResultIn,
const double *const __restrict__ nodes,
const double *const __restrict__ weights,
const double *const __restrict__ FLCoeff,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const double dt,
const int order,
const int unknowns,
const int auxiliaryVariables,
const bool callMaxEigenvalue )

Add cell-local contributions to new solution or update vector.

Note
This routine must be called AFTER corrector_addCellContributions_loop_AoS ! Otherwise calling adjustSolution and maxAbsoluteEigenvalue will not have the intended effect.
Returns
maximum eigenvalue (absolute value) computed at the quadrature nodes in the interior of the cell after the full DG update and after all solution adjustments.
Parameters
[in,out]UOut
[in]riemannResultIn
[in]weightsquadrature weights; size: (order+1)
[in]FLCoeffvalues of basis functions evaluated at x=0.0 (left); size: order+1
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]dttime step size
[in]orderthe DG approximation order, which corresponds to order+1 DG nodes/Lagrange basis functions per coordinate axis
[in]unknownsthe number of PDE unknowns that we evolve
[in]auxiliaryVariablesother quantities such as material parameters that we do not evolve

Definition at line 394 of file CorrectorAoS.cpp.

References corrector_addRiemannContributions_body_AoS(), corrector_adjustSolution_computeMaxEigenvalue_body_AoS(), dx, FLCoeff, and getNodesPerCell().

Referenced by exahype2::aderdg::tests::ADERDGTest::runADERDGStep().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ corrector_addSourceContributions_body_AoS()

void exahype2::aderdg::corrector_addSourceContributions_body_AoS ( std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, double *__restrict__ S) > algebraicSource,
double *__restrict__ UOut,
double *__restrict__ SAux,
const double *__restrict__ QIn,
const double *__restrict__ nodes,
const double *__restrict__ weights,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const double dt,
const int nodesPerAxis,
const int unknowns,
const int strideQ,
const int scalarIndex )

Add source contributions to the solution.

Writes to a single spatial degree of freedom.

Note
Directly run after Picard iterations in order to not store predictor.
Parameters
[in]algebraicSource
[in,out]UOut
[in,out]SAux
[in]QIn
[in]nodesquadrature nodes; size: (order+1)
[in]weightsquadrature weights; size: (order+1)
[in]cellCentre
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]ttime stamp
[in]dttime step size
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]unknownsthe number of PDE unknowns that we evolve
[in]strideQ
[in]scalarIndex

Definition at line 121 of file CorrectorAoS.cpp.

References delineariseIndex(), dx, getCoordinates(), getStrides(), and it.

Referenced by corrector_addCellContributions_loop_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ corrector_adjustSolution_computeMaxEigenvalue_body_AoS()

double exahype2::aderdg::corrector_adjustSolution_computeMaxEigenvalue_body_AoS ( std::function< void(double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t) > adjustSolution,
std::function< double(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, const int normal) > maxAbsoluteEigenvalue,
double *__restrict__ UOut,
const double *const __restrict__ nodes,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const int nodesPerAxis,
const int unknowns,
const int strideQ,
const int callMaxEigenvalue,
const int scalarIndex )

Allow user to modify solution at the given coordinates and time and compute the max eigenvalue afterwards.

Writes to a single spatial degree of freedom.

Parameters
[in]adjustSolution
[in,out]UOut
[in]nodesquadrature nodes; size: (order+1)
[in]cellCentre
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]ttime stamp
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]unknownsthe number of PDE unknowns that we evolve
[in]strideQ
[in]scalarIndex

Definition at line 16 of file CorrectorAoS.cpp.

References delineariseIndex(), dx, getCoordinates(), and getStrides().

Referenced by corrector_addRiemannContributions_loop_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ delineariseIndex()

◆ getCoordinates()

GPUCallableMethod tarch::la::Vector< Dimensions+1, double > exahype2::aderdg::getCoordinates ( const tarch::la::Vector< Dimensions+1, int > & index,
const tarch::la::Vector< Dimensions, double > & centre,
const double dx,
const double t,
const double dt,
const double *__restrict__ nodes )

Compute coordinates from cell geometry and quadrature index.

Returns
a tuple (t,x,y,z) where the t-component has index 0.
Parameters
[in]indext-,x-,y-, and z-direction (reference coordinates) components of descalar scalar index time stamp
Note
coords[0] = t if if t-direction component of index is negative.

Definition at line 116 of file KernelUtils.cpp.

References dx.

Referenced by corrector_addFluxContributions_body_AoS(), corrector_addNcpContributions_body_AoS(), corrector_addSourceContributions_body_AoS(), corrector_adjustSolution_computeMaxEigenvalue_body_AoS(), spaceTimePredictor_PicardLoop_addContributions_body_AoS(), spaceTimePredictor_PicardLoop_addFluxContributionsToRhs_body_AoS(), spaceTimePredictor_PicardLoop_addNcpContributionToRhs_body_AoS(), and spaceTimePredictor_PicardLoop_addSourceContributionToRhs_body_AoS().

Here is the caller graph for this function:

◆ getCoordinatesOnFace()

GPUCallableMethod tarch::la::Vector< Dimensions+1, double > exahype2::aderdg::getCoordinatesOnFace ( const tarch::la::Vector< Dimensions+1, int > & indexOnFace,
const tarch::la::Vector< Dimensions, double > & faceCentre,
const int direction,
const double dx,
const double t,
const double dt,
const double *__restrict__ nodes )

Compute coordinates from cell geometry and quadrature index.

Returns
a tuple (t,x,y,z) where the t-component has index 0.
Parameters
[in]indext-,x-,y-, and z-direction (reference coordinates) components of descalar scalar index time stamp
[in]directionencodes direction of face normal (x: 0, y: 1, z: 2)
Note
coords[0] = t if if t-direction component of index is negative.

Definition at line 145 of file KernelUtils.cpp.

References dx.

Referenced by riemann_maxAbsoluteEigenvalue_body_AoS(), riemann_setBoundaryState_body_AoS(), and rusanovNonlinear_body_AoS().

Here is the caller graph for this function:

◆ getNodesPerCell()

◆ getQuadraturePoint() [1/2]

tarch::la::Vector< 2, double > exahype2::aderdg::getQuadraturePoint ( const tarch::la::Vector< 2, double > & cellCentre,
const tarch::la::Vector< 2, double > & cellSize,
const tarch::la::Vector< 2, int > & index,
int polynomialOrder,
const double *__restrict__ quadraturePoints )

In ExaHyPE's Finite Volume setup, a cell hosts a patch of Finite Volumes.

When we iterate over these volumes, we typically have to know the centre of the volume.

Parameters
quadraturePoints1d quadrature points over unit interval. This array has polynomialOrder+1 entries.

Definition at line 6 of file AderUtils.cpp.

References tarch::la::multiplyComponents().

Referenced by computeCellErrorIntegral().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ getQuadraturePoint() [2/2]

tarch::la::Vector< 3, double > exahype2::aderdg::getQuadraturePoint ( const tarch::la::Vector< 3, double > & cellCentre,
const tarch::la::Vector< 3, double > & cellSize,
const tarch::la::Vector< 3, int > & index,
int polynomialOrder,
const double *__restrict__ quadraturePoints )

Definition at line 22 of file AderUtils.cpp.

References tarch::la::multiplyComponents().

Here is the call graph for this function:

◆ getSpaceTimeNodesPerCell()

GPUCallableMethod int exahype2::aderdg::getSpaceTimeNodesPerCell ( int nodesPerAxis)

Definition at line 25 of file KernelUtils.cpp.

Referenced by spaceTimePredictor_PicardLoop_loop_AoS().

Here is the caller graph for this function:

◆ getStrides()

GPUCallableMethod tarch::la::Vector< Dimensions+1, int > exahype2::aderdg::getStrides ( int nodesPerAxis,
bool strides4SpaceTimeQuantity = true )

If the stride is used for a space-time quantity, then this function generates:

(1,nodesPerAxis,nodesPerAxis^2,nodesPerAxis^3) for strides (t,x,y,z).

Otherwise, for time-invariant fields, it generates

(0,1,nodesPerAxis,nodesPerAxis^2) for strides (t,x,y,z).

Parameters
[in]strides4SpaceTimeQuantitygenerate strides for space-time quantity (default=true).
[in]nodesPerAxisnodes per coordinate axis nodes/Lagrange basis functions per coordinate axis (order+1)
Returns
Strides per direction (t,x,y,z), where t-direction stride has index 0.

Definition at line 41 of file KernelUtils.cpp.

Referenced by corrector_addFluxContributions_body_AoS(), corrector_addNcpContributions_body_AoS(), corrector_addRiemannContributions_body_AoS(), corrector_addSourceContributions_body_AoS(), corrector_adjustSolution_computeMaxEigenvalue_body_AoS(), gradient_AoS(), interpolate_AoS(), mapSpaceTimeFaceIndexToScalarCellIndex(), riemann_maxAbsoluteEigenvalue_body_AoS(), riemann_setBoundaryState_body_AoS(), rusanovNonlinear_body_AoS(), spaceTimePredictor_extrapolate_body_AoS(), spaceTimePredictor_extrapolate_Lobatto_body_AoS(), spaceTimePredictor_PicardLoop_addContributions_body_AoS(), spaceTimePredictor_PicardLoop_addFluxContributionsToRhs_body_AoS(), spaceTimePredictor_PicardLoop_addNcpContributionToRhs_body_AoS(), spaceTimePredictor_PicardLoop_addSourceContributionToRhs_body_AoS(), spaceTimePredictor_PicardLoop_initialiseRhs_AoS(), and spaceTimePredictor_PicardLoop_invert_body_AoS().

Here is the caller graph for this function:

◆ gradient_AoS()

GPUCallableMethod void exahype2::aderdg::gradient_AoS ( const double *__restrict__ const QIn,
const double *__restrict__ const dudx,
const double invDx,
const int nodesPerAxis,
const int strideQ,
const int scalarIndex,
double *__restrict__ gradQAux )

Computes the gradient at a given (space-time) quadrature point.

Computes

u(x,y,z) = sum_ijk u_j phi_i(x) phi_j(y) phi_k(z)

=> (d/dx u) (x,y,z) = sum_ijk u_ijk (d/dx phi_i)(x) phi_j(y) phi_k(z) (d/dy u) (x,y,z) = sum_ijk u_ijk phi_i(x) (d/dy phi_j)(y) phi_k(z) (d/dz u) (x,y,z) = sum_ijk u_ijk phi_i(x) phi_k(y) (d/dz phi_k)(z)

Lagrange basis property simplifies this to:

=> (d/dx u) (x_a,y_i,z_k) = sum_a u_ijk (d/dx phi_i)(x_a) (d/dy u) (x_i,y_a,z_k) = sum_a u_ijk (d/dy phi_j)(y_a) (d/dz u) (x_i,y_j,z_a) = sum_a u_ijk (d/dz phi_k)(z_a)

For elements that do not equal the unit cube, we write:

=> (d/dx u) (x_a,y_i,z_k) = sum_a u_ajk 1/lx * (d/dx phi_i)(x_a) (d/dy u) (x_i,y_a,z_k) = sum_a u_iak 1/ly * (d/dy phi_j)(y_a) (d/dz u) (x_i,y_j,z_a) = sum_a u_ija 1/lz * (d/dz phi_k)(z_a)

where lx,ly,lz are the lengths of the element.

Note
Assumes degrees of freedom are stored in order (iz,iy,ix,it)
This function assumes that the gradient zeroed out
Parameters
[in]QIn
[in]dudxderivative operator; size: (order+1)*(order+1)
[in]invDx
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]unknownsthe number of PDE unknowns that we evolve
[in]strideQ
[in]strideGradQ
[in]scalarIndex
[in]gradQAux

Definition at line 267 of file KernelUtils.cpp.

References delineariseIndex(), and getStrides().

Referenced by corrector_addNcpContributions_body_AoS(), spaceTimePredictor_PicardLoop_addContributions_body_AoS(), and spaceTimePredictor_PicardLoop_addNcpContributionToRhs_body_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ interpolate_AoS()

GPUCallableMethod void exahype2::aderdg::interpolate_AoS ( const double *__restrict__ const UIn,
const double *__restrict__ const nodes,
const double *__restrict__ const barycentricWeights,
const tarch::la::Vector< Dimensions, double > & referenceCoodinates,
const int nodesPerAxis,
const int strideQ,
double *__restrict__ pointwiseQOut )

Evalutes DG polynomial at arbitrary reference space coordinate in reference domain [0,1]^d.

See also
: mapToReferenceCoordinates

Definition at line 329 of file KernelUtils.cpp.

References delineariseIndex(), getNodesPerCell(), and getStrides().

Referenced by exahype2::aderdg::tests::ADERDGTest::testInterpolate().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ lineariseIndex()

GPUCallableMethod int exahype2::aderdg::lineariseIndex ( const tarch::la::Vector< Dimensions+1, int > & index,
const tarch::la::Vector< Dimensions+1, int > & strides )
Parameters
[in]linearisesthe index with the given strides.
[in]stridesstrides per direction (t,x,y,z). time stamp
Note
if stride[0] = 0 this implies that we have no time-index

Definition at line 66 of file KernelUtils.cpp.

Referenced by mapSpaceTimeFaceIndexToScalarCellIndex().

Here is the caller graph for this function:

◆ mapCellIndexToScalarFaceIndex()

GPUCallableMethod int exahype2::aderdg::mapCellIndexToScalarFaceIndex ( const tarch::la::Vector< Dimensions+1, int > & indexCell,
const int direction,
const int nodesPerAxis )

Map cell index to an index for a single face.

Parameters
[in]directioncoordinate direction of the (reference) element face normal (0:
Note
Result must be scaled additionally by nodesPerAxis if it used to access space-time quantities.

Definition at line 181 of file KernelUtils.cpp.

◆ mapCellIndexToScalarHullIndex()

GPUCallableMethod int exahype2::aderdg::mapCellIndexToScalarHullIndex ( const tarch::la::Vector< Dimensions+1, int > & indexCell,
const int direction,
const int orientation,
const int nodesPerAxis )

Map cell index to an index for the whole hull, i.e.

all 2*Dimensios faces.

Parameters
[in]directioncoordinate direction of the (reference) element face normal (0:
[in]orientationorientation of the (reference) element face normal (0: negative, 1: positive).
Note
Result must be scaled additionally by nodesPerAxis if it used to access space-time quantities.

Definition at line 208 of file KernelUtils.cpp.

Referenced by corrector_addRiemannContributions_body_AoS().

Here is the caller graph for this function:

◆ mapSpaceTimeFaceIndexToScalarCellIndex()

GPUCallableMethod int exahype2::aderdg::mapSpaceTimeFaceIndexToScalarCellIndex ( const tarch::la::Vector< Dimensions+1, int > & indexFace,
const int direction,
const int id,
const int nodesPerAxis )
Parameters
[in]directioncoordinate direction of the (reference) element face normal (0:
[in]orientationorientation of the (reference) element face normal (0: negative, 1: positive).
Returns
scalar cell index

Definition at line 239 of file KernelUtils.cpp.

References getStrides(), and lineariseIndex().

Referenced by spaceTimePredictor_extrapolate_body_AoS(), and spaceTimePredictor_extrapolate_Lobatto_body_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ mapToReferenceCoordinates()

GPUCallableMethod tarch::la::Vector< Dimensions, double > exahype2::aderdg::mapToReferenceCoordinates ( const tarch::la::Vector< Dimensions, double > & x,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx )

Map point in physical domain to corresponding coordinates in reference domain.

Parameters
[in]xcoordinate
[in]cellCentrecentre of the cell that contains x
[in]dxextent of the cell that contains x
Returns

Definition at line 306 of file KernelUtils.cpp.

References dx.

◆ plotNode()

std::string exahype2::aderdg::plotNode ( double Q[],
int unknowns )

Render a DG nodes components as a string.

Definition at line 11 of file Generic.cpp.

◆ riemann_maxAbsoluteEigenvalue_body_AoS()

double exahype2::aderdg::riemann_maxAbsoluteEigenvalue_body_AoS ( std::function< double(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, const int normal) > maxAbsoluteEigenvalue,
const double *const __restrict__ QLR,
const double *const __restrict__ nodes,
const tarch::la::Vector< Dimensions, double > & faceCentre,
const double dx,
const double t,
const double dt,
const int nodesPerAxis,
const int unknowns,
const int strideQ,
const int direction,
const int scalarIndexFace )

Determines the maximum absolute value among the eigenvalues computed at a space-time quadrature point on a face.

Note
: Different to the ExaHyPE1 implementation, we compute the maximum of the absolute values of the eigenvalues using all space-time quadrature points.
Returns
the largest absolute eigenvalue
Parameters
[in]maxAbsoluteEigenvalue
[in]QLR

param [in] t time stamp

Parameters
[in]dttime step size
[in]faceCentre
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]unknownsthe number of PDE unknowns that we evolve
[in]strideQ
[in]direction
[in]scalarIndexFace

Definition at line 50 of file RiemannAoS.cpp.

References delineariseIndex(), dx, getCoordinatesOnFace(), and getStrides().

Referenced by riemann_maxAbsoluteEigenvalue_loop_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ riemann_maxAbsoluteEigenvalue_loop_AoS()

void exahype2::aderdg::riemann_maxAbsoluteEigenvalue_loop_AoS ( std::function< double(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, const int direction) > maxAbsoluteEigenvalue,
double maxEigenvaluePerFaceOut[Dimensions *2],
double *const QHullIn[Dimensions *2],
const double *const __restrict__ nodes,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const double dt,
const int order,
const int unknowns,
const int auxiliaryVariables )

Per face, determine the eigenvalue with the largest absolute value among the boundary-extrapolated predictor values of the current cell and its face neighbour(s).

Returns
maxAbsoluteEigenvalue
Parameters
[in,out]maxEigenvaluePerFaceOut
[in]QHullIn
[in]nodesquadrature nodes; size: (order+1)
[in]ttime stamp
[in]dttime step size
[in]faceCentre
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]orderthe DG approximation order, which corresponds to order+1 DG nodes/Lagrange basis functions per coordinate axis
[in]unknownsthe number of PDE unknowns that we evolve
[in]auxiliaryVariablesother quantities such as material parameters that we do not evolve
[in]direction

Definition at line 146 of file RiemannAoS.cpp.

References dx, getNodesPerCell(), and riemann_maxAbsoluteEigenvalue_body_AoS().

Referenced by exahype2::aderdg::tests::ADERDGTest::runADERDGStep().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ riemann_setBoundaryState_body_AoS()

void exahype2::aderdg::riemann_setBoundaryState_body_AoS ( std::function< void(double *__restrict__ QInside, double *__restrict__ OOutside, const tarch::la::Vector< Dimensions, double > &x, double t, int normal) > boundaryState,
double *__restrict__ QOut,
double *__restrict__ QIn,
const double *const __restrict__ nodes,
const tarch::la::Vector< Dimensions, double > & faceCentre,
const double dx,
const double t,
const double dt,
const int nodesPerAxis,
const int unknowns,
const int strideQ,
const int direction,
const int scalarIndexFace )

Set boundary state at given coordinates and time.

Parameters
[in]boundaryState
[in,out]QOut
[in,out]QIn
[in]nodesquadrature nodes; size: (order+1)
[in]ttime stamp
[in]dttime step size
[in]faceCentre
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]unknownsthe number of PDE unknowns that we evolve
[in]strideQ
[in]direction
[in]orientation
[in]scalarIndexFace

Definition at line 12 of file RiemannAoS.cpp.

References delineariseIndex(), dx, getCoordinatesOnFace(), and getStrides().

Referenced by riemann_setBoundaryState_loop_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ riemann_setBoundaryState_loop_AoS()

void exahype2::aderdg::riemann_setBoundaryState_loop_AoS ( std::function< void(const double *const __restrict__ QInside, double *__restrict__ OOutside, const tarch::la::Vector< Dimensions, double > &x, double t, int normal) > boundaryState,
double * QHullOut[Dimensions *2],
const double *const __restrict__ nodes,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const double dt,
const int order,
const int unknowns,
const int auxiliaryVariables,
const bool atBoundary[Dimensions *2] )

Prescribe boundary values (QOut) and overwrite interior values (QInside).

Parameters
[in]boundaryState
[in,out]QHullOut
[in]nodesquadrature nodes; size: (order+1)
[in]ttime stamp
[in]dttime step size
[in]faceCentre
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]orderthe DG approximation order, which corresponds to order+1 DG nodes/Lagrange basis functions per coordinate axis
[in]unknownsthe number of PDE unknowns that we evolve
[in]auxiliaryVariablesother quantities such as material parameters that we do not evolve
[in]atBoundary

Definition at line 87 of file RiemannAoS.cpp.

References dx, getNodesPerCell(), and riemann_setBoundaryState_body_AoS().

Here is the call graph for this function:

◆ rusanovNonlinear_body_AoS()

void exahype2::aderdg::rusanovNonlinear_body_AoS ( std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux,
std::function< void(double *__restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct,
double *__restrict__ riemannResultOut,
double *__restrict__ FLAux,
double *__restrict__ FRAux,
double *__restrict__ QAvgAux,
double *__restrict__ dQAux,
double *__restrict__ SAux,
const double *const __restrict__ QLIn,
const double *const __restrict__ QRIn,
const double smax,
const double *const __restrict__ nodes,
const double *const __restrict__ weights,
const tarch::la::Vector< Dimensions, double > & faceCentre,
const double dx,
const double t,
const double dt,
const int nodesPerAxis,
const int unknowns,
const int strideQ,
const int strideF,
const int direction,
const bool orientationToCell,
const bool callFlux,
const bool callNonconservativeProduct,
const int scalarIndexFace )

Definition at line 15 of file RusanovNonlinearAoS.cpp.

References delineariseIndex(), dx, getCoordinatesOnFace(), getStrides(), and it.

Referenced by rusanovNonlinear_loop_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ rusanovNonlinear_loop_AoS()

void exahype2::aderdg::rusanovNonlinear_loop_AoS ( std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux,
std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > boundaryFlux,
std::function< void(double *__restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct,
std::function< void(double *__restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > boundaryNonconservativeProduct,
double *__restrict__ riemannResultOut,
const double *const __restrict__ QHullIn[Dimensions *2],
const double maxEigenvaluePerFace[Dimensions *2],
const double *const __restrict__ nodes,
const double *const __restrict__ weights,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const double dt,
const int order,
const int unknowns,
const int auxiliaryVariables,
const bool atBoundary[Dimensions *2],
const bool callFlux,
const bool callNonconservativeProduct )

Apply a Rusanov Riemann solver to the pair of extrapolated predictor values from two neighbouring cells.


Parameters
[in]flux

param [in] boundaryFlux,

Parameters
[in]nonconservativeProduct

param [in] boundaryNonconservativeProduct,

Parameters
[in,out]FLOut
[in,out]FROut
[in]QHullIn

param [in] maxEigenvaluePerFace,

Parameters
[in]smax
[in]nodesquadrature nodes; size: (order+1)
[in]weightsquadrature weights; size: (order+1)
[in]ttime stamp
[in]dttime step size
[in]faceCentre
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]orderthe DG approximation order, which corresponds to order+1 DG nodes/Lagrange basis functions per coordinate axis
[in]unknownsthe number of PDE unknowns that we evolve
[in]auxiliaryVariablesother quantities such as material parameters that we do not evolve
[in]direction
[in]callFlux
[in]callNonconservativeProduct

Definition at line 122 of file RusanovNonlinearAoS.cpp.

References dx, getNodesPerCell(), and rusanovNonlinear_body_AoS().

Referenced by exahype2::aderdg::tests::ADERDGTest::runADERDGStep().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ spaceTimePredictor_extrapolate_body_AoS()

void exahype2::aderdg::spaceTimePredictor_extrapolate_body_AoS ( double *__restrict__ QHullOut[Dimensions *2],
const double *const __restrict__ QIn,
const double *const __restrict__ FLRCoeff[2],
const int nodesPerAxis,
const int strideQLR,
const int strideQ,
const int scalarIndexHull )

Extrapolate the predictor onto DoF associated with the faces of a cell (cell hull).

These serve as Riemann solver inputs.

Parameters
[in,out]QHullOut
[in]QIn
[in]FLRCoeffBasis functions evaluated at reference coordinates 0.0 (L,component 0) and 1.0 (R, component 1).
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]strideQLR
[in]strideQ
[in]scalarIndexHull

Definition at line 377 of file PredictorAoS.cpp.

References delineariseIndex(), getStrides(), and mapSpaceTimeFaceIndexToScalarCellIndex().

Referenced by spaceTimePredictor_extrapolate_loop_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ spaceTimePredictor_extrapolate_Lobatto_body_AoS()

void exahype2::aderdg::spaceTimePredictor_extrapolate_Lobatto_body_AoS ( double *__restrict__ QHullOut[Dimensions *2],
const double *const __restrict__ QIn,
const int nodesPerAxis,
const int strideQLR,
const int strideQ,
const int scalarIndexHull )

Extrapolate the predictor onto DoF associated with the faces of a cell (cell hull).

These serve as Riemann solver inputs.

Simplified version of spaceTimePredictor_extrapolate_body_AoS specifically for Lobatto nodes.

Parameters
[in,out]QHullOut
[in]QIn
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]strideQLR
[in]strideQ
[in]scalarIndexHull

Definition at line 411 of file PredictorAoS.cpp.

References delineariseIndex(), getStrides(), and mapSpaceTimeFaceIndexToScalarCellIndex().

Referenced by spaceTimePredictor_extrapolate_Lobatto_loop_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ spaceTimePredictor_extrapolate_Lobatto_loop_AoS()

void exahype2::aderdg::spaceTimePredictor_extrapolate_Lobatto_loop_AoS ( double *__restrict__ QHullOut[Dimensions *2],
const double *const __restrict__ QIn,
const double *const __restrict__ FLCoeff,
const double *const __restrict__ FRCoeff,
const int order,
const int unknowns,
const int auxiliaryVariables )

Extrapolate all of cell's predictor DOF to the hull of the element.

Note
assumes that Gauss-Lobatto nodes are used as support points for the basis functions. In this case, the outermost nodes are located directly at the boundary of the reference element. Boundary-extrapolation thus becomes a simple copy of some predictor coefficients.
FLCoeff,FRCoeff are arguments in order to have the same signature as the generic routine. They are not actually needed.
See also
spaceTimePredictor_extrapolate_loop_AoS, spaceTimePredictor_extrapolate_body_AoS, spaceTimePredictor_extrapolate_Lobatto_body_AoS
Parameters
[in,out]QHullOut
[in]QIn
[in]FLCoeffvalues of basis functions evaluated at x=0.0 (left); size: order+1
[in]FRCoeffvalues of basis functions evaluated at x=1.0 (right); size: order+1
[in]orderthe DG approximation order, which corresponds to order+1 DG nodes/Lagrange basis functions per coordinate axis
[in]unknownsthe number of PDE unknowns that we evolve
[in]auxiliaryVariablesother quantities such as material parameters that we do not evolve

Definition at line 901 of file PredictorAoS.cpp.

References getNodesPerCell(), and spaceTimePredictor_extrapolate_Lobatto_body_AoS().

Here is the call graph for this function:

◆ spaceTimePredictor_extrapolate_loop_AoS()

void exahype2::aderdg::spaceTimePredictor_extrapolate_loop_AoS ( double *__restrict__ QHullOut[Dimensions *2],
const double *const __restrict__ QIn,
const double *const __restrict__ FLCoeff,
const double *const __restrict__ FRCoeff,
const int order,
const int unknowns,
const int auxiliaryVariables )

Extrapolate all of cell's predictor DOF to the hull of the element.

Note
This version works for Gauss-Legendre and Gauss-Lobatto nodes.
Parameters
[in,out]QHullOut
[in]QIn
[in]FLCoeffvalues of basis functions evaluated at x=0.0 (left); size: order+1
[in]FRCoeffvalues of basis functions evaluated at x=1.0 (right); size: order+1
[in]orderthe DG approximation order, which corresponds to order+1 DG nodes/Lagrange basis functions per coordinate axis
[in]unknownsthe number of PDE unknowns that we evolve
[in]auxiliaryVariablesother quantities such as material parameters that we do not evolve

Definition at line 868 of file PredictorAoS.cpp.

References FLCoeff, FRCoeff, getNodesPerCell(), and spaceTimePredictor_extrapolate_body_AoS().

Referenced by exahype2::aderdg::tests::ADERDGTest::runADERDGStep().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ spaceTimePredictor_extrapolateInTime_body_AoS()

void exahype2::aderdg::spaceTimePredictor_extrapolateInTime_body_AoS ( double *__restrict__ UOut,
const double *const __restrict__ QIn,
const double *const __restrict__ FRCoeff,
const double nodesPerAxis,
const int strideQ,
const int scalarIndexCell )

Extrapolates the predictor to time t+dt.

Parameters
[in,out]UOut
[in]QIn
[in]FRCoeffvalues of basis functions evaluated at x=1.0 (right); size: order+1
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]unknownsthe number of PDE unknowns that we evolve
[in]strideQ
[in]scalarIndexHull

Definition at line 444 of file PredictorAoS.cpp.

References FRCoeff, and it.

Referenced by spaceTimePredictor_extrapolateInTime_loop_AoS().

Here is the caller graph for this function:

◆ spaceTimePredictor_extrapolateInTime_Lobatto_body_AoS()

void exahype2::aderdg::spaceTimePredictor_extrapolateInTime_Lobatto_body_AoS ( double *__restrict__ UOut,
const double *const __restrict__ QIn,
const double *const __restrict__ FRCoeff,
const double nodesPerAxis,
const int strideQ,
const int scalarIndexCell )

Extrapolates the predictor to time t+dt.

Simplified version of spaceTimePredictor_extrapolate_body_AoS specifically for Lobatto nodes.

Parameters
[in,out]UOut
[in]QIn
[in]FRCoeffvalues of basis functions evaluated at x=1.0 (right); size: order+1
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]unknownsthe number of PDE unknowns that we evolve
[in]strideQ
[in]scalarIndexHull

Definition at line 473 of file PredictorAoS.cpp.

References it.

Referenced by spaceTimePredictor_extrapolateInTime_Lobatto_loop_AoS().

Here is the caller graph for this function:

◆ spaceTimePredictor_extrapolateInTime_Lobatto_loop_AoS()

void exahype2::aderdg::spaceTimePredictor_extrapolateInTime_Lobatto_loop_AoS ( double *__restrict__ UOut,
const double *const __restrict__ QIn,
const double *const __restrict__ FRCoeff,
const int order,
const int unknowns,
const int auxiliaryVariables )

Extrapolates the predictor to t+dt.

Simplified variant for Lobatto nodes.

Note
FRCoeff is passed as argument in order to have the same signature as the generic routine. This argument is not actually needed.
Parameters
[in,out]UOut
[in]QIn
[in]FRCoeffvalues of basis functions evaluated at x=1.0 (right); size: order+1
[in]dttime step size
[in]orderthe DG approximation order, which corresponds to order+1 DG nodes/Lagrange basis functions per coordinate axis
[in]unknownsthe number of PDE unknowns that we evolve
[in]auxiliaryVariablesother quantities such as material parameters that we do not evolve

Definition at line 955 of file PredictorAoS.cpp.

References FRCoeff, getNodesPerCell(), and spaceTimePredictor_extrapolateInTime_Lobatto_body_AoS().

Here is the call graph for this function:

◆ spaceTimePredictor_extrapolateInTime_loop_AoS()

void exahype2::aderdg::spaceTimePredictor_extrapolateInTime_loop_AoS ( double *__restrict__ UOut,
const double *const __restrict__ QIn,
const double *const __restrict__ FRCoeff,
const int order,
const int unknowns,
const int auxiliaryVariables )

Extrapolates the predictor to t+dt.

The output variable QIn will be completely overwritten. QIn is the space-time polynomial.

Parameters
[in,out]UOut
[in]QIn
[in]FRCoeffvalues of basis functions evaluated at x=1.0 (right); size: order+1
[in]dttime step size
[in]orderthe DG approximation order, which corresponds to order+1 DG nodes/Lagrange basis functions per coordinate axis
[in]unknownsthe number of PDE unknowns that we evolve
[in]auxiliaryVariablesother quantities such as material parameters that we do not evolve

Definition at line 931 of file PredictorAoS.cpp.

References FRCoeff, getNodesPerCell(), and spaceTimePredictor_extrapolateInTime_body_AoS().

Here is the call graph for this function:

◆ spaceTimePredictor_initialGuess_body_AoS()

void exahype2::aderdg::spaceTimePredictor_initialGuess_body_AoS ( double *__restrict__ Qout,
const double *__restrict__ UIn,
const int nodesPerAxis,
const int strideQ,
const int scalarIndex )

Initial guess for Q before STP Picard loop.

This body writes to a single space-time volume coefficient.

Parameters
[in]Qout
[in]UIn
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]unknownsthe number of PDE unknowns that we evolve

Definition at line 16 of file PredictorAoS.cpp.

Referenced by spaceTimePredictor_PicardLoop_loop_AoS().

Here is the caller graph for this function:

◆ spaceTimePredictor_PicardLoop_addContributions_body_AoS()

void exahype2::aderdg::spaceTimePredictor_PicardLoop_addContributions_body_AoS ( std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux,
std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, double *__restrict__ S) > algebraicSource,
std::function< void(const double *const __restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct,
double *__restrict__ rhsOut,
double *__restrict__ FAux,
double *__restrict__ gradQAux,
double *__restrict__ SAux,
const double *__restrict__ QIn,
const double *__restrict__ nodes,
const double *__restrict__ weights,
const double *__restrict__ Kxi,
const double *__restrict__ dudx,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const double dt,
const int nodesPerAxis,
const int unknowns,
const int strideQ,
const int strideRhs,
const bool callFlux,
const bool callSource,
const bool callNonconservativeProduct,
const int scalarIndex )

Add flux, NCP, and source contributions to RHS during STP Picard loop.

This body writes to a single space-time volume coefficient.

Parameters
[in]flux
[in]algebraicSource
[in]nonconservativeProduct
[in]flux
[in]rhsOut
[in]FAux
[in]QIn
[in]nodesquadrature nodes; size: (order+1)
[in]weightsquadrature weights; size: (order+1)
[in]Kxistiffness matrix; size: (order+1)*(order+1)
[in]ttime stamp
[in]dttime step size
[in]centre
[in]dx
[in]unknownsthe number of PDE unknowns that we evolve
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]strideQ
[in]strideRhs
[in]callFlux
[in]callAlgebraicSource
[in]callNonconservativeProduct
[in]scalarIndex

Definition at line 60 of file PredictorAoS.cpp.

References delineariseIndex(), dx, getCoordinates(), getStrides(), and gradient_AoS().

Referenced by spaceTimePredictor_PicardLoop_loop_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ spaceTimePredictor_PicardLoop_addFluxContributionsToRhs_body_AoS()

void exahype2::aderdg::spaceTimePredictor_PicardLoop_addFluxContributionsToRhs_body_AoS ( std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux,
double *__restrict__ rhsOut,
double *__restrict__ FAux,
const double *__restrict__ QIn,
const double *__restrict__ nodes,
const double *__restrict__ weights,
const double *__restrict__ Kxi,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const double dt,
const int nodesPerAxis,
const int unknowns,
const int strideQ,
const int strideRhs,
const int scalarIndex )

Add flux contributions to RHS during STP Picard loop.

This body writes to a single space-time volume coefficient.

Parameters
[in]flux
[in]rhsOut
[in]FAux
[in]QIn
[in]nodesquadrature nodes; size: (order+1)
[in]weightsquadrature weights; size: (order+1)
[in]Kxistiffness matrix; size: (order+1)*(order+1)
[in]ttime stamp
[in]dttime step size
[in]centre
[in]dx
[in]unknownsthe number of PDE unknowns that we evolve
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]strideQ
[in]strideRhs
[in]scalarIndex

Definition at line 158 of file PredictorAoS.cpp.

References delineariseIndex(), dx, getCoordinates(), and getStrides().

Here is the call graph for this function:

◆ spaceTimePredictor_PicardLoop_addNcpContributionToRhs_body_AoS()

void exahype2::aderdg::spaceTimePredictor_PicardLoop_addNcpContributionToRhs_body_AoS ( std::function< void(const double *const __restrict__ Q, double *__restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct,
double *__restrict__ rhsOut,
double *__restrict__ gradQAux,
double *__restrict__ SAux,
const double *const __restrict__ QIn,
const double *const __restrict__ nodes,
const double *const __restrict__ weights,
const double *const __restrict__ dudx,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const double dt,
const int nodesPerAxis,
const int unknowns,
const int strideQ,
const int strideRhs,
const int scalarIndex )

Definition at line 255 of file PredictorAoS.cpp.

References delineariseIndex(), dx, getCoordinates(), getStrides(), and gradient_AoS().

Here is the call graph for this function:

◆ spaceTimePredictor_PicardLoop_addSourceContributionToRhs_body_AoS()

void exahype2::aderdg::spaceTimePredictor_PicardLoop_addSourceContributionToRhs_body_AoS ( std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, double *__restrict__ S) > algebraicSource,
double *__restrict__ rhsOut,
double *__restrict__ SAux,
const double *__restrict__ QIn,
const double *__restrict__ nodes,
const double *__restrict__ weights,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const double dt,
const int nodesPerAxis,
const int unknowns,
const int strideQ,
const int strideRhs,
const int scalarIndex )

Add source contributions to RHS during STP Picard loop.

This body writes to a single space-time volume coefficient.

Parameters
[in]algebraicSource
[in]rhsOut
[in]SAux
[in]QIn
[in]nodesquadrature nodes; size: (order+1)
[in]weightsquadrature weights; size: (order+1)
[in]cellCentre
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]ttime stamp
[in]dttime step size
[in]unknownsthe number of PDE unknowns that we evolve
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]strideQ
[in]strideRhs
[in]scalarIndex

Definition at line 211 of file PredictorAoS.cpp.

References delineariseIndex(), dx, getCoordinates(), and getStrides().

Here is the call graph for this function:

◆ spaceTimePredictor_PicardLoop_initialiseRhs_AoS()

void exahype2::aderdg::spaceTimePredictor_PicardLoop_initialiseRhs_AoS ( double *__restrict__ rhsOut,
const double *const __restrict__ UIn,
const double *const __restrict__ FLCoeff,
const int nodesPerAxis,
const int strideQ,
const int strideRhs,
const int scalarIndex )

Initialise RHS of Picard iterations.

Initialise the right-hand of the equation system that is inverted in every step of the predictor's Picard loop.

This body writes to a single space-time volume coefficient.

Parameters
[in,out]rhsOut
[in]UIn
[in]FLCoeffvalues of basis functions evaluated at x=0.0 (left); size: order+1
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]strideQ

Definition at line 36 of file PredictorAoS.cpp.

References delineariseIndex(), FLCoeff, getStrides(), and it.

Referenced by spaceTimePredictor_PicardLoop_loop_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ spaceTimePredictor_PicardLoop_invert_body_AoS()

void exahype2::aderdg::spaceTimePredictor_PicardLoop_invert_body_AoS ( double *__restrict__ QOut,
double & squaredResiduumOut,
const double *const __restrict__ rhsIn,
const double *const __restrict__ iK1,
const int nodesPerAxis,
const int unknowns,
const int strideQ,
const int strideRhs,
const int scalarIndex )

Invert the Picard space-time system to compute the next solution.


Note
Unlike the kernel bodies that write to the RHS, this kernel cannot be fused with the other space-time predictor Picard loop kernel bodies due to the neighbour lookups of flux and gradient.
Parameters
[in,out]QOut
[in,out]squaredResiduumOut
[in]rhsIn
[in]iK1inverse of the predictor left-hand side operator; size: (order+1)*(order+1)
[in]nodesPerAxisnodes/Lagrange basis functions per coordinate axis (order+1)
[in]unknownsthe number of PDE unknowns that we evolve
[in]strideQ

Definition at line 309 of file PredictorAoS.cpp.

References delineariseIndex(), getStrides(), and it.

Referenced by spaceTimePredictor_PicardLoop_loop_AoS().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ spaceTimePredictor_PicardLoop_loop_AoS()

void exahype2::aderdg::spaceTimePredictor_PicardLoop_loop_AoS ( std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ F) > flux,
std::function< void(const double *const __restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t, double *__restrict__ S) > algebraicSource,
std::function< void(const double *const __restrict__ Q, const double *const __restrict__ dQ_or_deltaQ, const tarch::la::Vector< Dimensions, double > &x, double t, int normal, double *__restrict__ BgradQ) > nonconservativeProduct,
double *__restrict__ QOut,
const double *const __restrict__ UIn,
const double *const __restrict__ nodes,
const double *const __restrict__ weights,
const double *const __restrict__ Kxi,
const double *const __restrict__ iK1,
const double *const __restrict__ FLCoeff,
const double *const __restrict__ dudx,
const tarch::la::Vector< Dimensions, double > & cellCentre,
const double dx,
const double t,
const double dt,
const int order,
const int unknowns,
const int auxiliaryVariables,
const double atol,
const bool callFlux,
const bool callSource,
const bool callNonconservativeProduct )

Compute the space-time predictor (Qout) from the current solution (UIn).

Parameters
[in]flux
[in]algebraicSource
[in]nonconservativeProduct
[in,out]QOut
[in]UIn
[in]weightsquadrature weights; size: (order+1)
[in]nodesquadrature nodes; size: (order+1)
[in]Kxistiffness matrix; size: (order+1)*(order+1)
[in]iK1inverse of the predictor left-hand side operator; size: (order+1)*(order+1)
[in]FLCoeffvalues of basis functions evaluated at x=0.0 (left); size: order+1
[in]dudxderivative operator; size: (order+1)*(order+1)
[in]dxcell spacing (we assume the same spacing in all coordinate directions)
[in]cellCentre
[in]ttime stamp
[in]dttime step size
[in]orderthe DG approximation order, which corresponds to order+1 DG nodes/Lagrange basis functions per coordinate axis
[in]unknownsthe number of PDE unknowns that we evolve
[in]auxiliaryVariablesother quantities such as material parameters that we do not evolve
[in]atol
[in]callFluxcall the flux function
[in]callSourcecall the algebraicSource function
[in]callNonconservativeProductcall the nonconservativeProduct

Definition at line 712 of file PredictorAoS.cpp.

References _log, dx, FLCoeff, getSpaceTimeNodesPerCell(), logWarning, spaceTimePredictor_initialGuess_body_AoS(), spaceTimePredictor_PicardLoop_addContributions_body_AoS(), spaceTimePredictor_PicardLoop_initialiseRhs_AoS(), and spaceTimePredictor_PicardLoop_invert_body_AoS().

Referenced by exahype2::aderdg::tests::ADERDGTest::runADERDGStep().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ validatePatch()

void exahype2::aderdg::validatePatch ( const double *__restrict__ UIn,
int unknowns,
int auxiliaryVariables,
int order,
const std::string & location = "" )

Just runs over the patch and ensures that no entry is non or infinite.

Parameters
[in]UIn
[in]unknownsthe number of PDE unknowns that we evolve
[in]auxiliaryVariablesother quantities such as material parameters that we do not evolve
[in]orderthe DG approximation order, which corresponds to order+1 DG nodes/Lagrange basis functions per coordinate axis
[in]locationString that tells system from where this routine got called

Definition at line 21 of file Generic.cpp.

References dfor, peano4::utils::dLinearised(), k, and nonCriticalAssertion6.

Here is the call graph for this function:

◆ validateSpacetimePatch()

void exahype2::aderdg::validateSpacetimePatch ( const double *__restrict__ QIn,
int unknowns,
int auxiliaryVariables,
int order,
const std::string & location = "" )

Just runs over the space-time patch and ensures that no entry is non or infinite.

Parameters
[in]QIn
[in]unknownsthe number of PDE unknowns that we evolve
[in]auxiliaryVariablesother quantities such as material parameters that we do not evolve
[in]orderthe DG approximation order, which corresponds to order+1 DG nodes/Lagrange basis functions per coordinate axis
[in]locationString that tells system from where this routine got called

Definition at line 37 of file Generic.cpp.

References dfor, peano4::utils::dLinearised(), k, and nonCriticalAssertion6.

Here is the call graph for this function: