57 for (
size_t i = 0; i <
m_fields.size(); i++)
68 "No UPWINDTYPE defined in session.");
100 for (
size_t n = 0; n < (size_t)
m_fields[0]->GetBndConditions().size(); ++n)
102 std::string type =
m_fields[0]->GetBndConditions()[n]->GetUserDefined();
104 if (
m_fields[0]->GetBndConditions()[n]->GetBoundaryConditionType() ==
115 cnt +=
m_fields[0]->GetBndCondExpansions()[n]->GetExpSize();
137 " is not supported by this solver";
156 "Local time stepping requires CFL parameter.");
167 "Unsupported projection type.");
169 string advName, riemName;
170 m_session->LoadSolverInfo(
"AdvectionType", advName,
"WeakDG");
187 m_session->LoadSolverInfo(
"UpwindType", riemName,
"Average");
213 size_t nvariables = inarray.size();
230 for (
size_t i = 0; i < nvariables; ++i)
234 m_fields[i]->GetFwdBwdTracePhys(inarray[i], Fwd[i], Bwd[i]);
245 for (
size_t i = 0; i < nvariables; ++i)
259 x->Apply(
m_fields, inarray, outarray, time);
264 size_t nElements =
m_fields[0]->GetExpSize();
273 for (
size_t n = 0; n < nElements; ++n)
275 nq =
m_fields[0]->GetExp(n)->GetTotPoints();
276 offset =
m_fields[0]->GetPhys_Offset(n);
278 for (
size_t i = 0; i < nvariables; ++i)
281 tmp = outarray[i] + offset, 1);
300 size_t nvariables = inarray.size();
303 for (
size_t i = 0; i < nvariables; ++i)
320 "compressible Navier-Stokes equations");
338 int nvariables = inarray.size();
361 size_t nvariables = physarray.size();
364 for (
size_t i = 0; i < nvariables; ++i)
367 m_fields[i]->ExtractTracePhys(physarray[i], Fwd[i]);
375 x->Apply(Fwd, physarray, time);
405 size_t nVariables = physfield.size();
406 size_t nPts = physfield[0].size();
408 constexpr unsigned short maxVel = 3;
409 constexpr unsigned short maxFld = 5;
412 ASSERTL1(nVariables <= maxFld,
"GetFluxVector, hard coded max fields");
414 for (
size_t p = 0;
p < nPts; ++
p)
417 std::array<NekDouble, maxFld> fieldTmp;
418 std::array<NekDouble, maxVel>
velocity;
421 for (
size_t f = 0; f < nVariables; ++f)
423 fieldTmp[f] = physfield[f][
p];
432 flux[0][
d][
p] = fieldTmp[
d + 1];
468 size_t nq = physfield[0].size();
469 size_t nVariables =
m_fields.size();
473 nq =
m_fields[0]->Get1DScaledTotPoints(OneDptscale);
481 for (
size_t i = 0; i < nVariables; ++i)
485 m_fields[0]->PhysInterp1DScaled(OneDptscale, physfield[i],
486 physfield_interp[i]);
500 m_fields[0]->PhysGalerkinProjection1DScaled(
501 OneDptscale, physfield_interp[i + 1], flux[0][i]);
513 flux_interp[i + 1][j], 1);
518 flux_interp[i + 1][i], 1);
526 m_fields[0]->PhysGalerkinProjection1DScaled(
527 OneDptscale, flux_interp[i + 1][j], flux[i + 1][j]);
541 m_fields[0]->PhysGalerkinProjection1DScaled(
555 size_t nElements =
m_fields[0]->GetExpSize();
568 for (
size_t n = 0; n < nElements; ++n)
583 int nElements =
m_fields[0]->GetExpSize();
598 if (
m_session->DefinesSolverInfo(
"ICTYPE"))
600 if (boost::iequals(
m_session->GetSolverInfo(
"ICTYPE"),
605 else if (boost::iequals(
m_session->GetSolverInfo(
"ICTYPE"),
621 NekDouble initialtime,
bool dumpInitialConditions,
622 [[maybe_unused]]
const int domain)
624 if (
m_session->DefinesSolverInfo(
"ICTYPE") &&
625 boost::iequals(
m_session->GetSolverInfo(
"ICTYPE"),
"IsentropicVortex"))
628 for (
int i = 0; i <
m_fields.size(); ++i)
643 int phystot =
m_fields[0]->GetTotPoints();
646 m_session->LoadParameter(
"Noise", Noise, 0.0);
647 int m_nConvectiveFields =
m_fields.size();
651 int seed = -
m_comm->GetSpaceComm()->GetRank() * m_nConvectiveFields;
652 for (
int i = 0; i < m_nConvectiveFields; i++)
665 !
m_comm->IsParallelInTime())
669 else if (dumpInitialConditions &&
m_nchk == 0 &&
m_comm->IsParallelInTime())
672 if (!fs::is_directory(newdir))
674 fs::create_directory(newdir);
676 if (
m_comm->GetTimeComm()->GetRank() == 0)
688 if (!
m_session->DefinesFunction(
"ExactSolution") &&
691 if (boost::iequals(
m_session->GetSolverInfo(
"ICTYPE"),
696 else if (boost::iequals(
m_session->GetSolverInfo(
"ICTYPE"),
716 size_t n_element =
m_fields[0]->GetExpSize();
717 size_t expdim =
m_fields[0]->GetGraph()->GetMeshDimension();
723 for (
size_t i = 0; i < nfields; ++i)
725 physfields[i] =
m_fields[i]->GetPhys();
745 m_varConv->GetSoundSpeed(physfields, soundspeed);
747 for (
size_t el = 0; el < n_element; ++el)
749 ptsKeys =
m_fields[0]->GetExp(el)->GetPointsKeys();
750 offset =
m_fields[0]->GetPhys_Offset(el);
751 int nq =
m_fields[0]->GetExp(el)->GetTotPoints();
754 m_fields[0]->GetExp(el)->GetGeom()->GetMetricInfo();
760 ->GetDerivFactors(ptsKeys);
768 for (
size_t i = 0; i < expdim; ++i)
771 tmp = stdVelocity[i] + offset, 1);
772 Vmath::Vmul(nq, gmat[i], 1, soundspeed + offset, 1,
773 tmp = stdSoundSpeed[i] + offset, 1);
774 for (
size_t j = 1; j < expdim; ++j)
778 stdVelocity[i] + offset, 1,
779 tmp = stdVelocity[i] + offset, 1);
781 soundspeed + offset, 1,
782 stdSoundSpeed[i] + offset, 1,
783 tmp = stdSoundSpeed[i] + offset, 1);
789 for (
size_t i = 0; i < expdim; ++i)
792 tmp = stdVelocity[i] + offset, 1);
793 Vmath::Smul(nq, gmat[i][0], soundspeed + offset, 1,
794 tmp = stdSoundSpeed[i] + offset, 1);
795 for (
size_t j = 1; j < expdim; ++j)
799 stdVelocity[i] + offset, 1,
800 tmp = stdVelocity[i] + offset, 1);
802 soundspeed + offset, 1,
803 stdSoundSpeed[i] + offset, 1,
804 tmp = stdSoundSpeed[i] + offset, 1);
810 for (
size_t i = 0; i < nq; ++i)
813 for (
size_t j = 0; j < expdim; ++j)
816 vel =
std::abs(stdVelocity[j][offset + i]) +
817 SpeedSoundFactor *
std::abs(stdSoundSpeed[j][offset + i]);
818 pntVelocity += vel * vel;
820 pntVelocity =
sqrt(pntVelocity);
821 if (pntVelocity > stdV[el])
823 stdV[el] = pntVelocity;
840 ASSERTL0(n <= 20,
"Illegal modes dimension for CFL calculation "
841 "(P has to be less then 20)");
843 NekDouble CFLDG[21] = {2.0000, 6.0000, 11.8424, 19.1569, 27.8419,
844 37.8247, 49.0518, 61.4815, 75.0797, 89.8181,
845 105.6700, 122.6200, 140.6400, 159.7300, 179.8500,
846 201.0100, 223.1800, 246.3600, 270.5300, 295.6900,
857 "coefficients not introduced yet.");
875 for (i = 0; i <
m_fields[0]->GetExpSize(); i++)
884 std::vector<std::string> &variables)
887 m_session->MatchSolverInfo(
"OutputExtraFields",
"True", extraFields,
true);
890 const int nPhys =
m_fields[0]->GetNpoints();
891 const int nCoeffs =
m_fields[0]->GetNcoeffs();
894 for (
int i = 0; i <
m_fields.size(); ++i)
914 m_varConv->GetTemperature(tmp, temperature);
916 m_varConv->GetSoundSpeed(tmp, soundspeed);
917 m_varConv->GetMach(tmp, soundspeed, mach);
920 m_session->LoadParameter(
"SensorOffset", sensorOffset, 1);
929 string velNames[3] = {
"u",
"v",
"w"};
933 variables.push_back(velNames[i]);
934 fieldcoeffs.push_back(velFwd[i]);
938 m_fields[0]->FwdTransLocalElmt(temperature, TFwd);
939 m_fields[0]->FwdTransLocalElmt(entropy, sFwd);
940 m_fields[0]->FwdTransLocalElmt(soundspeed, aFwd);
941 m_fields[0]->FwdTransLocalElmt(mach, mFwd);
942 m_fields[0]->FwdTransLocalElmt(sensor, sensFwd);
944 variables.push_back(
"p");
945 variables.push_back(
"T");
946 variables.push_back(
"s");
947 variables.push_back(
"a");
948 variables.push_back(
"Mach");
949 variables.push_back(
"Sensor");
950 fieldcoeffs.push_back(pFwd);
951 fieldcoeffs.push_back(TFwd);
952 fieldcoeffs.push_back(sFwd);
953 fieldcoeffs.push_back(aFwd);
954 fieldcoeffs.push_back(mFwd);
955 fieldcoeffs.push_back(sensFwd);
964 variables.push_back(
"ArtificialVisc");
965 fieldcoeffs.push_back(sensorFwd);
987 density = physfield[0];
1004 const size_t nFields =
m_fields.size();
1007 for (
size_t i = 0; i < nFields; ++i)
1010 inarray[i] =
m_fields[i]->UpdatePhys();
1020 for (
size_t i = 0; i < nFields; ++i)
1031 for (
size_t i = 0; i < nFields; ++i)
1033 L2[i] =
sqrt(residual[i] * onPoints);
1059 m_session->LoadParameter(
"IsentropicU0", u0, 1.0);
1060 m_session->LoadParameter(
"IsentropicV0", v0, 0.5);
1061 m_session->LoadParameter(
"IsentropicX0", x0, 5.0);
1062 m_session->LoadParameter(
"IsentropicY0", y0, 0.0);
1075 for (
int i = 0; i < nTotQuadPoints; ++i)
1077 xbar = x[i] - u0 * time - x0;
1078 ybar = y[i] - v0 * time - y0;
1079 r =
sqrt(xbar * xbar + ybar * ybar);
1080 tmp =
beta * exp(1 - r * r);
1082 pow(1.0 - (
m_gamma - 1.0) * tmp * tmp * fac, 1.0 / (
m_gamma - 1.0));
1083 u[1][i + o] = u[0][i + o] * (u0 - tmp * ybar / (2 * M_PI));
1084 u[2][i + o] = u[0][i + o] * (v0 + tmp * xbar / (2 * M_PI));
1087 0.5 * (u[1][i + o] * u[1][i + o] + u[2][i + o] * u[2][i + o]) /
1090 Vmath::Vcopy(nTotQuadPoints, u[field].get(), 1, outfield.get(), 1);
1112 NekDouble c, k, phi, r, J, VV, pp, sint,
P, ss;
1126 NekDouble gamma_1_2 = (gamma - 1.0) / 2.0;
1128 for (
int i = 0; i < nTotQuadPoints; ++i)
1130 while ((
abs(errV) > toll) || (
abs(errTheta) > toll))
1134 c =
sqrt(1.0 - gamma_1_2 * VV);
1138 J = 1.0 / c + 1.0 / (3.0 * c * c * c) +
1139 1.0 / (5.0 * c * c * c * c * c) -
1140 0.5 *
log((1.0 + c) / (1.0 - c));
1142 r = pow(c, 1.0 / gamma_1_2);
1143 xi = 1.0 / (2.0 * r) * (1.0 / VV - 2.0 * pp) + J / 2.0;
1144 yi = phi / (r * V) *
sqrt(1.0 - VV * pp);
1145 par1 = 25.0 - 5.0 * VV;
1152 39062.5 / pow(par1, 3.5) * (1.0 / VV - 2.0 / VV * ss) * V +
1153 1562.5 / pow(par1, 2.5) *
1154 (-2.0 / (VV * V) + 4.0 / (VV * V) * ss) +
1155 12.5 / pow(par1, 1.5) * V + 312.5 / pow(par1, 2.5) * V +
1156 7812.5 / pow(par1, 3.5) * V -
1158 (-1.0 / pow(par1, 0.5) * V / (1.0 - 0.2 * pow(par1, 0.5)) -
1159 (1.0 + 0.2 * pow(par1, 0.5)) /
1160 pow((1.0 - 0.2 * pow(par1, 0.5)), 2.0) /
1161 pow(par1, 0.5) * V) /
1162 (1.0 + 0.2 * pow(par1, 0.5)) * (1.0 - 0.2 * pow(par1, 0.5));
1164 J12 = -6250.0 / pow(par1, 2.5) / VV * sint * cos(theta);
1165 J21 = -6250.0 / (VV * V) * sint / pow(par1, 2.5) *
1166 pow((1.0 - ss), 0.5) +
1167 78125.0 / V * sint / pow(par1, 3.5) * pow((1.0 - ss), 0.5);
1170 if (
abs(y[i]) < toll &&
abs(cos(theta)) < toll)
1172 J22 = -39062.5 / pow(par1, 3.5) / V +
1173 3125 / pow(par1, 2.5) / (VV * V) +
1174 12.5 / pow(par1, 1.5) * V + 312.5 / pow(par1, 2.5) * V +
1175 7812.5 / pow(par1, 3.5) * V -
1177 (-1.0 / pow(par1, 0.5) * V /
1178 (1.0 - 0.2 * pow(par1, 0.5)) -
1179 (1.0 + 0.2 * pow(par1, 0.5)) /
1180 pow((1.0 - 0.2 * pow(par1, 0.5)), 2.0) /
1181 pow(par1, 0.5) * V) /
1182 (1.0 + 0.2 * pow(par1, 0.5)) *
1183 (1.0 - 0.2 * pow(par1, 0.5));
1186 dV = -1.0 / J22 * Fx;
1192 J22 = 3125.0 / VV * cos(theta) / pow(par1, 2.5) *
1193 pow((1.0 - ss), 0.5) -
1194 3125.0 / VV * ss / pow(par1, 2.5) / pow((1.0 - ss), 0.5) *
1197 det = -1.0 / (J11 * J22 - J12 * J21);
1200 dV = det * (J22 * Fx - J12 * Fy);
1201 dtheta = det * (-J21 * Fx + J11 * Fy);
1205 theta = theta + dtheta;
1208 errTheta =
abs(dtheta);
1211 c =
sqrt(1.0 - gamma_1_2 * VV);
1212 r = pow(c, 1.0 / gamma_1_2);
1215 rhou[i] = rho[i] * V * cos(theta);
1216 rhov[i] = rho[i] * V * sin(theta);
1217 P = (c * c) * rho[i] / gamma;
1218 E[i] =
P / (gamma - 1.0) +
1219 0.5 * (rhou[i] * rhou[i] / rho[i] + rhov[i] * rhov[i] / rho[i]);
1225 V = kExt * sin(theta);
1243 ASSERTL0(
false,
"Error in variable number!");
#define ASSERTL0(condition, msg)
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
void GetFluxVector(const Array< OneD, const Array< OneD, NekDouble > > &physfield, TensorOfArray3D< NekDouble > &flux)
Return the flux vector for the compressible Euler equations.
CompressibleFlowSystem(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
void GetElmtTimeStep(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &tstep)
Calculate the maximum timestep on each element subject to CFL restrictions.
void SetBoundaryConditions(Array< OneD, Array< OneD, NekDouble > > &physarray, NekDouble time)
void DoOdeProjection(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
Compute the projection and call the method for imposing the boundary conditions in case of discontinu...
NekDouble m_bndEvaluateTime
virtual bool v_SupportsShockCaptType(const std::string type) const =0
std::string m_shockCaptureType
MultiRegions::ExpListSharedPtr v_GetPressure() override
Array< OneD, NekDouble > GetStabilityLimitVector(const Array< OneD, int > &ExpOrder)
Function to calculate the stability limit for DG/CG (a vector of them).
void InitAdvection()
Create advection and diffusion objects for CFS.
NekDouble v_GetTimeStep(const Array< OneD, const Array< OneD, NekDouble > > &inarray) override
Calculate the maximum timestep subject to CFL restrictions.
const Array< OneD, const Array< OneD, NekDouble > > & GetNormals()
void v_GenerateSummary(SolverUtils::SummaryList &s) override
Print a summary of time stepping parameters.
virtual void v_DoDiffusion(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd)=0
void v_SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0) override
Set up logic for residual calculation.
void GetExactRinglebFlow(int field, Array< OneD, NekDouble > &outarray)
Ringleb Flow Test Case.
NekDouble m_filterExponent
void v_InitObject(bool DeclareFields=true) override
Initialization object for CompressibleFlowSystem class.
void v_GetDensity(const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &density) override
void EvaluateIsentropicVortex(unsigned int field, Array< OneD, NekDouble > &outfield, NekDouble time, const int o=0)
Isentropic Vortex Test Case.
NekDouble GetStabilityLimit(int n)
Function to calculate the stability limit for DG/CG.
void GetFluxVectorDeAlias(const Array< OneD, const Array< OneD, NekDouble > > &physfield, TensorOfArray3D< NekDouble > &flux)
Return the flux vector for the compressible Euler equations by using the de-aliasing technique.
std::vector< CFSBndCondSharedPtr > m_bndConds
Array< OneD, Array< OneD, NekDouble > > m_vecLocs
void v_ExtraFldOutput(std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables) override
VariableConverterSharedPtr m_varConv
const Array< OneD, const Array< OneD, NekDouble > > & GetVecLocs()
void SetBoundaryConditionsBwdWeight()
Set up a weight on physical boundaries for boundary condition applications.
std::vector< SolverUtils::ForcingSharedPtr > m_forcing
ArtificialDiffusionSharedPtr m_artificialDiffusion
Array< OneD, NekDouble > v_GetMaxStdVelocity(const NekDouble SpeedSoundFactor) override
Compute the advection velocity in the standard space for each element of the expansion.
void DoAdvection(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd)
Compute the advection terms for the right-hand side.
void InitialiseParameters()
Load CFS parameters from the session file.
void DoDiffusion(const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd)
Add the diffusions terms to the right-hand side.
void DoOdeRhs(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
Compute the right-hand side.
void v_SteadyStateResidual(int step, Array< OneD, NekDouble > &L2) override
void v_EvaluateExactSolution(unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time=0.0) override
void v_GetVelocity(const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &velocity) override
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
void DefineProjection(FuncPointerT func, ObjectPointerT obj)
void DefineOdeRhs(FuncPointerT func, ObjectPointerT obj)
void AccumulateRegion(std::string, int iolevel=0)
Accumulate elapsed time for a region.
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
A base class for PDEs which include an advection component.
SolverUtils::AdvectionSharedPtr m_advObject
Advection term.
SOLVER_UTILS_EXPORT Array< OneD, NekDouble > GetElmtCFLVals(const bool FlagAcousticCFL=true)
SOLVER_UTILS_EXPORT void v_InitObject(bool DeclareField=true) override
Initialisation object for EquationSystem.
int m_spacedim
Spatial dimension (>= expansion dim).
virtual SOLVER_UTILS_EXPORT void v_SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
SpatialDomains::MeshGraphSharedPtr m_graph
Pointer to graph defining mesh.
LibUtilities::CommSharedPtr m_comm
Communicator.
NekDouble m_timestep
Time step size.
NekDouble m_time
Current time of simulation.
SOLVER_UTILS_EXPORT int GetTraceTotPoints()
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
SOLVER_UTILS_EXPORT void Checkpoint_Output(const int n)
Write checkpoint file of m_fields.
virtual SOLVER_UTILS_EXPORT void v_EvaluateExactSolution(unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time)
SOLVER_UTILS_EXPORT void WriteFld(const std::string &outname)
Write field data to the given filename.
SOLVER_UTILS_EXPORT int GetExpSize()
std::string m_sessionName
Name of the session.
bool m_specHP_dealiasing
Flag to determine if dealisising is usde for the Spectral/hp element discretisation.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Array holding trace normals for DG simulations in the forwards direction.
enum HomogeneousType m_HomogeneousType
SOLVER_UTILS_EXPORT int GetNpoints()
int m_nchk
Number of checkpoints written so far.
enum MultiRegions::ProjectionType m_projectionType
Type of projection; e.g continuous or discontinuous.
int m_checksteps
Number of steps between checkpoints.
SOLVER_UTILS_EXPORT int GetTotPoints()
static SOLVER_UTILS_EXPORT std::vector< ForcingSharedPtr > Load(const LibUtilities::SessionReaderSharedPtr &pSession, const std::weak_ptr< EquationSystem > &pEquation, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const unsigned int &pNumForcingFields=0)
Base class for unsteady solvers.
LibUtilities::TimeIntegrationSchemeOperators m_ode
The time integration scheme operators to use.
NekDouble m_cflSafetyFactor
CFL safety factor (comprise between 0 to 1).
SOLVER_UTILS_EXPORT NekDouble MaxTimeStepEstimator()
Get the maximum timestep estimator for cfl control.
SOLVER_UTILS_EXPORT void v_GenerateSummary(SummaryList &s) override
Print a summary of time stepping parameters.
bool m_homoInitialFwd
Flag to determine if simulation should start in homogeneous forward transformed state.
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::vector< PointsKey > PointsKeyVector
@ beta
Gauss Radau pinned at x=-1,.
@ P
Monomial polynomials .
const std::vector< NekDouble > velocity
@ eMixed_CG_Discontinuous
std::shared_ptr< RiemannSolver > RiemannSolverSharedPtr
A shared pointer to an EquationSystem object.
AdvectionFactory & GetAdvectionFactory()
Gets the factory for initialising advection objects.
std::vector< std::pair< std::string, std::string > > SummaryList
void AddSummaryItem(SummaryList &l, const std::string &name, const std::string &value)
Adds a summary item to the summary info list.
RiemannSolverFactory & GetRiemannSolverFactory()
std::shared_ptr< GeomFactors > GeomFactorsSharedPtr
Pointer to a GeomFactors object.
@ eDeformed
Geometry is curved or has non-constant factors.
std::shared_ptr< MeshGraph > MeshGraphSharedPtr
std::vector< double > z(NPUPPER)
std::vector< double > d(NPUPPER *NPUPPER)
CFSBndCondFactory & GetCFSBndCondFactory()
Declaration of the boundary condition factory singleton.
static Array< OneD, Array< OneD, NekDouble > > NullNekDoubleArrayOfArray
ArtificialDiffusionFactory & GetArtificialDiffusionFactory()
Declaration of the artificial diffusion factory singleton.
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
void Svtvp(int n, const T alpha, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Svtvp (scalar times vector plus vector): z = alpha*x + y.
void Neg(int n, T *x, const int incx)
Negate x = -x.
T Vmin(int n, const T *x, const int incx)
Return the minimum element in x - called vmin to avoid conflict with min.
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
T Vsum(int n, const T *x, const int incx)
Subtract return sum(x)
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
void Zero(int n, T *x, const int incx)
Zero vector.
void FillWhiteNoise(int n, const T eps, T *x, const int incx, int outseed)
Fills a vector with white noise.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
scalarT< T > abs(scalarT< T > in)
scalarT< T > log(scalarT< T > in)
scalarT< T > sqrt(scalarT< T > in)