Nektar++
Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | List of all members
Nektar::SolverUtils::RiemannSolver Class Referenceabstract

The RiemannSolver class provides an abstract interface under which solvers for various Riemann problems can be implemented. More...

#include <RiemannSolver.h>

Inheritance diagram for Nektar::SolverUtils::RiemannSolver:
[legend]

Public Member Functions

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

int m_spacedim
 

Protected Member Functions

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)
 

Protected Attributes

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

The RiemannSolver class provides an abstract interface under which solvers for various Riemann problems can be implemented.

Definition at line 61 of file RiemannSolver.h.

Constructor & Destructor Documentation

◆ RiemannSolver() [1/2]

Nektar::SolverUtils::RiemannSolver::RiemannSolver ( )
protected

Definition at line 76 of file RiemannSolver.cpp.

77  : m_requiresRotation(false), m_rotStorage(3) {}
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_rotStorage
Rotation storage.
bool m_requiresRotation
Indicates whether the Riemann solver requires a rotation to be applied to the velocity fields.

◆ RiemannSolver() [2/2]

Nektar::SolverUtils::RiemannSolver::RiemannSolver ( const LibUtilities::SessionReaderSharedPtr pSession)
protected

Definition at line 79 of file RiemannSolver.cpp.

81  : m_requiresRotation(false), m_rotStorage(3)
82  {
83  boost::ignore_unused(pSession);
84  }

◆ ~RiemannSolver()

virtual SOLVER_UTILS_EXPORT Nektar::SolverUtils::RiemannSolver::~RiemannSolver ( )
inlineprotectedvirtual

Definition at line 177 of file RiemannSolver.h.

178  {};

Member Function Documentation

◆ CalcFluxJacobian()

void Nektar::SolverUtils::RiemannSolver::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.

Parameters
FwdForwards trace space.
BwdBackwards trace space.
fluxResultant flux along trace space.

Definition at line 525 of file RiemannSolver.cpp.

531  {
532  int nPts = Fwd[0].size();
533 
534  if (m_requiresRotation)
535  {
536  ASSERTL1(CheckVectors("N"), "N not defined.");
537  ASSERTL1(CheckAuxVec("vecLocs"), "vecLocs not defined.");
538  const Array<OneD, const Array<OneD, NekDouble> > normals =
539  m_vectors["N"]();
540  const Array<OneD, const Array<OneD, NekDouble> > vecLocs =
541  m_auxVec["vecLocs"]();
542 
543  v_CalcFluxJacobian(nDim, Fwd, Bwd,normals, FJac, BJac);
544  }
545  else
546  {
547  Array<OneD, Array<OneD, NekDouble> > normals(nDim);
548  for(int i=0;i< nDim;i++)
549  {
550  normals[i] = Array<OneD, NekDouble> (nPts,0.0);
551  }
552  Vmath::Fill(nPts, 1.0, normals[0],1);
553 
554  v_CalcFluxJacobian(nDim, Fwd, Bwd,normals, FJac, BJac);
555  }
556  }
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:250
SOLVER_UTILS_EXPORT bool CheckAuxVec(std::string name)
Determine whether a vector has been defined in m_auxVec.
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)
std::map< std::string, RSVecFuncType > m_vectors
Map of vector function types.
SOLVER_UTILS_EXPORT bool CheckVectors(std::string name)
Determine whether a vector has been defined in m_vectors.
std::map< std::string, RSVecFuncType > m_auxVec
Map of auxiliary vector function types.
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
Definition: Vmath.cpp:45

References ASSERTL1, CheckAuxVec(), CheckVectors(), Vmath::Fill(), m_auxVec, m_requiresRotation, m_vectors, and v_CalcFluxJacobian().

◆ CheckAuxScal()

bool Nektar::SolverUtils::RiemannSolver::CheckAuxScal ( std::string  name)
protected

Determine whether a scalar has been defined in m_auxScal.

Parameters
nameScalar name.

Definition at line 369 of file RiemannSolver.cpp.

370  {
371  return m_auxScal.find(name) != m_auxScal.end();
372  }
std::map< std::string, RSScalarFuncType > m_auxScal
Map of auxiliary scalar function types.

References m_auxScal, and CellMLToNektar.pycml::name.

◆ CheckAuxVec()

bool Nektar::SolverUtils::RiemannSolver::CheckAuxVec ( std::string  name)
protected

Determine whether a vector has been defined in m_auxVec.

Parameters
nameVector name.

Definition at line 379 of file RiemannSolver.cpp.

380  {
381  return m_auxVec.find(name) != m_auxVec.end();
382  }

References m_auxVec, and CellMLToNektar.pycml::name.

Referenced by CalcFluxJacobian(), and Solve().

◆ CheckParams()

bool Nektar::SolverUtils::RiemannSolver::CheckParams ( std::string  name)
protected

Determine whether a parameter has been defined in m_params.

Parameters
nameParameter name.

Definition at line 359 of file RiemannSolver.cpp.

360  {
361  return m_params.find(name) != m_params.end();
362  }
std::map< std::string, RSParamFuncType > m_params
Map of parameter function types.

References m_params, and CellMLToNektar.pycml::name.

Referenced by Nektar::UpwindPulseSolver::RiemannSolverUpwind().

◆ CheckScalars()

bool Nektar::SolverUtils::RiemannSolver::CheckScalars ( std::string  name)
protected

Determine whether a scalar has been defined in m_scalars.

Parameters
nameScalar name.

Definition at line 339 of file RiemannSolver.cpp.

340  {
341  return m_scalars.find(name) != m_scalars.end();
342  }
std::map< std::string, RSScalarFuncType > m_scalars
Map of scalar function types.

References m_scalars, and CellMLToNektar.pycml::name.

Referenced by Nektar::SolverUtils::UpwindLDGSolver::v_Solve(), Nektar::SolverUtils::UpwindSolver::v_Solve(), and Nektar::UpwindPulseSolver::v_Solve().

◆ CheckVectors()

bool Nektar::SolverUtils::RiemannSolver::CheckVectors ( std::string  name)
protected

Determine whether a vector has been defined in m_vectors.

Parameters
nameVector name.

Definition at line 349 of file RiemannSolver.cpp.

350  {
351  return m_vectors.find(name) != m_vectors.end();
352  }

References m_vectors, and CellMLToNektar.pycml::name.

Referenced by CalcFluxJacobian(), Nektar::AcousticSolver::GetRotBasefield(), Solve(), and Nektar::RoeSolverSIMD::v_Solve().

◆ FromToRotation()

void Nektar::SolverUtils::RiemannSolver::FromToRotation ( Array< OneD, const NekDouble > &  from,
Array< OneD, const NekDouble > &  to,
NekDouble mat 
)
protected

A function for creating a rotation matrix that rotates a vector from into another vector to.

Authors: Tomas Möller, John Hughes "Efficiently Building a Matrix to Rotate One Vector to Another" Journal of Graphics Tools, 4(4):1-4, 1999

Parameters
fromNormalised 3-vector to rotate from.
toNormalised 3-vector to rotate to.
outResulting 3x3 rotation matrix (row-major order).

Definition at line 432 of file RiemannSolver.cpp.

436  {
437  NekDouble v[3];
438  NekDouble e, h, f;
439 
440  CROSS(v, from, to);
441  e = DOT(from, to);
442  f = (e < 0)? -e:e;
443  if (f > 1.0 - EPSILON)
444  {
445  NekDouble u[3], v[3];
446  NekDouble x[3];
447  NekDouble c1, c2, c3;
448  int i, j;
449 
450  x[0] = (from[0] > 0.0)? from[0] : -from[0];
451  x[1] = (from[1] > 0.0)? from[1] : -from[1];
452  x[2] = (from[2] > 0.0)? from[2] : -from[2];
453 
454  if (x[0] < x[1])
455  {
456  if (x[0] < x[2])
457  {
458  x[0] = 1.0; x[1] = x[2] = 0.0;
459  }
460  else
461  {
462  x[2] = 1.0; x[0] = x[1] = 0.0;
463  }
464  }
465  else
466  {
467  if (x[1] < x[2])
468  {
469  x[1] = 1.0; x[0] = x[2] = 0.0;
470  }
471  else
472  {
473  x[2] = 1.0; x[0] = x[1] = 0.0;
474  }
475  }
476 
477  u[0] = x[0] - from[0];
478  u[1] = x[1] - from[1];
479  u[2] = x[2] - from[2];
480  v[0] = x[0] - to [0];
481  v[1] = x[1] - to [1];
482  v[2] = x[2] - to [2];
483 
484  c1 = 2.0 / DOT(u, u);
485  c2 = 2.0 / DOT(v, v);
486  c3 = c1 * c2 * DOT(u, v);
487 
488  for (i = 0; i < 3; i++) {
489  for (j = 0; j < 3; j++) {
490  mat[3*i+j] = - c1 * u[i] * u[j]
491  - c2 * v[i] * v[j]
492  + c3 * v[i] * u[j];
493  }
494  mat[i+3*i] += 1.0;
495  }
496  }
497  else
498  {
499  NekDouble hvx, hvz, hvxy, hvxz, hvyz;
500  h = 1.0/(1.0 + e);
501  hvx = h * v[0];
502  hvz = h * v[2];
503  hvxy = hvx * v[1];
504  hvxz = hvx * v[2];
505  hvyz = hvz * v[1];
506  mat[0] = e + hvx * v[0];
507  mat[1] = hvxy - v[2];
508  mat[2] = hvxz + v[1];
509  mat[3] = hvxy + v[2];
510  mat[4] = e + h * v[1] * v[1];
511  mat[5] = hvyz - v[0];
512  mat[6] = hvxz - v[1];
513  mat[7] = hvyz + v[0];
514  mat[8] = e + hvz * v[2];
515  }
516  }
#define EPSILON
#define CROSS(dest, v1, v2)
#define DOT(v1, v2)
double NekDouble

References CROSS, DOT, and EPSILON.

Referenced by GenerateRotationMatrices().

◆ GenerateRotationMatrices()

void Nektar::SolverUtils::RiemannSolver::GenerateRotationMatrices ( const Array< OneD, const Array< OneD, NekDouble > > &  normals)
protected

Generate rotation matrices for 3D expansions.

Definition at line 387 of file RiemannSolver.cpp.

389  {
390  Array<OneD, NekDouble> xdir(3,0.0);
391  Array<OneD, NekDouble> tn (3);
392  NekDouble tmp[9];
393  const int nq = normals[0].size();
394  int i, j;
395  xdir[0] = 1.0;
396 
397  // Allocate storage for rotation matrices.
398  m_rotMat = Array<OneD, Array<OneD, NekDouble> >(9);
399 
400  for (i = 0; i < 9; ++i)
401  {
402  m_rotMat[i] = Array<OneD, NekDouble>(nq);
403  }
404  for (i = 0; i < normals[0].size(); ++i)
405  {
406  // Generate matrix which takes us from (1,0,0) vector to trace
407  // normal.
408  tn[0] = normals[0][i];
409  tn[1] = normals[1][i];
410  tn[2] = normals[2][i];
411  FromToRotation(tn, xdir, tmp);
412 
413  for (j = 0; j < 9; ++j)
414  {
415  m_rotMat[j][i] = tmp[j];
416  }
417  }
418  }
Array< OneD, Array< OneD, NekDouble > > m_rotMat
Rotation matrices for each trace quadrature point.
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.

References FromToRotation(), and m_rotMat.

Referenced by rotateToNormal(), and Nektar::RoeSolverSIMD::v_Solve().

◆ GetParams()

std::map<std::string, RSParamFuncType>& Nektar::SolverUtils::RiemannSolver::GetParams ( )
inline

Definition at line 140 of file RiemannSolver.h.

141  {
142  return m_params;
143  }

References m_params.

◆ GetScalars()

std::map<std::string, RSScalarFuncType>& Nektar::SolverUtils::RiemannSolver::GetScalars ( )
inline

Definition at line 130 of file RiemannSolver.h.

131  {
132  return m_scalars;
133  }

References m_scalars.

◆ GetVectors()

std::map<std::string, RSVecFuncType>& Nektar::SolverUtils::RiemannSolver::GetVectors ( )
inline

Definition at line 135 of file RiemannSolver.h.

136  {
137  return m_vectors;
138  }

References m_vectors.

◆ rotateFromNormal()

void Nektar::SolverUtils::RiemannSolver::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 
)
protected

Rotate a vector field from trace normal.

This function performs a rotation of the triad of vector components provided in inarray so that the first component aligns with the Cartesian components; it performs the inverse operation of RiemannSolver::rotateToNormal.

Definition at line 258 of file RiemannSolver.cpp.

263  {
264  for (int i = 0; i < inarray.size(); ++i)
265  {
266  Vmath::Vcopy(inarray[i].size(), inarray[i], 1,
267  outarray[i], 1);
268  }
269 
270  for (int i = 0; i < vecLocs.size(); i++)
271  {
272  ASSERTL1(vecLocs[i].size() == normals.size(),
273  "vecLocs[i] element count mismatch");
274 
275  switch (normals.size())
276  {
277  case 1:
278  { // do nothing
279  const int nq = normals[0].size();
280  const int vx = (int)vecLocs[i][0];
281  Vmath::Vmul (nq, inarray [vx], 1, normals [0], 1,
282  outarray[vx], 1);
283  break;
284  }
285  case 2:
286  {
287  const int nq = normals[0].size();
288  const int vx = (int)vecLocs[i][0];
289  const int vy = (int)vecLocs[i][1];
290 
291  Vmath::Vmul (nq, inarray [vy], 1, normals [1], 1,
292  outarray[vx], 1);
293  Vmath::Vvtvm(nq, inarray [vx], 1, normals [0], 1,
294  outarray[vx], 1, outarray[vx], 1);
295  Vmath::Vmul (nq, inarray [vx], 1, normals [1], 1,
296  outarray[vy], 1);
297  Vmath::Vvtvp(nq, inarray [vy], 1, normals [0], 1,
298  outarray[vy], 1, outarray[vy], 1);
299  break;
300  }
301 
302  case 3:
303  {
304  const int nq = normals[0].size();
305  const int vx = (int)vecLocs[i][0];
306  const int vy = (int)vecLocs[i][1];
307  const int vz = (int)vecLocs[i][2];
308 
309  Vmath::Vvtvvtp(nq, inarray [vx], 1, m_rotMat[0], 1,
310  inarray [vy], 1, m_rotMat[3], 1,
311  outarray[vx], 1);
312  Vmath::Vvtvp (nq, inarray [vz], 1, m_rotMat[6], 1,
313  outarray[vx], 1, outarray[vx], 1);
314  Vmath::Vvtvvtp(nq, inarray [vx], 1, m_rotMat[1], 1,
315  inarray [vy], 1, m_rotMat[4], 1,
316  outarray[vy], 1);
317  Vmath::Vvtvp (nq, inarray [vz], 1, m_rotMat[7], 1,
318  outarray[vy], 1, outarray[vy], 1);
319  Vmath::Vvtvvtp(nq, inarray [vx], 1, m_rotMat[2], 1,
320  inarray [vy], 1, m_rotMat[5], 1,
321  outarray[vz], 1);
322  Vmath::Vvtvp (nq, inarray [vz], 1, m_rotMat[8], 1,
323  outarray[vz], 1, outarray[vz], 1);
324  break;
325  }
326 
327  default:
328  ASSERTL1(false, "Invalid space dimension.");
329  break;
330  }
331  }
332  }
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.
Definition: Vmath.cpp:192
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
Definition: Vmath.cpp:513
void Vvtvm(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)
vvtvm (vector times vector plus vector): z = w*x - y
Definition: Vmath.cpp:541
void Vvtvvtp(int n, const T *v, int incv, const T *w, int incw, const T *x, int incx, const T *y, int incy, T *z, int incz)
vvtvvtp (vector times vector plus vector times vector):
Definition: Vmath.cpp:625
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1199

References ASSERTL1, m_rotMat, Vmath::Vcopy(), Vmath::Vmul(), Vmath::Vvtvm(), Vmath::Vvtvp(), and Vmath::Vvtvvtp().

Referenced by Solve().

◆ rotateToNormal()

void Nektar::SolverUtils::RiemannSolver::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 
)
protected

Rotate a vector field to trace normal.

This function performs a rotation of a vector so that the first component aligns with the trace normal direction.

The vectors components are stored in inarray. Their locations must be specified in the "vecLocs" array. vecLocs[0] contains the locations of the first vectors components, vecLocs[1] those of the second and so on.

In 2D, this is accomplished through the transform:

\[ (u_x, u_y) = (n_x u_x + n_y u_y, -n_x v_x + n_y v_y) \]

In 3D, we generate a (non-unique) transformation using RiemannSolver::fromToRotation.

Definition at line 167 of file RiemannSolver.cpp.

172  {
173  for (int i = 0; i < inarray.size(); ++i)
174  {
175  Vmath::Vcopy(inarray[i].size(), inarray[i], 1,
176  outarray[i], 1);
177  }
178 
179  for (int i = 0; i < vecLocs.size(); i++)
180  {
181  ASSERTL1(vecLocs[i].size() == normals.size(),
182  "vecLocs[i] element count mismatch");
183 
184  switch (normals.size())
185  {
186  case 1:
187  { // do nothing
188  const int nq = inarray[0].size();
189  const int vx = (int)vecLocs[i][0];
190  Vmath::Vmul (nq, inarray [vx], 1, normals [0], 1,
191  outarray[vx], 1);
192  break;
193  }
194  case 2:
195  {
196  const int nq = inarray[0].size();
197  const int vx = (int)vecLocs[i][0];
198  const int vy = (int)vecLocs[i][1];
199 
200  Vmath::Vmul (nq, inarray [vx], 1, normals [0], 1,
201  outarray[vx], 1);
202  Vmath::Vvtvp(nq, inarray [vy], 1, normals [1], 1,
203  outarray[vx], 1, outarray[vx], 1);
204  Vmath::Vmul (nq, inarray [vx], 1, normals [1], 1,
205  outarray[vy], 1);
206  Vmath::Vvtvm(nq, inarray [vy], 1, normals [0], 1,
207  outarray[vy], 1, outarray[vy], 1);
208  break;
209  }
210 
211  case 3:
212  {
213  const int nq = inarray[0].size();
214  const int vx = (int)vecLocs[i][0];
215  const int vy = (int)vecLocs[i][1];
216  const int vz = (int)vecLocs[i][2];
217 
218  // Generate matrices if they don't already exist.
219  if (m_rotMat.size() == 0)
220  {
221  GenerateRotationMatrices(normals);
222  }
223 
224  // Apply rotation matrices.
225  Vmath::Vvtvvtp(nq, inarray [vx], 1, m_rotMat[0], 1,
226  inarray [vy], 1, m_rotMat[1], 1,
227  outarray[vx], 1);
228  Vmath::Vvtvp (nq, inarray [vz], 1, m_rotMat[2], 1,
229  outarray[vx], 1, outarray[vx], 1);
230  Vmath::Vvtvvtp(nq, inarray [vx], 1, m_rotMat[3], 1,
231  inarray [vy], 1, m_rotMat[4], 1,
232  outarray[vy], 1);
233  Vmath::Vvtvp (nq, inarray [vz], 1, m_rotMat[5], 1,
234  outarray[vy], 1, outarray[vy], 1);
235  Vmath::Vvtvvtp(nq, inarray [vx], 1, m_rotMat[6], 1,
236  inarray [vy], 1, m_rotMat[7], 1,
237  outarray[vz], 1);
238  Vmath::Vvtvp (nq, inarray [vz], 1, m_rotMat[8], 1,
239  outarray[vz], 1, outarray[vz], 1);
240  break;
241  }
242 
243  default:
244  ASSERTL1(false, "Invalid space dimension.");
245  break;
246  }
247  }
248  }
SOLVER_UTILS_EXPORT void GenerateRotationMatrices(const Array< OneD, const Array< OneD, NekDouble > > &normals)
Generate rotation matrices for 3D expansions.

References ASSERTL1, GenerateRotationMatrices(), m_rotMat, Vmath::Vcopy(), Vmath::Vmul(), Vmath::Vvtvm(), Vmath::Vvtvp(), and Vmath::Vvtvvtp().

Referenced by Nektar::AcousticSolver::GetRotBasefield(), and Solve().

◆ SetAuxScal()

template<typename FuncPointerT , typename ObjectPointerT >
void Nektar::SolverUtils::RiemannSolver::SetAuxScal ( std::string  name,
FuncPointerT  func,
ObjectPointerT  obj 
)
inline

Definition at line 110 of file RiemannSolver.h.

113  {
114  m_auxScal[name] = std::bind(func, obj);
115  }

References m_auxScal, and CellMLToNektar.pycml::name.

◆ SetAuxVec() [1/2]

template<typename FuncPointerT , typename ObjectPointerT >
void Nektar::SolverUtils::RiemannSolver::SetAuxVec ( std::string  name,
FuncPointerT  func,
ObjectPointerT  obj 
)
inline

Definition at line 118 of file RiemannSolver.h.

121  {
122  m_auxVec[name] = std::bind(func, obj);
123  }

References m_auxVec, and CellMLToNektar.pycml::name.

◆ SetAuxVec() [2/2]

void Nektar::SolverUtils::RiemannSolver::SetAuxVec ( std::string  name,
RSVecFuncType  fp 
)
inline

Definition at line 125 of file RiemannSolver.h.

126  {
127  m_auxVec[name] = fp;
128  }

References m_auxVec, and CellMLToNektar.pycml::name.

◆ SetParam() [1/2]

template<typename FuncPointerT , typename ObjectPointerT >
void Nektar::SolverUtils::RiemannSolver::SetParam ( std::string  name,
FuncPointerT  func,
ObjectPointerT  obj 
)
inline

Definition at line 97 of file RiemannSolver.h.

100  {
101  m_params[name] = std::bind(func, obj);
102  }

References m_params, and CellMLToNektar.pycml::name.

◆ SetParam() [2/2]

void Nektar::SolverUtils::RiemannSolver::SetParam ( std::string  name,
RSParamFuncType  fp 
)
inline

Definition at line 104 of file RiemannSolver.h.

105  {
106  m_params[name] = fp;
107  }

References m_params, and CellMLToNektar.pycml::name.

◆ SetScalar() [1/2]

template<typename FuncPointerT , typename ObjectPointerT >
void Nektar::SolverUtils::RiemannSolver::SetScalar ( std::string  name,
FuncPointerT  func,
ObjectPointerT  obj 
)
inline

Definition at line 71 of file RiemannSolver.h.

74  {
75  m_scalars[name] = std::bind(func, obj);
76  }

References m_scalars, and CellMLToNektar.pycml::name.

◆ SetScalar() [2/2]

void Nektar::SolverUtils::RiemannSolver::SetScalar ( std::string  name,
RSScalarFuncType  fp 
)
inline

Definition at line 78 of file RiemannSolver.h.

79  {
80  m_scalars[name] = fp;
81  }

References m_scalars, and CellMLToNektar.pycml::name.

◆ SetVector() [1/2]

template<typename FuncPointerT , typename ObjectPointerT >
void Nektar::SolverUtils::RiemannSolver::SetVector ( std::string  name,
FuncPointerT  func,
ObjectPointerT  obj 
)
inline

Definition at line 84 of file RiemannSolver.h.

87  {
88  m_vectors[name] = std::bind(func, obj);
89  }

References m_vectors, and CellMLToNektar.pycml::name.

◆ SetVector() [2/2]

void Nektar::SolverUtils::RiemannSolver::SetVector ( std::string  name,
RSVecFuncType  fp 
)
inline

Definition at line 91 of file RiemannSolver.h.

92  {
93  m_vectors[name] = fp;
94  }

References m_vectors, and CellMLToNektar.pycml::name.

◆ Solve()

void Nektar::SolverUtils::RiemannSolver::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.

This routine calls the virtual function v_Solve to perform the Riemann solve. If the flag m_requiresRotation is set, then the velocity field is rotated to the normal direction to perform dimensional splitting, and the resulting fluxes are rotated back to the Cartesian directions before being returned. For the Rotation to work, the normal vectors "N" and the location of the vector components in Fwd "vecLocs"must be set via the SetAuxVec() method.

Parameters
FwdForwards trace space.
BwdBackwards trace space.
fluxResultant flux along trace space.

Definition at line 102 of file RiemannSolver.cpp.

107  {
108  if (m_requiresRotation)
109  {
110  ASSERTL1(CheckVectors("N"), "N not defined.");
111  ASSERTL1(CheckAuxVec("vecLocs"), "vecLocs not defined.");
112  const Array<OneD, const Array<OneD, NekDouble> > normals =
113  m_vectors["N"]();
114  const Array<OneD, const Array<OneD, NekDouble> > vecLocs =
115  m_auxVec["vecLocs"]();
116 
117  int nFields = Fwd .size();
118  int nPts = Fwd[0].size();
119 
120  if (m_rotStorage[0].size() != nFields ||
121  m_rotStorage[0][0].size() != nPts)
122  {
123  for (int i = 0; i < 3; ++i)
124  {
125  m_rotStorage[i] =
126  Array<OneD, Array<OneD, NekDouble> >(nFields);
127  for (int j = 0; j < nFields; ++j)
128  {
129  m_rotStorage[i][j] = Array<OneD, NekDouble>(nPts);
130  }
131  }
132  }
133 
134  rotateToNormal (Fwd, normals, vecLocs, m_rotStorage[0]);
135  rotateToNormal (Bwd, normals, vecLocs, m_rotStorage[1]);
136  v_Solve (nDim, m_rotStorage[0], m_rotStorage[1],
137  m_rotStorage[2]);
138  rotateFromNormal(m_rotStorage[2], normals, vecLocs, flux);
139  }
140  else
141  {
142  v_Solve(nDim, Fwd, Bwd, flux);
143  }
144  }
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.
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 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.

References ASSERTL1, CheckAuxVec(), CheckVectors(), m_auxVec, m_requiresRotation, m_rotStorage, m_vectors, rotateFromNormal(), rotateToNormal(), and v_Solve().

◆ v_CalcFluxJacobian()

void Nektar::SolverUtils::RiemannSolver::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 
)
protectedvirtual

Definition at line 559 of file RiemannSolver.cpp.

566  {
567  boost::ignore_unused(nDim,Fwd,Bwd,normals,FJac,BJac);
568  NEKERROR(ErrorUtil::efatal, "v_CalcFluxJacobian not specified.");
569  }
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:209

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by CalcFluxJacobian().

◆ v_Solve()

virtual void Nektar::SolverUtils::RiemannSolver::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 
)
protectedpure virtual

Member Data Documentation

◆ m_auxScal

std::map<std::string, RSScalarFuncType> Nektar::SolverUtils::RiemannSolver::m_auxScal
protected

Map of auxiliary scalar function types.

Definition at line 165 of file RiemannSolver.h.

Referenced by CheckAuxScal(), and SetAuxScal().

◆ m_auxVec

std::map<std::string, RSVecFuncType> Nektar::SolverUtils::RiemannSolver::m_auxVec
protected

Map of auxiliary vector function types.

Definition at line 167 of file RiemannSolver.h.

Referenced by CalcFluxJacobian(), CheckAuxVec(), SetAuxVec(), and Solve().

◆ m_params

std::map<std::string, RSParamFuncType > Nektar::SolverUtils::RiemannSolver::m_params
protected

◆ m_requiresRotation

bool Nektar::SolverUtils::RiemannSolver::m_requiresRotation
protected

◆ m_rotMat

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::RiemannSolver::m_rotMat
protected

Rotation matrices for each trace quadrature point.

Definition at line 169 of file RiemannSolver.h.

Referenced by GenerateRotationMatrices(), rotateFromNormal(), rotateToNormal(), and Nektar::RoeSolverSIMD::v_Solve().

◆ m_rotStorage

Array<OneD, Array<OneD, Array<OneD, NekDouble> > > Nektar::SolverUtils::RiemannSolver::m_rotStorage
protected

Rotation storage.

Definition at line 171 of file RiemannSolver.h.

Referenced by Solve().

◆ m_scalars

std::map<std::string, RSScalarFuncType> Nektar::SolverUtils::RiemannSolver::m_scalars
protected

◆ m_spacedim

int Nektar::SolverUtils::RiemannSolver::m_spacedim

Definition at line 145 of file RiemannSolver.h.

◆ m_vectors

std::map<std::string, RSVecFuncType> Nektar::SolverUtils::RiemannSolver::m_vectors
protected