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:
Inheritance graph
[legend]
Collaboration diagram for Nektar::ExactSolverToro:
Collaboration graph
[legend]

Static Public Member Functions

static RiemannSolverSharedPtr create ()
 

Static Public Attributes

static std::string solverName
 

Protected Member Functions

 ExactSolverToro ()
 
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 ()
 
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)
 
- Protected Member Functions inherited from Nektar::SolverUtils::RiemannSolver
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
 
- 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 43 of file ExactSolverToro.h.

Constructor & Destructor Documentation

Nektar::ExactSolverToro::ExactSolverToro ( )
protected

Definition at line 49 of file ExactSolverToro.cpp.

Referenced by create().

50  {
51 
52  }

Member Function Documentation

static RiemannSolverSharedPtr Nektar::ExactSolverToro::create ( )
inlinestatic

Definition at line 46 of file ExactSolverToro.h.

References ExactSolverToro().

47  {
49  new ExactSolverToro());
50  }
boost::shared_ptr< RiemannSolver > RiemannSolverSharedPtr
A shared pointer to an EquationSystem object.
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 170 of file ExactSolverToro.cpp.

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

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

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

Definition at line 52 of file ExactSolverToro.h.