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:
Inheritance graph
[legend]
Collaboration diagram for Nektar::SolverUtils::RiemannSolver:
Collaboration graph
[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)
 
std::map< std::string, RSScalarFuncType > & GetScalars ()
 
std::map< std::string, RSVecFuncType > & GetVectors ()
 
std::map< std::string, RSParamFuncType > & GetParams ()
 

Public Attributes

int m_spacedim
 

Protected Member Functions

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

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 62 of file RiemannSolver.h.

Constructor & Destructor Documentation

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

Definition at line 79 of file RiemannSolver.cpp.

79  : m_requiresRotation(false),
80  m_rotStorage (3)
81  {
82 
83  }
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...

Member Function Documentation

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 367 of file RiemannSolver.cpp.

References Nektar::iterator, and m_auxScal.

368  {
370  m_auxScal.find(name);
371 
372  return it != m_auxScal.end();
373  }
std::map< std::string, RSScalarFuncType > m_auxScal
Map of auxiliary scalar function types.
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
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 380 of file RiemannSolver.cpp.

References Nektar::iterator, and m_auxVec.

Referenced by Solve().

381  {
383  m_auxVec.find(name);
384 
385  return it != m_auxVec.end();
386  }
std::map< std::string, RSVecFuncType > m_auxVec
Map of auxiliary vector function types.
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
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 354 of file RiemannSolver.cpp.

References Nektar::iterator, and m_params.

Referenced by Nektar::LaxFriedrichsSolver::v_PointSolve(), and Nektar::UpwindSolver::v_PointSolve().

355  {
357  m_params.find(name);
358 
359  return it != m_params.end();
360  }
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
std::map< std::string, RSParamFuncType > m_params
Map of parameter function types.
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 328 of file RiemannSolver.cpp.

References Nektar::iterator, and m_scalars.

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

329  {
331  m_scalars.find(name);
332 
333  return it != m_scalars.end();
334  }
std::map< std::string, RSScalarFuncType > m_scalars
Map of scalar function types.
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
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 341 of file RiemannSolver.cpp.

References Nektar::iterator, and m_vectors.

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

342  {
344  m_vectors.find(name);
345 
346  return it != m_vectors.end();
347  }
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
std::map< std::string, RSVecFuncType > m_vectors
Map of vector function types.
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 436 of file RiemannSolver.cpp.

References CROSS, DOT, and EPSILON.

Referenced by GenerateRotationMatrices().

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

Generate rotation matrices for 3D expansions.

Definition at line 391 of file RiemannSolver.cpp.

References FromToRotation(), and m_rotMat.

Referenced by rotateToNormal().

393  {
394  Array<OneD, NekDouble> xdir(3,0.0);
395  Array<OneD, NekDouble> tn (3);
396  NekDouble tmp[9];
397  const int nq = normals[0].num_elements();
398  int i, j;
399  xdir[0] = 1.0;
400 
401  // Allocate storage for rotation matrices.
402  m_rotMat = Array<OneD, Array<OneD, NekDouble> >(9);
403 
404  for (i = 0; i < 9; ++i)
405  {
406  m_rotMat[i] = Array<OneD, NekDouble>(nq);
407  }
408  for (i = 0; i < normals[0].num_elements(); ++i)
409  {
410  // Generate matrix which takes us from (1,0,0) vector to trace
411  // normal.
412  tn[0] = normals[0][i];
413  tn[1] = normals[1][i];
414  tn[2] = normals[2][i];
415  FromToRotation(tn, xdir, tmp);
416 
417  for (j = 0; j < 9; ++j)
418  {
419  m_rotMat[j][i] = tmp[j];
420  }
421  }
422  }
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...
double NekDouble
Array< OneD, Array< OneD, NekDouble > > m_rotMat
Rotation matrices for each trace quadrature point.
std::map<std::string, RSParamFuncType>& Nektar::SolverUtils::RiemannSolver::GetParams ( )
inline

Definition at line 136 of file RiemannSolver.h.

References m_params.

137  {
138  return m_params;
139  }
std::map< std::string, RSParamFuncType > m_params
Map of parameter function types.
std::map<std::string, RSScalarFuncType>& Nektar::SolverUtils::RiemannSolver::GetScalars ( )
inline

Definition at line 126 of file RiemannSolver.h.

References m_scalars.

127  {
128  return m_scalars;
129  }
std::map< std::string, RSScalarFuncType > m_scalars
Map of scalar function types.
std::map<std::string, RSVecFuncType>& Nektar::SolverUtils::RiemannSolver::GetVectors ( )
inline

Definition at line 131 of file RiemannSolver.h.

References m_vectors.

132  {
133  return m_vectors;
134  }
std::map< std::string, RSVecFuncType > m_vectors
Map of vector function types.
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 251 of file RiemannSolver.cpp.

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

Referenced by Solve().

256  {
257  for (int i = 0; i < inarray.num_elements(); ++i)
258  {
259  Vmath::Vcopy(inarray[i].num_elements(), inarray[i], 1,
260  outarray[i], 1);
261  }
262 
263  for (int i = 0; i < vecLocs.num_elements(); i++)
264  {
265  ASSERTL1(vecLocs[i].num_elements() == normals.num_elements(),
266  "vecLocs[i] element count mismatch");
267 
268  switch (normals.num_elements())
269  {
270  case 1:
271  // do nothing
272  break;
273 
274  case 2:
275  {
276  const int nq = normals[0].num_elements();
277  const int vx = (int)vecLocs[i][0];
278  const int vy = (int)vecLocs[i][1];
279 
280  Vmath::Vmul (nq, inarray [vy], 1, normals [1], 1,
281  outarray[vx], 1);
282  Vmath::Vvtvm(nq, inarray [vx], 1, normals [0], 1,
283  outarray[vx], 1, outarray[vx], 1);
284  Vmath::Vmul (nq, inarray [vx], 1, normals [1], 1,
285  outarray[vy], 1);
286  Vmath::Vvtvp(nq, inarray [vy], 1, normals [0], 1,
287  outarray[vy], 1, outarray[vy], 1);
288  break;
289  }
290 
291  case 3:
292  {
293  const int nq = normals[0].num_elements();
294  const int vx = (int)vecLocs[i][0];
295  const int vy = (int)vecLocs[i][1];
296  const int vz = (int)vecLocs[i][2];
297 
298  Vmath::Vvtvvtp(nq, inarray [vx], 1, m_rotMat[0], 1,
299  inarray [vy], 1, m_rotMat[3], 1,
300  outarray[vx], 1);
301  Vmath::Vvtvp (nq, inarray [vz], 1, m_rotMat[6], 1,
302  outarray[vx], 1, outarray[vx], 1);
303  Vmath::Vvtvvtp(nq, inarray [vx], 1, m_rotMat[1], 1,
304  inarray [vy], 1, m_rotMat[4], 1,
305  outarray[vy], 1);
306  Vmath::Vvtvp (nq, inarray [vz], 1, m_rotMat[7], 1,
307  outarray[vy], 1, outarray[vy], 1);
308  Vmath::Vvtvvtp(nq, inarray [vx], 1, m_rotMat[2], 1,
309  inarray [vy], 1, m_rotMat[5], 1,
310  outarray[vz], 1);
311  Vmath::Vvtvp (nq, inarray [vz], 1, m_rotMat[8], 1,
312  outarray[vz], 1, outarray[vz], 1);
313  break;
314  }
315 
316  default:
317  ASSERTL1(false, "Invalid space dimension.");
318  break;
319  }
320  }
321  }
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:428
Array< OneD, Array< OneD, NekDouble > > m_rotMat
Rotation matrices for each trace quadrature point.
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:523
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:451
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1038
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:169
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 164 of file RiemannSolver.cpp.

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

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

169  {
170  for (int i = 0; i < inarray.num_elements(); ++i)
171  {
172  Vmath::Vcopy(inarray[i].num_elements(), inarray[i], 1,
173  outarray[i], 1);
174  }
175 
176  for (int i = 0; i < vecLocs.num_elements(); i++)
177  {
178  ASSERTL1(vecLocs[i].num_elements() == normals.num_elements(),
179  "vecLocs[i] element count mismatch");
180 
181  switch (normals.num_elements())
182  {
183  case 1:
184  // do nothing
185  break;
186 
187  case 2:
188  {
189  const int nq = inarray[0].num_elements();
190  const int vx = (int)vecLocs[i][0];
191  const int vy = (int)vecLocs[i][1];
192 
193  Vmath::Vmul (nq, inarray [vx], 1, normals [0], 1,
194  outarray[vx], 1);
195  Vmath::Vvtvp(nq, inarray [vy], 1, normals [1], 1,
196  outarray[vx], 1, outarray[vx], 1);
197  Vmath::Vmul (nq, inarray [vx], 1, normals [1], 1,
198  outarray[vy], 1);
199  Vmath::Vvtvm(nq, inarray [vy], 1, normals [0], 1,
200  outarray[vy], 1, outarray[vy], 1);
201  break;
202  }
203 
204  case 3:
205  {
206  const int nq = inarray[0].num_elements();
207  const int vx = (int)vecLocs[i][0];
208  const int vy = (int)vecLocs[i][1];
209  const int vz = (int)vecLocs[i][2];
210 
211  // Generate matrices if they don't already exist.
212  if (m_rotMat.num_elements() == 0)
213  {
214  GenerateRotationMatrices(normals);
215  }
216 
217  // Apply rotation matrices.
218  Vmath::Vvtvvtp(nq, inarray [vx], 1, m_rotMat[0], 1,
219  inarray [vy], 1, m_rotMat[1], 1,
220  outarray[vx], 1);
221  Vmath::Vvtvp (nq, inarray [vz], 1, m_rotMat[2], 1,
222  outarray[vx], 1, outarray[vx], 1);
223  Vmath::Vvtvvtp(nq, inarray [vx], 1, m_rotMat[3], 1,
224  inarray [vy], 1, m_rotMat[4], 1,
225  outarray[vy], 1);
226  Vmath::Vvtvp (nq, inarray [vz], 1, m_rotMat[5], 1,
227  outarray[vy], 1, outarray[vy], 1);
228  Vmath::Vvtvvtp(nq, inarray [vx], 1, m_rotMat[6], 1,
229  inarray [vy], 1, m_rotMat[7], 1,
230  outarray[vz], 1);
231  Vmath::Vvtvp (nq, inarray [vz], 1, m_rotMat[8], 1,
232  outarray[vz], 1, outarray[vz], 1);
233  break;
234  }
235 
236  default:
237  ASSERTL1(false, "Invalid space dimension.");
238  break;
239  }
240  }
241  }
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:428
void GenerateRotationMatrices(const Array< OneD, const Array< OneD, NekDouble > > &normals)
Generate rotation matrices for 3D expansions.
Array< OneD, Array< OneD, NekDouble > > m_rotMat
Rotation matrices for each trace quadrature point.
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:523
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:451
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1038
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:169
template<typename FuncPointerT , typename ObjectPointerT >
void Nektar::SolverUtils::RiemannSolver::SetAuxScal ( std::string  name,
FuncPointerT  func,
ObjectPointerT  obj 
)
inline

Definition at line 111 of file RiemannSolver.h.

References m_auxScal.

114  {
115  m_auxScal[name] = boost::bind(func, obj);
116  }
std::map< std::string, RSScalarFuncType > m_auxScal
Map of auxiliary scalar function types.
template<typename FuncPointerT , typename ObjectPointerT >
void Nektar::SolverUtils::RiemannSolver::SetAuxVec ( std::string  name,
FuncPointerT  func,
ObjectPointerT  obj 
)
inline

Definition at line 119 of file RiemannSolver.h.

References m_auxVec.

122  {
123  m_auxVec[name] = boost::bind(func, obj);
124  }
std::map< std::string, RSVecFuncType > m_auxVec
Map of auxiliary vector function types.
template<typename FuncPointerT , typename ObjectPointerT >
void Nektar::SolverUtils::RiemannSolver::SetParam ( std::string  name,
FuncPointerT  func,
ObjectPointerT  obj 
)
inline

Definition at line 98 of file RiemannSolver.h.

References m_params.

101  {
102  m_params[name] = boost::bind(func, obj);
103  }
std::map< std::string, RSParamFuncType > m_params
Map of parameter function types.
void Nektar::SolverUtils::RiemannSolver::SetParam ( std::string  name,
RSParamFuncType  fp 
)
inline

Definition at line 105 of file RiemannSolver.h.

References m_params.

106  {
107  m_params[name] = fp;
108  }
std::map< std::string, RSParamFuncType > m_params
Map of parameter function types.
template<typename FuncPointerT , typename ObjectPointerT >
void Nektar::SolverUtils::RiemannSolver::SetScalar ( std::string  name,
FuncPointerT  func,
ObjectPointerT  obj 
)
inline

Definition at line 72 of file RiemannSolver.h.

References m_scalars.

75  {
76  m_scalars[name] = boost::bind(func, obj);
77  }
std::map< std::string, RSScalarFuncType > m_scalars
Map of scalar function types.
void Nektar::SolverUtils::RiemannSolver::SetScalar ( std::string  name,
RSScalarFuncType  fp 
)
inline

Definition at line 79 of file RiemannSolver.h.

References m_scalars.

80  {
81  m_scalars[name] = fp;
82  }
std::map< std::string, RSScalarFuncType > m_scalars
Map of scalar function types.
template<typename FuncPointerT , typename ObjectPointerT >
void Nektar::SolverUtils::RiemannSolver::SetVector ( std::string  name,
FuncPointerT  func,
ObjectPointerT  obj 
)
inline

Definition at line 85 of file RiemannSolver.h.

References m_vectors.

88  {
89  m_vectors[name] = boost::bind(func, obj);
90  }
std::map< std::string, RSVecFuncType > m_vectors
Map of vector function types.
void Nektar::SolverUtils::RiemannSolver::SetVector ( std::string  name,
RSVecFuncType  fp 
)
inline

Definition at line 92 of file RiemannSolver.h.

References m_vectors.

93  {
94  m_vectors[name] = fp;
95  }
std::map< std::string, RSVecFuncType > m_vectors
Map of vector function types.
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 101 of file RiemannSolver.cpp.

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

106  {
107  if (m_requiresRotation)
108  {
109  ASSERTL1(CheckVectors("N"), "N not defined.");
110  ASSERTL1(CheckAuxVec("vecLocs"), "vecLocs not defined.");
111  const Array<OneD, const Array<OneD, NekDouble> > normals =
112  m_vectors["N"]();
113  const Array<OneD, const Array<OneD, NekDouble> > vecLocs =
114  m_auxVec["vecLocs"]();
115 
116  int nFields = Fwd .num_elements();
117  int nPts = Fwd[0].num_elements();
118 
119  if (m_rotStorage[0].num_elements() != nFields ||
120  m_rotStorage[0][0].num_elements() != nPts)
121  {
122  for (int i = 0; i < 3; ++i)
123  {
124  m_rotStorage[i] =
125  Array<OneD, Array<OneD, NekDouble> >(nFields);
126  for (int j = 0; j < nFields; ++j)
127  {
128  m_rotStorage[i][j] = Array<OneD, NekDouble>(nPts);
129  }
130  }
131  }
132 
133  rotateToNormal (Fwd, normals, vecLocs, m_rotStorage[0]);
134  rotateToNormal (Bwd, normals, vecLocs, m_rotStorage[1]);
135  v_Solve (nDim, m_rotStorage[0], m_rotStorage[1],
136  m_rotStorage[2]);
137  rotateFromNormal(m_rotStorage[2], normals, vecLocs, flux);
138  }
139  else
140  {
141  v_Solve(nDim, Fwd, Bwd, flux);
142  }
143  }
SOLVER_UTILS_EXPORT bool CheckVectors(std::string name)
Determine whether a vector has been defined in m_vectors.
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.
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
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_rotStorage
Rotation storage.
SOLVER_UTILS_EXPORT bool CheckAuxVec(std::string name)
Determine whether a vector has been defined in m_auxVec.
std::map< std::string, RSVecFuncType > m_auxVec
Map of auxiliary vector function types.
bool m_requiresRotation
Indicates whether the Riemann solver requires a rotation to be applied to the velocity fields...
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
std::map< std::string, RSVecFuncType > m_vectors
Map of vector function types.
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 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

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

Map of auxiliary scalar function types.

Definition at line 154 of file RiemannSolver.h.

Referenced by CheckAuxScal(), and SetAuxScal().

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

Map of auxiliary vector function types.

Definition at line 156 of file RiemannSolver.h.

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

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

Indicates whether the Riemann solver requires a rotation to be applied to the velocity fields.

Definition at line 146 of file RiemannSolver.h.

Referenced by Nektar::APESolver::APESolver(), Nektar::CompressibleSolver::CompressibleSolver(), Nektar::LinearSWESolver::LinearSWESolver(), Nektar::NonlinearSWESolver::NonlinearSWESolver(), and Solve().

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

Rotation matrices for each trace quadrature point.

Definition at line 158 of file RiemannSolver.h.

Referenced by GenerateRotationMatrices(), rotateFromNormal(), and rotateToNormal().

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

Rotation storage.

Definition at line 160 of file RiemannSolver.h.

Referenced by Solve().

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

Definition at line 141 of file RiemannSolver.h.

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

Map of vector function types.

Definition at line 150 of file RiemannSolver.h.

Referenced by CheckVectors(), Nektar::APESolver::GetRotBasefield(), GetVectors(), SetVector(), and Solve().