Nektar++
Advection.h
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File: Advection.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 // Permission is hereby granted, free of charge, to any person obtaining a
14 // copy of this software and associated documentation files (the "Software"),
15 // to deal in the Software without restriction, including without limitation
16 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
17 // and/or sell copies of the Software, and to permit persons to whom the
18 // Software is furnished to do so, subject to the following conditions:
19 //
20 // The above copyright notice and this permission notice shall be included
21 // in all copies or substantial portions of the Software.
22 //
23 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
24 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
25 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
26 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
27 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
28 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
29 // DEALINGS IN THE SOFTWARE.
30 //
31 // Description: Abstract base class for advection.
32 //
33 ///////////////////////////////////////////////////////////////////////////////
34 
35 #ifndef NEKTAR_SOLVERUTILS_ADVECTION
36 #define NEKTAR_SOLVERUTILS_ADVECTION
37 
38 #include <functional>
39 #include <string>
40 
43 #include <MultiRegions/ExpList.h>
46 #include <iomanip>
47 
48 using namespace std;
49 
50 namespace Nektar
51 {
52 namespace SolverUtils
53 {
54 
55 /**
56  * Defines a callback function type which evaluates the flux vector \f$ F(u) \f$
57  * in a conservative advection of the form \f$ \nabla\cdot F(u) \f$.
58  */
59 typedef std::function<void(const Array<OneD, Array<OneD, NekDouble>> &,
60  Array<OneD, Array<OneD, Array<OneD, NekDouble>>> &)>
62 
63 /**
64  * @brief An abstract base class encapsulating the concept of advection
65  * of a vector field.
66  *
67  * Subclasses override the Advection::v_InitObject function to
68  * initialise the object and the Advection::v_Advect function to
69  * evaluate the advection of the vector field.
70  */
71 class Advection
72 {
73 public:
75 
76  /// Interface function to initialise the advection object.
77  SOLVER_UTILS_EXPORT void InitObject(
80 
81  /// Interface function to advect the vector field.
82  SOLVER_UTILS_EXPORT void Advect(
83  const int nConvectiveFields,
85  const Array<OneD, Array<OneD, NekDouble>> &advVel,
86  const Array<OneD, Array<OneD, NekDouble>> &inarray,
87  Array<OneD, Array<OneD, NekDouble>> &outarray, const NekDouble &time,
88  const Array<OneD, Array<OneD, NekDouble>> &pFwd =
90  const Array<OneD, Array<OneD, NekDouble>> &pBwd =
92 
93  /// Interface function to advect the Volume field.
95  const int nConvectiveFields,
97  const Array<OneD, Array<OneD, NekDouble>> &pAdvVel,
98  const Array<OneD, Array<OneD, NekDouble>> &pInarray,
99  TensorOfArray3D<NekDouble> &pVolumeFlux, const NekDouble &pTime)
100  {
101  v_AdvectVolumeFlux(nConvectiveFields, pFields, pAdvVel, pInarray,
102  pVolumeFlux, pTime);
103  }
104 
105  /// Interface function to advect the Trace field.
107  const int nConvectiveFields,
109  const Array<OneD, Array<OneD, NekDouble>> &pAdvVel,
110  const Array<OneD, Array<OneD, NekDouble>> &pInarray,
111  Array<OneD, Array<OneD, NekDouble>> &pTraceFlux, const NekDouble &pTime,
112  const Array<OneD, Array<OneD, NekDouble>> &pFwd =
114  const Array<OneD, Array<OneD, NekDouble>> &pBwd =
116  {
117  v_AdvectTraceFlux(nConvectiveFields, pFields, pAdvVel, pInarray,
118  pTraceFlux, pTime, pFwd, pBwd);
119  }
120 
121  /**
122  * @brief Similar with Advection::Advect(): calculate the advection flux
123  * The difference is in the outarray:
124  * it is the coefficients of basis for AdvectCoeffs()
125  * it is the physics on quadrature points for Advect()
126  */
127  SOLVER_UTILS_EXPORT void AdvectCoeffs(
128  const int nConvectiveFields,
130  const Array<OneD, Array<OneD, NekDouble>> &advVel,
131  const Array<OneD, Array<OneD, NekDouble>> &inarray,
132  Array<OneD, Array<OneD, NekDouble>> &outarray, const NekDouble &time,
133  const Array<OneD, Array<OneD, NekDouble>> &pFwd =
135  const Array<OneD, Array<OneD, NekDouble>> &pBwd =
137 
138  /**
139  * @brief Set the flux vector callback function.
140  *
141  * This routine is a utility function to avoid the explicit use of
142  * std::bind. A function and object can be passed to this function
143  * instead.
144  */
145  template <typename FuncPointerT, typename ObjectPointerT>
146  void SetFluxVector(FuncPointerT func, ObjectPointerT obj)
147  {
148  m_fluxVector =
149  std::bind(func, obj, std::placeholders::_1, std::placeholders::_2);
150  }
151 
152  /**
153  * @brief Set a Riemann solver object for this advection object.
154  *
155  * @param riemann The RiemannSolver object.
156  */
158  {
159  m_riemann = riemann;
160  }
161 
162  /**
163  * @brief Set the flux vector callback function.
164  *
165  * @param fluxVector The callback function to override.
166  */
167  inline void SetFluxVector(AdvectionFluxVecCB fluxVector)
168  {
169  m_fluxVector = fluxVector;
170  }
171 
172  /**
173  * @brief Set the base flow used for linearised advection objects.
174  *
175  * @param inarray Vector to use as baseflow
176  */
177  inline void SetBaseFlow(
178  const Array<OneD, Array<OneD, NekDouble>> &inarray,
180  {
181  v_SetBaseFlow(inarray, fields);
182  }
183 
184  template <typename DataType, typename TypeNekBlkMatSharedPtr>
186  const int nConvectiveFields, const int nSpaceDim,
188  const Array<OneD, TypeNekBlkMatSharedPtr> &TracePntJacCons,
190  const Array<OneD, TypeNekBlkMatSharedPtr> &TracePntJacGrad,
191  const Array<OneD, Array<OneD, DataType>> &TracePntJacGradSign);
192 
193  template <typename DataType, typename TypeNekBlkMatSharedPtr>
194  void CalcJacobTraceInteg(
195  const Array<OneD, MultiRegions::ExpListSharedPtr> &pFields, const int m,
196  const int n, const Array<OneD, const TypeNekBlkMatSharedPtr> &PntJac,
197  const Array<OneD, const Array<OneD, DataType>> &PntJacSign,
198  Array<OneD, DNekMatSharedPtr> &TraceJacFwd,
199  Array<OneD, DNekMatSharedPtr> &TraceJacBwd);
200 
203  const int &nConvectiveFields,
204  const TensorOfArray5D<NekDouble> &ElmtJacArray,
206  {
207  v_AddVolumJacToMat(pFields, nConvectiveFields, ElmtJacArray, gmtxarray);
208  }
209 
210 protected:
211  /// Callback function to the flux vector (set when advection is in
212  /// conservative form).
214  /// Riemann solver for DG-type schemes.
216  /// Storage for space dimension. Used for homogeneous extension.
218 
219  /// Initialises the advection object.
220  SOLVER_UTILS_EXPORT virtual void v_InitObject(
223 
224  /// Advects a vector field.
226  const int nConvectiveFields,
228  const Array<OneD, Array<OneD, NekDouble>> &advVel,
229  const Array<OneD, Array<OneD, NekDouble>> &inarray,
230  Array<OneD, Array<OneD, NekDouble>> &outarray, const NekDouble &time,
231  const Array<OneD, Array<OneD, NekDouble>> &pFwd =
233  const Array<OneD, Array<OneD, NekDouble>> &pBwd =
235 
236  /// Advects Volume Flux.
237  SOLVER_UTILS_EXPORT virtual void v_AdvectVolumeFlux(
238  const int nConvectiveFields,
240  const Array<OneD, Array<OneD, NekDouble>> &advVel,
241  const Array<OneD, Array<OneD, NekDouble>> &inarray,
242  TensorOfArray3D<NekDouble> &pVolumeFlux, const NekDouble &time);
243 
244  /// Advects Trace Flux.
245  SOLVER_UTILS_EXPORT virtual void v_AdvectTraceFlux(
246  const int nConvectiveFields,
248  const Array<OneD, Array<OneD, NekDouble>> &advVel,
249  const Array<OneD, Array<OneD, NekDouble>> &inarray,
250  Array<OneD, Array<OneD, NekDouble>> &pTraceFlux, const NekDouble &time,
251  const Array<OneD, Array<OneD, NekDouble>> &pFwd =
253  const Array<OneD, Array<OneD, NekDouble>> &pBwd =
255 
256  SOLVER_UTILS_EXPORT virtual void v_AdvectCoeffs(
257  const int nConvectiveFields,
259  const Array<OneD, Array<OneD, NekDouble>> &advVel,
260  const Array<OneD, Array<OneD, NekDouble>> &inarray,
261  Array<OneD, Array<OneD, NekDouble>> &outarray, const NekDouble &time,
262  const Array<OneD, Array<OneD, NekDouble>> &pFwd =
264  const Array<OneD, Array<OneD, NekDouble>> &pBwd =
266 
267  /// Overrides the base flow used during linearised advection
268  SOLVER_UTILS_EXPORT virtual void v_SetBaseFlow(
269  const Array<OneD, Array<OneD, NekDouble>> &inarray,
271 
272  SOLVER_UTILS_EXPORT virtual void v_AddVolumJacToMat(
274  const int &nConvectiveFields,
275  const TensorOfArray5D<NekDouble> &ElmtJacArray,
277 };
278 
279 /// A shared pointer to an Advection object.
280 typedef std::shared_ptr<Advection> AdvectionSharedPtr;
281 
282 /// Datatype of the NekFactory used to instantiate classes derived
283 /// from the Advection class.
286 
287 /// Gets the factory for initialising advection objects.
289 
290 } // namespace SolverUtils
291 } // namespace Nektar
292 
293 #endif
#define SOLVER_UTILS_EXPORT
Provides a generic Factory class.
Definition: NekFactory.hpp:105
An abstract base class encapsulating the concept of advection of a vector field.
Definition: Advection.h:72
void SetBaseFlow(const Array< OneD, Array< OneD, NekDouble >> &inarray, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields)
Set the base flow used for linearised advection objects.
Definition: Advection.h:177
int m_spaceDim
Storage for space dimension. Used for homogeneous extension.
Definition: Advection.h:217
virtual SOLVER_UTILS_EXPORT ~Advection()
Definition: Advection.h:74
AdvectionFluxVecCB m_fluxVector
Callback function to the flux vector (set when advection is in conservative form).
Definition: Advection.h:213
virtual SOLVER_UTILS_EXPORT void v_Advect(const int nConvectiveFields, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble >> &advVel, const Array< OneD, Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray, const NekDouble &time, const Array< OneD, Array< OneD, NekDouble >> &pFwd=NullNekDoubleArrayOfArray, const Array< OneD, Array< OneD, NekDouble >> &pBwd=NullNekDoubleArrayOfArray)=0
Advects a vector field.
SOLVER_UTILS_EXPORT void AdvectTraceFlux(const int nConvectiveFields, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, Array< OneD, NekDouble >> &pAdvVel, const Array< OneD, Array< OneD, NekDouble >> &pInarray, Array< OneD, Array< OneD, NekDouble >> &pTraceFlux, const NekDouble &pTime, const Array< OneD, Array< OneD, NekDouble >> &pFwd=NullNekDoubleArrayOfArray, const Array< OneD, Array< OneD, NekDouble >> &pBwd=NullNekDoubleArrayOfArray)
Interface function to advect the Trace field.
Definition: Advection.h:106
SOLVER_UTILS_EXPORT void AddVolumJacToMat(const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int &nConvectiveFields, const TensorOfArray5D< NekDouble > &ElmtJacArray, Array< OneD, Array< OneD, SNekBlkMatSharedPtr >> &gmtxarray)
Definition: Advection.h:201
SOLVER_UTILS_EXPORT void AdvectVolumeFlux(const int nConvectiveFields, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, Array< OneD, NekDouble >> &pAdvVel, const Array< OneD, Array< OneD, NekDouble >> &pInarray, TensorOfArray3D< NekDouble > &pVolumeFlux, const NekDouble &pTime)
Interface function to advect the Volume field.
Definition: Advection.h:94
SOLVER_UTILS_EXPORT void AddTraceJacToMat(const int nConvectiveFields, const int nSpaceDim, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, TypeNekBlkMatSharedPtr > &TracePntJacCons, Array< OneD, Array< OneD, TypeNekBlkMatSharedPtr >> &gmtxarray, const Array< OneD, TypeNekBlkMatSharedPtr > &TracePntJacGrad, const Array< OneD, Array< OneD, DataType >> &TracePntJacGradSign)
RiemannSolverSharedPtr m_riemann
Riemann solver for DG-type schemes.
Definition: Advection.h:215
void SetFluxVector(AdvectionFluxVecCB fluxVector)
Set the flux vector callback function.
Definition: Advection.h:167
void SetRiemannSolver(RiemannSolverSharedPtr riemann)
Set a Riemann solver object for this advection object.
Definition: Advection.h:157
void SetFluxVector(FuncPointerT func, ObjectPointerT obj)
Set the flux vector callback function.
Definition: Advection.h:146
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::function< void(const Array< OneD, Array< OneD, NekDouble >> &, Array< OneD, Array< OneD, Array< OneD, NekDouble >>> &)> AdvectionFluxVecCB
Definition: Advection.h:61
std::shared_ptr< RiemannSolver > RiemannSolverSharedPtr
A shared pointer to an EquationSystem object.
AdvectionFactory & GetAdvectionFactory()
Gets the factory for initialising advection objects.
Definition: Advection.cpp:47
LibUtilities::NekFactory< std::string, Advection, std::string > AdvectionFactory
Datatype of the NekFactory used to instantiate classes derived from the Advection class.
Definition: Advection.h:285
std::shared_ptr< Advection > AdvectionSharedPtr
A shared pointer to an Advection object.
Definition: Advection.h:280
The above copyright notice and this permission notice shall be included.
Definition: CoupledSolver.h:1
static Array< OneD, Array< OneD, NekDouble > > NullNekDoubleArrayOfArray
double NekDouble