Peano 4
Loading...
Searching...
No Matches
exahype2::dg::rusanov Namespace Reference

Namespaces

namespace  omp
 
namespace  tests
 

Functions

void solveRiemannProblem_pointwise_in_situ (::exahype2::dg::Flux flux, ::exahype2::dg::NonConservativeProduct ncp, ::exahype2::dg::MaxEigenvalue maxEigenvalue, const tarch::la::Vector< Dimensions, double > &faceCentre, const tarch::la::Vector< Dimensions, double > &cellSize, double t, double dt, int order, int unknowns, int auxiliaryVariables, int faceNumber, const double *__restrict__ quadraturePoints, bool useFlux, bool useNcp, const double *__restrict__ projectedValues, double *__restrict__ solution)
 Solve Riemann problem on face.
 
void solveRiemannProblem_pointwise_in_situ_with_gradient_projection (::exahype2::dg::Flux flux, ::exahype2::dg::NonConservativeProduct ncp, ::exahype2::dg::MaxEigenvalue maxEigenvalue, const tarch::la::Vector< Dimensions, double > &faceCentre, const tarch::la::Vector< Dimensions, double > &cellSize, double t, double dt, int order, int unknowns, int auxiliaryVariables, int faceNumber, const double *__restrict__ quadraturePoints, bool useFlux, bool useNcp, const double *__restrict__ projectedValues, double *__restrict__ solution)
 
static void cellIntegral_patchwise_in_situ_GaussLegendre_functors (::exahype2::CellData &cellData, const int order, const int unknowns, const int auxiliaryVariables, Flux flux, NonConservativeProduct nonConservativeProduct, Source source, PointSources pointSources, const double *__restrict__ QuadratureNodes1d, const double *__restrict__ MassMatrixDiagonal1d, const double *__restrict__ StiffnessMatrix1d, const double *__restrict__ DerivativeOperator1d, bool evaluateFlux, bool evaluateNonconservativeProduct, bool evaluateSource, bool evaluatePointSources)
 Delegate to generic implementation in parent directory.
 
static void multiplyWithInvertedMassMatrix_GaussLegendre (::exahype2::CellData &cellData, const int order, const int unknowns, const int auxiliaryVariables, const double *__restrict__ MassMatrixDiagonal1d)
 Delegate to generic implementation.
 
static void integrateOverRiemannSolutionsAndAddToVolume_GaussLegendre (const double *const __restrict__ faceQLeft, const double *const __restrict__ faceQRight, const double *const __restrict__ faceQBottom, const double *const __restrict__ faceQUp, int order, int unknowns, const int auxiliaryVariables, const tarch::la::Vector< 2, double > &cellSize, const double *const __restrict__ BasisFunctionValuesLeft, const double *__restrict__ MassMatrixDiagonal1d, double *__restrict__ cellQ)
 Delegate to generic implementation.
 
static void integrateOverRiemannSolutionsAndAddToVolume_GaussLegendre (const double *const __restrict__ faceQLeft, const double *const __restrict__ faceQRight, const double *const __restrict__ faceQBottom, const double *const __restrict__ faceQUp, const double *const __restrict__ faceQFront, const double *const __restrict__ faceQBack, int order, int unknowns, const int auxiliaryVariables, const tarch::la::Vector< 3, double > &cellSize, const double *const __restrict__ BasisFunctionValuesLeft, const double *__restrict__ MassMatrixDiagonal1d, double *__restrict__ cellQ)
 Delegate to generic implementation.
 
template<typename Solver , int order, int unknowns, int auxiliaryVariables>
void cellIntegral_patchwise_in_situ_GaussLegendre (::exahype2::CellData &cellData, bool evaluateFlux, bool evaluateNonconservativeProduct, bool evaluateSource, bool evaluatePointSources)
 

Function Documentation

◆ cellIntegral_patchwise_in_situ_GaussLegendre()

template<typename Solver , int order, int unknowns, int auxiliaryVariables>
void exahype2::dg::rusanov::cellIntegral_patchwise_in_situ_GaussLegendre ( ::exahype2::CellData & cellData,
bool evaluateFlux,
bool evaluateNonconservativeProduct,
bool evaluateSource,
bool evaluatePointSources )

Definition at line 265 of file Rusanov.h.

◆ cellIntegral_patchwise_in_situ_GaussLegendre_functors()

static void exahype2::dg::rusanov::cellIntegral_patchwise_in_situ_GaussLegendre_functors ( ::exahype2::CellData & cellData,
const int order,
const int unknowns,
const int auxiliaryVariables,
Flux flux,
NonConservativeProduct nonConservativeProduct,
Source source,
PointSources pointSources,
const double *__restrict__ QuadratureNodes1d,
const double *__restrict__ MassMatrixDiagonal1d,
const double *__restrict__ StiffnessMatrix1d,
const double *__restrict__ DerivativeOperator1d,
bool evaluateFlux,
bool evaluateNonconservativeProduct,
bool evaluateSource,
bool evaluatePointSources )
static

Delegate to generic implementation in parent directory.

static keyword required to avoid multiple definition linker errors.

Definition at line 149 of file Rusanov.h.

References exahype2::dg::cellIntegral_patchwise_in_situ_GaussLegendre_functors().

Here is the call graph for this function:

◆ integrateOverRiemannSolutionsAndAddToVolume_GaussLegendre() [1/2]

static void exahype2::dg::rusanov::integrateOverRiemannSolutionsAndAddToVolume_GaussLegendre ( const double *const __restrict__ faceQLeft,
const double *const __restrict__ faceQRight,
const double *const __restrict__ faceQBottom,
const double *const __restrict__ faceQUp,
const double *const __restrict__ faceQFront,
const double *const __restrict__ faceQBack,
int order,
int unknowns,
const int auxiliaryVariables,
const tarch::la::Vector< 3, double > & cellSize,
const double *const __restrict__ BasisFunctionValuesLeft,
const double *__restrict__ MassMatrixDiagonal1d,
double *__restrict__ cellQ )
static

Delegate to generic implementation.

static keyword required to avoid multiple definition linker errors.

Definition at line 233 of file Rusanov.h.

References exahype2::dg::integrateOverRiemannSolutionsAndAddToVolume_GaussLegendre().

Here is the call graph for this function:

◆ integrateOverRiemannSolutionsAndAddToVolume_GaussLegendre() [2/2]

static void exahype2::dg::rusanov::integrateOverRiemannSolutionsAndAddToVolume_GaussLegendre ( const double *const __restrict__ faceQLeft,
const double *const __restrict__ faceQRight,
const double *const __restrict__ faceQBottom,
const double *const __restrict__ faceQUp,
int order,
int unknowns,
const int auxiliaryVariables,
const tarch::la::Vector< 2, double > & cellSize,
const double *const __restrict__ BasisFunctionValuesLeft,
const double *__restrict__ MassMatrixDiagonal1d,
double *__restrict__ cellQ )
static

Delegate to generic implementation.

static keyword required to avoid multiple definition linker errors.

Definition at line 204 of file Rusanov.h.

References exahype2::dg::integrateOverRiemannSolutionsAndAddToVolume_GaussLegendre().

Here is the call graph for this function:

◆ multiplyWithInvertedMassMatrix_GaussLegendre()

static void exahype2::dg::rusanov::multiplyWithInvertedMassMatrix_GaussLegendre ( ::exahype2::CellData & cellData,
const int order,
const int unknowns,
const int auxiliaryVariables,
const double *__restrict__ MassMatrixDiagonal1d )
static

Delegate to generic implementation.

static keyword required to avoid multiple definition linker errors.

Definition at line 185 of file Rusanov.h.

References exahype2::dg::multiplyWithInvertedMassMatrix_GaussLegendre().

Here is the call graph for this function:

◆ solveRiemannProblem_pointwise_in_situ()

void exahype2::dg::rusanov::solveRiemannProblem_pointwise_in_situ ( ::exahype2::dg::Flux flux,
::exahype2::dg::NonConservativeProduct ncp,
::exahype2::dg::MaxEigenvalue maxEigenvalue,
const tarch::la::Vector< Dimensions, double > & faceCentre,
const tarch::la::Vector< Dimensions, double > & cellSize,
double t,
double dt,
int order,
int unknowns,
int auxiliaryVariables,
int faceNumber,
const double *__restrict__ quadraturePoints,
bool useFlux,
bool useNcp,
const double *__restrict__ projectedValues,
double *__restrict__ solution )

Solve Riemann problem on face.

This routine solves the Riemann problem on one face through a Rusanov solver. It takes the average of the flux left and right, reduces it by the jump subject to a scaling with the maximum eigenvalue, and finally adds/substracts the ncp terms.

The implementation is purely point-wise, i.e. we run over the quadrature points and use Rusanov per point. It is not clear if this is reasonable in all cases. It might be better for example to determine the max eigenvalue per face once and then update all quadrature point subject to this one eigenvalue.

Implementation details

The DoFs per face are ordered in standard Peano order:

    2 3
    0 1

2 3 [ o o ] 2 3 0 1 [ o o ] 0 1 2 3 0 1

The lrEnumerator is used for quantities that we evaluate right and left of the face. The faceEnumerator is for data that we evaluate directly on the face.

Linearity in h

The ncp - as so often - deserves special attention: As we have a non-conservative term here, we don't write the same flux to the left and the right face. Indeed, they carry exactly switched signs. We furthermore note that the original PDE formulation speaks of a term \( B(Q) \nabla B \) whereas the Finite Volume formulation's Rusanov solver uses a \( 0.5 B( 0.5(Q^++Q^-) ) (Q^+-Q^-) \) term. Therefore, we need an additional factor of 0.5. However, we do not need any additional calibration, as we simply pass in the Q differences into the DG ncp term which expects the gradient. Whatever comes out of it, should be by a factor of h too big already, i.e. have the right cardinality of the FV expression.

In theory, it might be that B is not linear in the vector, i.e. not linear in an h scaling. However, we always can assume that the B that is passed in here is a linearisation along the face. Therefore, we can assume that the scaling is by definition fine with our observation.

Sign discussion

We rely on the formulation as discussed on the page Discontinuous Galerkin. See dg.h or the doxygen html pages for details. I also recommend to read through the documentation of cellIntegral_patchwise_in_situ_GaussLegendre_functors() before we dive into the discussion.

All the comments below assume that both the left and the right side of the flux compute the flux in the normal direction. It is then the backpropagation of the Riemann solution which takes into the account that the right'' faces denote an outflow, while theleft'' faces encode an inflow.

See also
integrateOverRiemannSolutionsAndAddToVolume_GaussLegendre()

Bringing all these facts together, we conclude that we end up with a 1:1 translation of Finite Volumes into the DG world.

See also
exahype2/fv/rusanov/LoopBody.cpph for further information.

Definition at line 6 of file Rusanov.cpp.

References dfore.

◆ solveRiemannProblem_pointwise_in_situ_with_gradient_projection()

void exahype2::dg::rusanov::solveRiemannProblem_pointwise_in_situ_with_gradient_projection ( ::exahype2::dg::Flux flux,
::exahype2::dg::NonConservativeProduct ncp,
::exahype2::dg::MaxEigenvalue maxEigenvalue,
const tarch::la::Vector< Dimensions, double > & faceCentre,
const tarch::la::Vector< Dimensions, double > & cellSize,
double t,
double dt,
int order,
int unknowns,
int auxiliaryVariables,
int faceNumber,
const double *__restrict__ quadraturePoints,
bool useFlux,
bool useNcp,
const double *__restrict__ projectedValues,
double *__restrict__ solution )
See also
solveRiemannProblem_pointwise_in_situ()

Definition at line 127 of file Rusanov.cpp.

References assertionMsg.