Nektar++
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Static Public Member Functions | Static Public Attributes | Protected Member Functions | List of all members
Nektar::ExactSolverToro Class Reference

#include <ExactSolverToro.h>

Inheritance diagram for Nektar::ExactSolverToro:
[legend]

Static Public Member Functions

static RiemannSolverSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession)
 

Static Public Attributes

static std::string solverName
 

Protected Member Functions

 ExactSolverToro (const LibUtilities::SessionReaderSharedPtr &pSession)
 
void v_PointSolve (NekDouble rhoL, NekDouble rhouL, NekDouble rhovL, NekDouble rhowL, NekDouble EL, NekDouble rhoR, NekDouble rhouR, NekDouble rhovR, NekDouble rhowR, NekDouble ER, NekDouble &rhof, NekDouble &rhouf, NekDouble &rhovf, NekDouble &rhowf, NekDouble &Ef) override
 Exact Riemann solver for the Euler equations. More...
 
- Protected Member Functions inherited from Nektar::CompressibleSolver
 CompressibleSolver (const LibUtilities::SessionReaderSharedPtr &pSession)
 Session ctor. More...
 
 CompressibleSolver ()
 Programmatic ctor. More...
 
void v_Solve (const int nDim, const Array< OneD, const Array< OneD, ND > > &Fwd, const Array< OneD, const Array< OneD, ND > > &Bwd, Array< OneD, Array< OneD, ND > > &flux) override
 
virtual void v_ArraySolve (const Array< OneD, const Array< OneD, ND > > &Fwd, const Array< OneD, const Array< OneD, ND > > &Bwd, Array< OneD, Array< OneD, ND > > &flux)
 
virtual void v_PointSolve (ND rhoL, ND rhouL, ND rhovL, ND rhowL, ND EL, ND rhoR, ND rhouR, ND rhovR, ND rhowR, ND ER, ND &rhof, ND &rhouf, ND &rhovf, ND &rhowf, ND &Ef)
 
ND GetRoeSoundSpeed (ND rhoL, ND pL, ND eL, ND HL, ND srL, ND rhoR, ND pR, ND eR, ND HR, ND srR, ND HRoe, ND URoe2, ND srLR)
 
- Protected Member Functions inherited from Nektar::SolverUtils::RiemannSolver
SOLVER_UTILS_EXPORT RiemannSolver ()
 
SOLVER_UTILS_EXPORT RiemannSolver (const LibUtilities::SessionReaderSharedPtr &pSession)
 
virtual SOLVER_UTILS_EXPORT ~RiemannSolver ()
 
virtual void v_Solve (const int nDim, const Array< OneD, const Array< OneD, NekDouble > > &Fwd, const Array< OneD, const Array< OneD, NekDouble > > &Bwd, Array< OneD, Array< OneD, NekDouble > > &flux)=0
 
SOLVER_UTILS_EXPORT void GenerateRotationMatrices (const Array< OneD, const Array< OneD, NekDouble > > &normals)
 Generate rotation matrices for 3D expansions. More...
 
void FromToRotation (Array< OneD, const NekDouble > &from, Array< OneD, const NekDouble > &to, NekDouble *mat)
 A function for creating a rotation matrix that rotates a vector from into another vector to. More...
 
SOLVER_UTILS_EXPORT void rotateToNormal (const Array< OneD, const Array< OneD, NekDouble > > &inarray, const Array< OneD, const Array< OneD, NekDouble > > &normals, const Array< OneD, const Array< OneD, NekDouble > > &vecLocs, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Rotate a vector field to trace normal. More...
 
SOLVER_UTILS_EXPORT void rotateFromNormal (const Array< OneD, const Array< OneD, NekDouble > > &inarray, const Array< OneD, const Array< OneD, NekDouble > > &normals, const Array< OneD, const Array< OneD, NekDouble > > &vecLocs, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Rotate a vector field from trace normal. More...
 
SOLVER_UTILS_EXPORT bool CheckScalars (std::string name)
 Determine whether a scalar has been defined in m_scalars. More...
 
SOLVER_UTILS_EXPORT bool CheckVectors (std::string name)
 Determine whether a vector has been defined in m_vectors. More...
 
SOLVER_UTILS_EXPORT bool CheckParams (std::string name)
 Determine whether a parameter has been defined in m_params. More...
 
SOLVER_UTILS_EXPORT bool CheckAuxScal (std::string name)
 Determine whether a scalar has been defined in m_auxScal. More...
 
SOLVER_UTILS_EXPORT bool CheckAuxVec (std::string name)
 Determine whether a vector has been defined in m_auxVec. More...
 
virtual SOLVER_UTILS_EXPORT void v_CalcFluxJacobian (const int nDim, const Array< OneD, const Array< OneD, NekDouble > > &Fwd, const Array< OneD, const Array< OneD, NekDouble > > &Bwd, const Array< OneD, const Array< OneD, NekDouble > > &normals, DNekBlkMatSharedPtr &FJac, DNekBlkMatSharedPtr &BJac)
 

Additional Inherited Members

- Public Member Functions inherited from Nektar::SolverUtils::RiemannSolver
SOLVER_UTILS_EXPORT void Solve (const int nDim, const Array< OneD, const Array< OneD, NekDouble > > &Fwd, const Array< OneD, const Array< OneD, NekDouble > > &Bwd, Array< OneD, Array< OneD, NekDouble > > &flux)
 Perform the Riemann solve given the forwards and backwards spaces. More...
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetScalar (std::string name, FuncPointerT func, ObjectPointerT obj)
 
void SetScalar (std::string name, RSScalarFuncType fp)
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetVector (std::string name, FuncPointerT func, ObjectPointerT obj)
 
void SetVector (std::string name, RSVecFuncType fp)
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetParam (std::string name, FuncPointerT func, ObjectPointerT obj)
 
void SetALEFlag (bool &ALE)
 
void SetParam (std::string name, RSParamFuncType fp)
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetAuxScal (std::string name, FuncPointerT func, ObjectPointerT obj)
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetAuxVec (std::string name, FuncPointerT func, ObjectPointerT obj)
 
void SetAuxVec (std::string name, RSVecFuncType fp)
 
std::map< std::string, RSScalarFuncType > & GetScalars ()
 
std::map< std::string, RSVecFuncType > & GetVectors ()
 
std::map< std::string, RSParamFuncType > & GetParams ()
 
SOLVER_UTILS_EXPORT void CalcFluxJacobian (const int nDim, const Array< OneD, const Array< OneD, NekDouble > > &Fwd, const Array< OneD, const Array< OneD, NekDouble > > &Bwd, DNekBlkMatSharedPtr &FJac, DNekBlkMatSharedPtr &BJac)
 Calculate the flux jacobian of Fwd and Bwd. More...
 
- Public Attributes inherited from Nektar::SolverUtils::RiemannSolver
int m_spacedim
 
- Protected Types inherited from Nektar::CompressibleSolver
using ND = NekDouble
 
- Protected Attributes inherited from Nektar::CompressibleSolver
bool m_pointSolve
 
EquationOfStateSharedPtr m_eos
 
bool m_idealGas
 
- Protected Attributes inherited from Nektar::SolverUtils::RiemannSolver
bool m_requiresRotation
 Indicates whether the Riemann solver requires a rotation to be applied to the velocity fields. More...
 
std::map< std::string, RSScalarFuncTypem_scalars
 Map of scalar function types. More...
 
std::map< std::string, RSVecFuncTypem_vectors
 Map of vector function types. More...
 
std::map< std::string, RSParamFuncTypem_params
 Map of parameter function types. More...
 
std::map< std::string, RSScalarFuncTypem_auxScal
 Map of auxiliary scalar function types. More...
 
std::map< std::string, RSVecFuncTypem_auxVec
 Map of auxiliary vector function types. More...
 
Array< OneD, Array< OneD, NekDouble > > m_rotMat
 Rotation matrices for each trace quadrature point. More...
 
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_rotStorage
 Rotation storage. More...
 
bool m_ALESolver = false
 Flag if using the ALE formulation. More...
 

Detailed Description

Definition at line 43 of file ExactSolverToro.h.

Constructor & Destructor Documentation

◆ ExactSolverToro()

Nektar::ExactSolverToro::ExactSolverToro ( const LibUtilities::SessionReaderSharedPtr pSession)
protected

Definition at line 47 of file ExactSolverToro.cpp.

49 : CompressibleSolver(pSession)
50{
51}
CompressibleSolver()
Programmatic ctor.

Referenced by create().

Member Function Documentation

◆ create()

static RiemannSolverSharedPtr Nektar::ExactSolverToro::create ( const LibUtilities::SessionReaderSharedPtr pSession)
inlinestatic

Definition at line 46 of file ExactSolverToro.h.

48 {
49 return RiemannSolverSharedPtr(new ExactSolverToro(pSession));
50 }
ExactSolverToro(const LibUtilities::SessionReaderSharedPtr &pSession)
std::shared_ptr< RiemannSolver > RiemannSolverSharedPtr
A shared pointer to an EquationSystem object.

References ExactSolverToro().

◆ v_PointSolve()

void Nektar::ExactSolverToro::v_PointSolve ( NekDouble  rhoL,
NekDouble  rhouL,
NekDouble  rhovL,
NekDouble  rhowL,
NekDouble  EL,
NekDouble  rhoR,
NekDouble  rhouR,
NekDouble  rhovR,
NekDouble  rhowR,
NekDouble  ER,
NekDouble rhof,
NekDouble rhouf,
NekDouble rhovf,
NekDouble rhowf,
NekDouble Ef 
)
overrideprotectedvirtual

Exact Riemann solver for the Euler equations.

This algorithm is transcribed from:

"Riemann Solvers and Numerical Methods for Fluid Dynamics: A Practical Introduction", E. F. Toro (3rd edition, 2009).

The full Fortran 77 routine can be found at the end of chapter 4 (section 4.9). This transcription is essentially the functions STARPU and SAMPLE glued together, and variable names are kept mostly the same. See the preceding chapter which explains the derivation of the solver. The routines PREFUN and GUESSP are kept separate and are reproduced above.

Parameters
rhoLDensity left state.
rhoRDensity right state.
rhouLx-momentum component left state.
rhouRx-momentum component right state.
rhovLy-momentum component left state.
rhovRy-momentum component right state.
rhowLz-momentum component left state.
rhowRz-momentum component right state.
ELEnergy left state.
EREnergy right state.
rhofComputed Riemann flux for density.
rhoufComputed Riemann flux for x-momentum component
rhovfComputed Riemann flux for y-momentum component
rhowfComputed Riemann flux for z-momentum component
EfComputed Riemann flux for energy.

Reimplemented from Nektar::CompressibleSolver.

Definition at line 169 of file ExactSolverToro.cpp.

174{
175 static NekDouble gamma = m_params["gamma"]();
176
177 // Left and right variables.
178 NekDouble uL = rhouL / rhoL;
179 NekDouble vL = rhovL / rhoL;
180 NekDouble wL = rhowL / rhoL;
181 NekDouble uR = rhouR / rhoR;
182 NekDouble vR = rhovR / rhoR;
183 NekDouble wR = rhowR / rhoR;
184
185 // Left and right pressure.
186 NekDouble pL =
187 (gamma - 1.0) * (EL - 0.5 * (rhouL * uL + rhovL * vL + rhowL * wL));
188 NekDouble pR =
189 (gamma - 1.0) * (ER - 0.5 * (rhouR * uR + rhovR * vR + rhowR * wR));
190
191 // Compute gammas.
192 NekDouble g[] = {gamma,
193 (gamma - 1.0) / (2.0 * gamma),
194 (gamma + 1.0) / (2.0 * gamma),
195 2.0 * gamma / (gamma - 1.0),
196 2.0 / (gamma - 1.0),
197 2.0 / (gamma + 1.0),
198 (gamma - 1.0) / (gamma + 1.0),
199 0.5 * (gamma - 1.0),
200 gamma - 1.0};
201
202 // Compute sound speeds.
203 NekDouble cL = sqrt(gamma * pL / rhoL);
204 NekDouble cR = sqrt(gamma * pR / rhoR);
205
206 ASSERTL0(g[4] * (cL + cR) > (uR - uL),
207 "Vacuum is generated by given data.");
208
209 // Guess initial pressure.
210 NekDouble pOld = guessp(g, rhoL, uL, pL, cL, rhoR, uR, pR, cR);
211 NekDouble uDiff = uR - uL;
212 NekDouble p, fL, fR, fLd, fRd, change;
213 int k;
214
215 // Newton-Raphson iteration for pressure in star region.
216 for (k = 0; k < NRITER; ++k)
217 {
218 prefun(g, pOld, rhoL, pL, cL, fL, fLd);
219 prefun(g, pOld, rhoR, pR, cR, fR, fRd);
220 p = pOld - (fL + fR + uDiff) / (fLd + fRd);
221 change = 2 * fabs((p - pOld) / (p + pOld));
222
223 if (change <= TOL)
224 {
225 break;
226 }
227
228 if (p < 0.0)
229 {
230 p = TOL;
231 }
232
233 pOld = p;
234 }
235
236 ASSERTL0(k < NRITER, "Divergence in Newton-Raphson scheme");
237
238 // Compute velocity in star region.
239 NekDouble u = 0.5 * (uL + uR + fR - fL);
240
241 // -- SAMPLE ROUTINE --
242 // The variable S aligns with the Fortran parameter S of the SAMPLE
243 // routine, but is hard-coded as 0 (and should be optimised out by the
244 // compiler). Since we are using a Godunov scheme we pick this as 0 (see
245 // chapter 6 of Toro 2009).
246 const NekDouble S = 0.0;
247
248 // Computed primitive variables.
249 NekDouble outRho, outU, outV, outW, outP;
250
251 if (S <= u)
252 {
253 if (p <= pL)
254 {
255 // Left rarefaction
256 NekDouble shL = uL - cL;
257 if (S <= shL)
258 {
259 // Sampled point is left data state.
260 outRho = rhoL;
261 outU = uL;
262 outV = vL;
263 outW = wL;
264 outP = pL;
265 }
266 else
267 {
268 NekDouble cmL = cL * pow(p / pL, g[1]);
269 NekDouble stL = u - cmL;
270
271 if (S > stL)
272 {
273 // Sampled point is star left state
274 outRho = rhoL * pow(p / pL, 1.0 / gamma);
275 outU = u;
276 outV = vL;
277 outW = wL;
278 outP = p;
279 }
280 else
281 {
282 // Sampled point is inside left fan
283 NekDouble c = g[5] * (cL + g[7] * (uL - S));
284 outRho = rhoL * pow(c / cL, g[4]);
285 outU = g[5] * (cL + g[7] * uL + S);
286 outV = vL;
287 outW = wL;
288 outP = pL * pow(c / cL, g[3]);
289 }
290 }
291 }
292 else
293 {
294 // Left shock
295 NekDouble pmL = p / pL;
296 NekDouble SL = uL - cL * sqrt(g[2] * pmL + g[1]);
297 if (S <= SL)
298 {
299 // Sampled point is left data state
300 outRho = rhoL;
301 outU = uL;
302 outV = vL;
303 outW = wL;
304 outP = pL;
305 }
306 else
307 {
308 // Sampled point is star left state
309 outRho = rhoL * (pmL + g[6]) / (pmL * g[6] + 1.0);
310 outU = u;
311 outV = vL;
312 outW = wL;
313 outP = p;
314 }
315 }
316 }
317 else
318 {
319 if (p > pR)
320 {
321 // Right shock
322 NekDouble pmR = p / pR;
323 NekDouble SR = uR + cR * sqrt(g[2] * pmR + g[1]);
324 if (S >= SR)
325 {
326 // Sampled point is right data state
327 outRho = rhoR;
328 outU = uR;
329 outV = vR;
330 outW = wR;
331 outP = pR;
332 }
333 else
334 {
335 // Sampled point is star right state
336 outRho = rhoR * (pmR + g[6]) / (pmR * g[6] + 1.0);
337 outU = u;
338 outV = vR;
339 outW = wR;
340 outP = p;
341 }
342 }
343 else
344 {
345 // Right rarefaction
346 NekDouble shR = uR + cR;
347
348 if (S >= shR)
349 {
350 // Sampled point is right data state
351 outRho = rhoR;
352 outU = uR;
353 outV = vR;
354 outW = wR;
355 outP = pR;
356 }
357 else
358 {
359 NekDouble cmR = cR * pow(p / pR, g[1]);
360 NekDouble stR = u + cmR;
361
362 if (S <= stR)
363 {
364 // Sampled point is star right state
365 outRho = rhoR * pow(p / pR, 1.0 / gamma);
366 outU = u;
367 outV = vR;
368 outW = wR;
369 outP = p;
370 }
371 else
372 {
373 // Sampled point is inside left fan
374 NekDouble c = g[5] * (cR - g[7] * (uR - S));
375 outRho = rhoR * pow(c / cR, g[4]);
376 outU = g[5] * (-cR + g[7] * uR + S);
377 outV = vR;
378 outW = wR;
379 outP = pR * pow(c / cR, g[3]);
380 }
381 }
382 }
383 }
384
385 // Transform computed primitive variables to fluxes.
386 rhof = outRho * outU;
387 rhouf = outP + outRho * outU * outU;
388 rhovf = outRho * outU * outV;
389 rhowf = outRho * outU * outW;
390 Ef = outU *
391 (outP / (gamma - 1.0) +
392 0.5 * outRho * (outU * outU + outV * outV + outW * outW) + outP);
393}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
#define TOL
#define NRITER
std::map< std::string, RSParamFuncType > m_params
Map of parameter function types.
NekDouble guessp(NekDouble g[], NekDouble rhoL, NekDouble uL, NekDouble pL, NekDouble cL, NekDouble rhoR, NekDouble uR, NekDouble pR, NekDouble cR)
Use either PVRS, two-rarefaction or two-shock Riemann solvers to calculate an initial pressure for th...
void prefun(NekDouble *g, NekDouble p, NekDouble dk, NekDouble pk, NekDouble ck, NekDouble &f, NekDouble &fd)
Evaluate pressure functions fL and fR in Newton iteration of Riemann solver (see equation 4....
double NekDouble
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:285

References ASSERTL0, Nektar::guessp(), Nektar::SolverUtils::RiemannSolver::m_params, NRITER, CellMLToNektar.cellml_metadata::p, Nektar::prefun(), tinysimd::sqrt(), and TOL.

Member Data Documentation

◆ solverName

std::string Nektar::ExactSolverToro::solverName
static
Initial value:
=
"ExactToro", ExactSolverToro::create, "Exact Riemann solver")
static RiemannSolverSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession)
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
RiemannSolverFactory & GetRiemannSolverFactory()

Definition at line 52 of file ExactSolverToro.h.