25 const char* name =
"GaugeWave";
26 int length = strlen(name);
30 std::cerr <<
"initial scenario " <<
Scenario <<
" is not supported" << std::endl << std::endl << std::endl;
36 double * __restrict__ Q,
54 logError(
"adjustSolution(...)",
"initial scenario " <<
Scenario <<
" is not supported" );
57 for (
int i=0; i<NumberOfUnknowns; i++) {
58 assertion3( std::isfinite(Q[i]), volumeX, t, i );
61 for (
int i=NumberOfUnknowns; i<NumberOfUnknowns+NumberOfAuxiliaryVariables; i++) {
73 const double * __restrict__ Q,
78 double * __restrict__ S
81 for(
int i=0; i<NumberOfUnknowns; i++){
85 memset(S, 0, NumberOfUnknowns*
sizeof(
double));
88 for(
int i=0; i<NumberOfUnknowns; i++){
97 const double * __restrict__ Qinside,
98 double * __restrict__ Qoutside,
105 for(
int i=0; i<NumberOfUnknowns+NumberOfAuxiliaryVariables; i++) {
106 assertion4( Qinside[i]==Qinside[i], faceCentre, t, normal, i );
107 Qoutside[i]=Qinside[i];
114 const double * __restrict__ Q,
121#if defined(MGCCZ4EINSTEIN)
122 const double qmin = std::min({Q[0],Q[3],Q[5]});
123 const double alpha = std::max({1.0, std::exp(Q[16])}) * std::max({1.0, std::exp(Q[54])}) / std::sqrt(qmin);
125 const double alpha = 1.0;
128 constexpr double sqrtwo = 1.4142135623730951;
130 const double tempA = alpha * std::max({sqrtwo, MGCCZ4e, MGCCZ4ds, MGCCZ4GLMc/alpha, MGCCZ4GLMd/alpha});
131 const double tempB = Q[17+normal];
132 return std::max({1.0, std::abs(-tempA-tempB), std::abs(tempA-tempB)});
166 const double * __restrict__ Q,
167 const double * __restrict__ deltaQ,
172 double * __restrict__ BgradQ
174#if !defined(GPUOffloadingOMP)
179 double gradQSerialised[NumberOfUnknowns*3];
180 for (
int i=0; i<NumberOfUnknowns; i++) {
181 gradQSerialised[i+0*NumberOfUnknowns] = 0.0;
182 gradQSerialised[i+1*NumberOfUnknowns] = 0.0;
183 gradQSerialised[i+2*NumberOfUnknowns] = 0.0;
185 gradQSerialised[i+normal*NumberOfUnknowns] = deltaQ[i];
188 pdencp_(BgradQ, Q, gradQSerialised);
189#if !defined(GPUOffloadingOMP)
190 for (
int i=0; i<NumberOfUnknowns; i++) {
#define assertion4(expr, param0, param1, param2, param3)
#define assertion3(expr, param0, param1, param2)
#define logError(methodName, logMacroMessageStream)
Wrapper macro around tarch::tarch::logging::Log to improve logging.
#define logTraceOut(methodName)
#define logTraceInWith4Arguments(methodName, argument0, argument1, argument2, argument3)
#define nonCriticalAssertion3(expr, param0, param1, param2)
#define nonCriticalAssertion4(expr, param0, param1, param2, param3)
virtual void boundaryConditions(const double *__restrict__ Qinside, double *__restrict__ Qoutside, const tarch::la::Vector< Dimensions, double > &faceCentre, const tarch::la::Vector< Dimensions, double > &volumeH, double t, int normal) override
void adjustSolution(double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &volumeCentre, const tarch::la::Vector< Dimensions, double > &volumeH, double t, double dt) override
double maxEigenvalue(const double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &faceCentre, const tarch::la::Vector< Dimensions, double > &volumeH, double t, int normal)
void sourceTerm(const double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &volumeCentre, const tarch::la::Vector< Dimensions, double > &volumeH, double t, double dt, double *__restrict__ S)
static tarch::logging::Log _log
void nonconservativeProduct(const double *__restrict__ Q, const double *__restrict__ deltaQ, const tarch::la::Vector< Dimensions, double > &faceCentre, const tarch::la::Vector< Dimensions, double > &volumeH, double t, int normal, double *__restrict__ BgradQ)
void pdesource_(double *S, const double *const Q)
void pdencp_(double *BgradQ, const double *const Q, const double *const gradQ)
void enforceMGCCZ4constraints(double *luh)
void linearWave(double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &X, double t)
void forcedflat(double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &X, double t)
void gaugeWave(double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t)
bool equals(const Matrix< Rows, Cols, Scalar > &lhs, const Matrix< Rows, Cols, Scalar > &rhs, const Scalar &tolerance=NUMERICAL_ZERO_DIFFERENCE)
Compares to matrices on equality by means of a numerical accuracy.