Nektar++
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)
 
virtual 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)
 Exact Riemann solver for the Euler equations. More...
 
- Protected Member Functions inherited from Nektar::CompressibleSolver
 CompressibleSolver (const LibUtilities::SessionReaderSharedPtr &pSession)
 
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)
 
virtual void v_ArraySolve (const Array< OneD, const Array< OneD, NekDouble > > &Fwd, const Array< OneD, const Array< OneD, NekDouble > > &Bwd, Array< OneD, Array< OneD, NekDouble > > &flux)
 
virtual void v_PointSolveVisc (NekDouble rhoL, NekDouble rhouL, NekDouble rhovL, NekDouble rhowL, NekDouble EL, NekDouble EpsL, NekDouble rhoR, NekDouble rhouR, NekDouble rhovR, NekDouble rhowR, NekDouble ER, NekDouble EpsR, NekDouble &rhof, NekDouble &rhouf, NekDouble &rhovf, NekDouble &rhowf, NekDouble &Ef, NekDouble &Epsf)
 
NekDouble GetRoeSoundSpeed (NekDouble rhoL, NekDouble pL, NekDouble eL, NekDouble HL, NekDouble srL, NekDouble rhoR, NekDouble pR, NekDouble eR, NekDouble HR, NekDouble srR, NekDouble HRoe, NekDouble URoe2, NekDouble srLR)
 
- Protected Member Functions inherited from Nektar::SolverUtils::RiemannSolver
SOLVER_UTILS_EXPORT RiemannSolver (const LibUtilities::SessionReaderSharedPtr &pSession)
 
virtual SOLVER_UTILS_EXPORT ~RiemannSolver ()
 
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...
 

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 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)
 
std::map< std::string, RSScalarFuncType > & GetScalars ()
 
std::map< std::string, RSVecFuncType > & GetVectors ()
 
std::map< std::string, RSParamFuncType > & GetParams ()
 
- Public Attributes inherited from Nektar::SolverUtils::RiemannSolver
int m_spacedim
 
- 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...
 

Detailed Description

Definition at line 42 of file ExactSolverToro.h.

Constructor & Destructor Documentation

◆ ExactSolverToro()

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

Definition at line 48 of file ExactSolverToro.cpp.

Referenced by create().

50  : CompressibleSolver(pSession)
51  {
52 
53  }
CompressibleSolver(const LibUtilities::SessionReaderSharedPtr &pSession)

Member Function Documentation

◆ create()

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

Definition at line 45 of file ExactSolverToro.h.

References ExactSolverToro().

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

◆ 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 
)
protectedvirtual

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 171 of file ExactSolverToro.cpp.

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

175  {
176  static NekDouble gamma = m_params["gamma"]();
177 
178  // Left and right variables.
179  NekDouble uL = rhouL / rhoL;
180  NekDouble vL = rhovL / rhoL;
181  NekDouble wL = rhowL / rhoL;
182  NekDouble uR = rhouR / rhoR;
183  NekDouble vR = rhovR / rhoR;
184  NekDouble wR = rhowR / rhoR;
185 
186  // Left and right pressure.
187  NekDouble pL = (gamma - 1.0) *
188  (EL - 0.5 * (rhouL*uL + rhovL*vL + rhowL*wL));
189  NekDouble pR = (gamma - 1.0) *
190  (ER - 0.5 * (rhouR*uR + rhovR*vR + rhowR*wR));
191 
192  // Compute gammas.
193  NekDouble g[] = { gamma,
194  (gamma-1.0)/(2.0*gamma),
195  (gamma+1.0)/(2.0*gamma),
196  2.0*gamma/(gamma-1.0),
197  2.0/(gamma-1.0),
198  2.0/(gamma+1.0),
199  (gamma-1.0)/(gamma+1.0),
200  0.5*(gamma-1.0),
201  gamma-1.0 };
202 
203  // Compute sound speeds.
204  NekDouble cL = sqrt(gamma * pL / rhoL);
205  NekDouble cR = sqrt(gamma * pR / rhoR);
206 
207  ASSERTL0(g[4]*(cL+cR) > (uR-uL), "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*(outP/(gamma-1.0) + 0.5*outRho*
391  (outU*outU + outV*outV + outW*outW) + outP);
392  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:216
#define TOL
double NekDouble
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...
#define NRITER
std::map< std::string, RSParamFuncType > m_params
Map of parameter function types.
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...

Member Data Documentation

◆ solverName

std::string Nektar::ExactSolverToro::solverName
static
Initial value:

Definition at line 52 of file ExactSolverToro.h.