Nektar++
RiemannSolver.h
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File: RiemannSolver.h
4 //
5 // For more information, please see: http://www.nektar.info
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description: Abstract base class for Riemann solvers with factory.
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 #ifndef NEKTAR_SOLVERUTILS_RIEMANNSOLVER
37 #define NEKTAR_SOLVERUTILS_RIEMANNSOLVER
38 
42 
43 #include <string>
44 
45 #include <boost/function.hpp>
46 #include <boost/call_traits.hpp>
47 
48 namespace Nektar
49 {
50  template <typename Dim, typename DataType>
51  class Array;
52 
53  namespace SolverUtils
54  {
55  typedef boost::function<
57  typedef boost::function<
59  typedef boost::function<
61 
63  {
64  public:
66  const int nDim,
67  const Array<OneD, const Array<OneD, NekDouble> > &Fwd,
68  const Array<OneD, const Array<OneD, NekDouble> > &Bwd,
70 
71  template<typename FuncPointerT, typename ObjectPointerT>
72  void SetScalar(std::string name,
73  FuncPointerT func,
74  ObjectPointerT obj)
75  {
76  m_scalars[name] = boost::bind(func, obj);
77  }
78 
79  void SetScalar(std::string name, RSScalarFuncType fp)
80  {
81  m_scalars[name] = fp;
82  }
83 
84  template<typename FuncPointerT, typename ObjectPointerT>
85  void SetVector(std::string name,
86  FuncPointerT func,
87  ObjectPointerT obj)
88  {
89  m_vectors[name] = boost::bind(func, obj);
90  }
91 
92  void SetVector(std::string name, RSVecFuncType fp)
93  {
94  m_vectors[name] = fp;
95  }
96 
97  template<typename FuncPointerT, typename ObjectPointerT>
98  void SetParam(std::string name,
99  FuncPointerT func,
100  ObjectPointerT obj)
101  {
102  m_params[name] = boost::bind(func, obj);
103  }
104 
105  void SetParam(std::string name, RSParamFuncType fp)
106  {
107  m_params[name] = fp;
108  }
109 
110  template<typename FuncPointerT, typename ObjectPointerT>
111  void SetAuxScal(std::string name,
112  FuncPointerT func,
113  ObjectPointerT obj)
114  {
115  m_auxScal[name] = boost::bind(func, obj);
116  }
117 
118  template<typename FuncPointerT, typename ObjectPointerT>
119  void SetAuxVec(std::string name,
120  FuncPointerT func,
121  ObjectPointerT obj)
122  {
123  m_auxVec[name] = boost::bind(func, obj);
124  }
125 
126  std::map<std::string, RSScalarFuncType> &GetScalars()
127  {
128  return m_scalars;
129  }
130 
131  std::map<std::string, RSVecFuncType> &GetVectors()
132  {
133  return m_vectors;
134  }
135 
136  std::map<std::string, RSParamFuncType> &GetParams()
137  {
138  return m_params;
139  }
140 
142 
143  protected:
144  /// Indicates whether the Riemann solver requires a rotation to be
145  /// applied to the velocity fields.
147  /// Map of scalar function types.
148  std::map<std::string, RSScalarFuncType> m_scalars;
149  /// Map of vector function types.
150  std::map<std::string, RSVecFuncType> m_vectors;
151  /// Map of parameter function types.
152  std::map<std::string, RSParamFuncType > m_params;
153  /// Map of auxiliary scalar function types.
154  std::map<std::string, RSScalarFuncType> m_auxScal;
155  /// Map of auxiliary vector function types.
156  std::map<std::string, RSVecFuncType> m_auxVec;
157  /// Rotation matrices for each trace quadrature point.
159  /// Rotation storage
161 
163 
164  virtual void v_Solve(
165  const int nDim,
166  const Array<OneD, const Array<OneD, NekDouble> > &Fwd,
167  const Array<OneD, const Array<OneD, NekDouble> > &Bwd,
168  Array<OneD, Array<OneD, NekDouble> > &flux) = 0;
169 
171  const Array<OneD, const Array<OneD, NekDouble> > &normals);
172  void FromToRotation(
175  NekDouble *mat);
177  const Array<OneD, const Array<OneD, NekDouble> > &inarray,
178  const Array<OneD, const Array<OneD, NekDouble> > &normals,
179  const Array<OneD, const Array<OneD, NekDouble> > &vecLocs,
180  Array<OneD, Array<OneD, NekDouble> > &outarray);
182  const Array<OneD, const Array<OneD, NekDouble> > &inarray,
183  const Array<OneD, const Array<OneD, NekDouble> > &normals,
184  const Array<OneD, const Array<OneD, NekDouble> > &vecLocs,
185  Array<OneD, Array<OneD, NekDouble> > &outarray);
186  SOLVER_UTILS_EXPORT bool CheckScalars (std::string name);
187  SOLVER_UTILS_EXPORT bool CheckVectors (std::string name);
188  SOLVER_UTILS_EXPORT bool CheckParams (std::string name);
189  SOLVER_UTILS_EXPORT bool CheckAuxScal (std::string name);
190  SOLVER_UTILS_EXPORT bool CheckAuxVec (std::string name);
191  };
192 
193  /// A shared pointer to an EquationSystem object
194  typedef boost::shared_ptr<RiemannSolver> RiemannSolverSharedPtr;
195  /// Datatype of the NekFactory used to instantiate classes derived
196  /// from the RiemannSolver class.
200  }
201 }
202 
203 #endif
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...
boost::function< const Array< OneD, const NekDouble > &()> RSScalarFuncType
Definition: RiemannSolver.h:56
SOLVER_UTILS_EXPORT bool CheckScalars(std::string name)
Determine whether a scalar has been defined in m_scalars.
boost::function< const Array< OneD, const Array< OneD, NekDouble > > &()> RSVecFuncType
Definition: RiemannSolver.h:58
void SetParam(std::string name, FuncPointerT func, ObjectPointerT obj)
Definition: RiemannSolver.h:98
SOLVER_UTILS_EXPORT bool CheckVectors(std::string name)
Determine whether a vector has been defined in m_vectors.
boost::shared_ptr< RiemannSolver > RiemannSolverSharedPtr
A shared pointer to an EquationSystem object.
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.
void SetAuxVec(std::string name, FuncPointerT func, ObjectPointerT obj)
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
LibUtilities::NekFactory< std::string, RiemannSolver > RiemannSolverFactory
Datatype of the NekFactory used to instantiate classes derived from the RiemannSolver class...
std::map< std::string, RSParamFuncType > & GetParams()
void SetVector(std::string name, FuncPointerT func, ObjectPointerT obj)
Definition: RiemannSolver.h:85
std::map< std::string, RSScalarFuncType > m_scalars
Map of scalar function types.
void GenerateRotationMatrices(const Array< OneD, const Array< OneD, NekDouble > > &normals)
Generate rotation matrices for 3D expansions.
void SetParam(std::string name, RSParamFuncType fp)
void SetScalar(std::string name, RSScalarFuncType fp)
Definition: RiemannSolver.h:79
RiemannSolverFactory & GetRiemannSolverFactory()
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.
std::map< std::string, RSScalarFuncType > m_auxScal
Map of auxiliary scalar function types.
void SetAuxScal(std::string name, FuncPointerT func, ObjectPointerT obj)
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_rotStorage
Rotation storage.
double NekDouble
SOLVER_UTILS_EXPORT bool CheckAuxScal(std::string name)
Determine whether a scalar has been defined in m_auxScal.
Array< OneD, Array< OneD, NekDouble > > m_rotMat
Rotation matrices for each trace quadrature point.
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.
SOLVER_UTILS_EXPORT RiemannSolver()
bool m_requiresRotation
Indicates whether the Riemann solver requires a rotation to be applied to the velocity fields...
The RiemannSolver class provides an abstract interface under which solvers for various Riemann proble...
Definition: RiemannSolver.h:62
#define SOLVER_UTILS_EXPORT
void SetScalar(std::string name, FuncPointerT func, ObjectPointerT obj)
Definition: RiemannSolver.h:72
std::map< std::string, RSParamFuncType > m_params
Map of parameter function types.
SOLVER_UTILS_EXPORT bool CheckParams(std::string name)
Determine whether a parameter has been defined in m_params.
void SetVector(std::string name, RSVecFuncType fp)
Definition: RiemannSolver.h:92
std::map< std::string, RSVecFuncType > & GetVectors()
std::map< std::string, RSScalarFuncType > & GetScalars()
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.
boost::function< NekDouble()> RSParamFuncType
Definition: RiemannSolver.h:60
Provides a generic Factory class.
Definition: NekFactory.hpp:116