44 "Upwind solver for the APE equation");
108 pF = pF + 0.5 * (c0L - u0L) *
109 (-rho0L * v0L * vL * (c0L * u0L + c0sqL) -
110 rho0L * w0L * wL * (c0L * u0L + c0sqL) +
111 (c0sqL - pow(u0L, 2)) * (rho0L * c0L * uL - pL)) /
112 (c0sqL - pow(u0L, 2));
116 (rho0L * c0L * (c0L * u0L + c0sqL) * (v0L * vL + w0L * wL) -
117 rho0L * c0sqL * uL * (c0sqL - pow(u0L, 2)) +
118 c0L * pL * (c0sqL - pow(u0L, 2))) /
119 (rho0L * c0sqL * (c0sqL - pow(u0L, 2)));
123 pF = pF + 0.5 * (c0R - u0R) *
124 (-rho0R * v0R * vR * (c0R * u0R + c0sqR) -
125 rho0R * w0R * wR * (c0R * u0R + c0sqR) +
126 (c0sqR - pow(u0R, 2)) * (rho0R * c0R * uR - pR)) /
127 (c0sqR - pow(u0R, 2));
131 (rho0R * c0R * (c0R * u0R + c0sqR) * (v0R * vR + w0R * wR) -
132 rho0R * c0sqR * uR * (c0sqR - pow(u0R, 2)) +
133 c0R * pR * (c0sqR - pow(u0R, 2))) /
134 (rho0R * c0sqR * (c0sqR - pow(u0R, 2)));
140 pF = pF + 0.5 * (c0L + u0L) *
141 (-rho0L * v0L * vL * (c0L * u0L - c0sqL) -
142 rho0L * w0L * wL * (c0L * u0L - c0sqL) +
143 (c0sqL - pow(u0L, 2)) * (rho0L * c0L * uL + pL)) /
144 (c0sqL - pow(u0L, 2));
148 (-rho0L * c0L * (c0L * u0L - c0sqL) * (v0L * vL + w0L * wL) +
149 rho0L * c0sqL * uL * (c0sqL - pow(u0L, 2)) +
150 c0L * pL * (c0sqL - pow(u0L, 2))) /
151 (rho0L * c0sqL * (c0sqL - pow(u0L, 2)));
155 pF = pF + 0.5 * (c0R + u0R) *
156 (-rho0R * v0R * vR * (c0R * u0R - c0sqR) -
157 rho0R * w0R * wR * (c0R * u0R - c0sqR) +
158 (c0sqR - pow(u0R, 2)) * (rho0R * c0R * uR + pR)) /
159 (c0sqR - pow(u0R, 2));
163 (-rho0R * c0R * (c0R * u0R - c0sqR) * (v0R * vR + w0R * wR) +
164 rho0R * c0sqR * uR * (c0sqR - pow(u0R, 2)) +
165 c0R * pR * (c0sqR - pow(u0R, 2))) /
166 (rho0R * c0sqR * (c0sqR - pow(u0R, 2)));
APEUpwindSolver(const LibUtilities::SessionReaderSharedPtr &pSession)
static std::string solverName
void v_PointSolve(NekDouble pL, NekDouble rhoL, NekDouble uL, NekDouble vL, NekDouble wL, NekDouble pR, NekDouble rhoR, NekDouble uR, NekDouble vR, NekDouble wR, NekDouble c0sqL, NekDouble rho0L, NekDouble u0L, NekDouble v0L, NekDouble w0L, NekDouble c0sqR, NekDouble rho0R, NekDouble u0R, NekDouble v0R, NekDouble w0R, NekDouble &pF, NekDouble &rhoF, NekDouble &uF, NekDouble &vF, NekDouble &wF) override
Upwind Riemann solver.
static RiemannSolverSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession)
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
std::shared_ptr< SessionReader > SessionReaderSharedPtr
RiemannSolverFactory & GetRiemannSolverFactory()
scalarT< T > sqrt(scalarT< T > in)