Nektar++
Public Member Functions | Static Public Member Functions | Static Public Attributes | Protected Member Functions | Static Protected Attributes | Private Attributes | Friends | List of all members
Nektar::NavierStokesAdvection Class Reference

#include <NavierStokesAdvection.h>

Inheritance diagram for Nektar::NavierStokesAdvection:
[legend]

Public Member Functions

void SetSpecHPDealiasing (bool value)
 
- Public Member Functions inherited from Nektar::SolverUtils::Advection
virtual SOLVER_UTILS_EXPORT ~Advection ()
 
SOLVER_UTILS_EXPORT void InitObject (LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)
 Interface function to initialise the advection object. More...
 
SOLVER_UTILS_EXPORT void 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)
 Interface function to advect the vector field. More...
 
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. More...
 
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. More...
 
SOLVER_UTILS_EXPORT void AdvectCoeffs (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)
 Similar with Advection::Advect(): calculate the advection flux The difference is in the outarray: it is the coefficients of basis for AdvectCoeffs() it is the physics on quadrature points for Advect() More...
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetFluxVector (FuncPointerT func, ObjectPointerT obj)
 Set the flux vector callback function. More...
 
void SetRiemannSolver (RiemannSolverSharedPtr riemann)
 Set a Riemann solver object for this advection object. More...
 
void SetFluxVector (AdvectionFluxVecCB fluxVector)
 Set the flux vector callback function. More...
 
void SetBaseFlow (const Array< OneD, Array< OneD, NekDouble >> &inarray, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields)
 Set the base flow used for linearised advection objects. More...
 
template<typename DataType , typename TypeNekBlkMatSharedPtr >
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)
 
template<typename DataType , typename TypeNekBlkMatSharedPtr >
void CalcJacobTraceInteg (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int m, const int n, const Array< OneD, const TypeNekBlkMatSharedPtr > &PntJac, const Array< OneD, const Array< OneD, DataType >> &PntJacSign, Array< OneD, DNekMatSharedPtr > &TraceJacFwd, Array< OneD, DNekMatSharedPtr > &TraceJacBwd)
 
SOLVER_UTILS_EXPORT void AddVolumJacToMat (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int &nConvectiveFields, const TensorOfArray5D< NekDouble > &ElmtJacArray, Array< OneD, Array< OneD, SNekBlkMatSharedPtr >> &gmtxarray)
 
template<typename DataType , typename TypeNekBlkMatSharedPtr >
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)
 

Static Public Member Functions

static SolverUtils::AdvectionSharedPtr create (std::string)
 Creates an instance of this class. More...
 

Static Public Attributes

static std::string className
 Name of class. More...
 
static std::string className2
 

Protected Member Functions

 NavierStokesAdvection ()
 
virtual ~NavierStokesAdvection ()
 
virtual void v_InitObject (LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields) override
 Initialises the advection object. More...
 
virtual 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) override
 Advects a vector field. More...
 
- Protected Member Functions inherited from Nektar::SolverUtils::Advection
virtual SOLVER_UTILS_EXPORT void v_AdvectVolumeFlux (const int nConvectiveFields, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble >> &advVel, const Array< OneD, Array< OneD, NekDouble >> &inarray, TensorOfArray3D< NekDouble > &pVolumeFlux, const NekDouble &time)
 Advects Volume Flux. More...
 
virtual SOLVER_UTILS_EXPORT void v_AdvectTraceFlux (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 >> &pTraceFlux, const NekDouble &time, const Array< OneD, Array< OneD, NekDouble >> &pFwd=NullNekDoubleArrayOfArray, const Array< OneD, Array< OneD, NekDouble >> &pBwd=NullNekDoubleArrayOfArray)
 Advects Trace Flux. More...
 
virtual SOLVER_UTILS_EXPORT void v_AdvectCoeffs (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)
 
virtual SOLVER_UTILS_EXPORT void v_SetBaseFlow (const Array< OneD, Array< OneD, NekDouble >> &inarray, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields)
 Overrides the base flow used during linearised advection. More...
 
virtual SOLVER_UTILS_EXPORT void v_AddVolumJacToMat (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int &nConvectiveFields, const TensorOfArray5D< NekDouble > &ElmtJacArray, Array< OneD, Array< OneD, SNekBlkMatSharedPtr >> &gmtxarray)
 

Static Protected Attributes

static std::string navierStokesAdvectionTypeLookupIds []
 

Private Attributes

bool m_specHP_dealiasing
 
bool m_homogen_dealiasing
 
bool m_SingleMode
 
bool m_HalfMode
 

Friends

class MemoryManager< NavierStokesAdvection >
 

Additional Inherited Members

- Protected Attributes inherited from Nektar::SolverUtils::Advection
AdvectionFluxVecCB m_fluxVector
 Callback function to the flux vector (set when advection is in conservative form). More...
 
RiemannSolverSharedPtr m_riemann
 Riemann solver for DG-type schemes. More...
 
int m_spaceDim
 Storage for space dimension. Used for homogeneous extension. More...
 

Detailed Description

Definition at line 43 of file NavierStokesAdvection.h.

Constructor & Destructor Documentation

◆ NavierStokesAdvection()

Nektar::NavierStokesAdvection::NavierStokesAdvection ( )
protected

Constructor. Creates ...

Parameters

Definition at line 63 of file NavierStokesAdvection.cpp.

63  : Advection()
64 
65 {
66 }

◆ ~NavierStokesAdvection()

Nektar::NavierStokesAdvection::~NavierStokesAdvection ( )
protectedvirtual

Definition at line 68 of file NavierStokesAdvection.cpp.

69 {
70 }

Member Function Documentation

◆ create()

static SolverUtils::AdvectionSharedPtr Nektar::NavierStokesAdvection::create ( std::string  )
inlinestatic

Creates an instance of this class.

Definition at line 50 of file NavierStokesAdvection.h.

51  {
53  }
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr().

◆ SetSpecHPDealiasing()

void Nektar::NavierStokesAdvection::SetSpecHPDealiasing ( bool  value)
inline

Definition at line 59 of file NavierStokesAdvection.h.

60  {
61  m_specHP_dealiasing = value;
62  }

References m_specHP_dealiasing.

◆ v_Advect()

void Nektar::NavierStokesAdvection::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 
)
overrideprotectedvirtual

Advects a vector field.

Implements Nektar::SolverUtils::Advection.

Definition at line 86 of file NavierStokesAdvection.cpp.

94 {
95  boost::ignore_unused(time, pFwd, pBwd);
96 
97  int nqtot = fields[0]->GetTotPoints();
98  ASSERTL1(nConvectiveFields == inarray.size(),
99  "Number of convective fields and Inarray are not compatible");
100 
101  // use dimension of Velocity vector to dictate dimension of operation
102  int ndim = advVel.size();
103  Array<OneD, Array<OneD, NekDouble>> AdvVel(advVel.size());
104 
105  Array<OneD, Array<OneD, NekDouble>> velocity(ndim);
106 
107  LibUtilities::Timer timer;
108 
109  for (int i = 0; i < ndim; ++i)
110  {
111  if (fields[i]->GetWaveSpace() && !m_SingleMode && !m_HalfMode &&
113  {
114  velocity[i] = Array<OneD, NekDouble>(nqtot, 0.0);
115  fields[i]->HomogeneousBwdTrans(nqtot, advVel[i], velocity[i]);
116  }
117  else
118  {
119  velocity[i] = advVel[i];
120  }
121  }
122 
123  int nPointsTot = fields[0]->GetNpoints();
124  Array<OneD, NekDouble> grad0, grad1, grad2, wkSp;
125 
126  NekDouble OneDptscale = 1.5; // factor to rescale 1d points in dealiasing
127 
129  {
130  // Get number of points to dealias a quadratic non-linearity
131  nPointsTot = fields[0]->Get1DScaledTotPoints(OneDptscale);
132  }
133 
134  // interpolate Advection velocity
135  if (m_specHP_dealiasing) // interpolate advection field to higher space.
136  {
137  for (int i = 0; i < ndim; ++i)
138  {
139  AdvVel[i] = Array<OneD, NekDouble>(nPointsTot);
140  // interpolate infield to 3/2 dimension
141  timer.Start();
142  fields[0]->PhysInterp1DScaled(OneDptscale, velocity[i], AdvVel[i]);
143  timer.Stop();
144  timer.AccumulateRegion("Interp1DScaled");
145  }
146  }
147  else
148  {
149  for (int i = 0; i < ndim; ++i)
150  {
151  AdvVel[i] = velocity[i];
152  }
153  }
154 
155  wkSp = Array<OneD, NekDouble>(nPointsTot);
156 
157  // Evaluate V\cdot Grad(u)
158  switch (ndim)
159  {
160  case 1:
161  grad0 = Array<OneD, NekDouble>(fields[0]->GetNpoints());
162  for (int n = 0; n < nConvectiveFields; ++n)
163  {
164  fields[0]->PhysDeriv(inarray[n], grad0);
165  if (m_specHP_dealiasing) // interpolate gradient field
166  {
167  Array<OneD, NekDouble> Outarray(nPointsTot);
168  fields[0]->PhysInterp1DScaled(OneDptscale, grad0, wkSp);
169  Vmath::Vmul(nPointsTot, wkSp, 1, AdvVel[0], 1, Outarray, 1);
170  // Galerkin project solution back to origianl spac
171  timer.Start();
172  fields[0]->PhysGalerkinProjection1DScaled(
173  OneDptscale, Outarray, outarray[n]);
174  timer.Stop();
175  timer.AccumulateRegion("GalerinProject");
176  }
177  else
178  {
179  Vmath::Vmul(nPointsTot, grad0, 1, AdvVel[0], 1, outarray[n],
180  1);
181  }
182  }
183  break;
184  case 2:
185  grad0 = Array<OneD, NekDouble>(nqtot);
186  grad1 = Array<OneD, NekDouble>(nqtot);
187  for (int n = 0; n < nConvectiveFields; ++n)
188  {
189  fields[0]->PhysDeriv(inarray[n], grad0, grad1);
190 
191  if (m_specHP_dealiasing) // interpolate gradient field
192  {
193  Array<OneD, NekDouble> Outarray(nPointsTot);
194  fields[0]->PhysInterp1DScaled(OneDptscale, grad0, wkSp);
195  Vmath::Vmul(nPointsTot, wkSp, 1, AdvVel[0], 1, Outarray, 1);
196  timer.Start();
197  fields[0]->PhysInterp1DScaled(OneDptscale, grad1, wkSp);
198  timer.Stop();
199  timer.AccumulateRegion("Interp1DScaled");
200  Vmath::Vvtvp(nPointsTot, wkSp, 1, AdvVel[1], 1, Outarray, 1,
201  Outarray, 1);
202  // Galerkin project solution back to original space
203  timer.Start();
204  fields[0]->PhysGalerkinProjection1DScaled(
205  OneDptscale, Outarray, outarray[n]);
206  timer.Stop();
207  timer.AccumulateRegion("GalerinProject");
208  }
209  else
210  {
211  Vmath::Vmul(nPointsTot, grad0, 1, AdvVel[0], 1, outarray[n],
212  1);
213  Vmath::Vvtvp(nPointsTot, grad1, 1, AdvVel[1], 1,
214  outarray[n], 1, outarray[n], 1);
215  }
216  }
217  break;
218  case 3:
219  if (m_homogen_dealiasing == true && m_specHP_dealiasing == true)
220  {
221  Array<OneD, Array<OneD, NekDouble>> grad(ndim);
222  Array<OneD, Array<OneD, NekDouble>> gradScaled(
223  ndim * nConvectiveFields);
224  Array<OneD, Array<OneD, NekDouble>> Outarray(nConvectiveFields);
225  for (int i = 0; i < ndim; i++)
226  {
227  grad[i] = Array<OneD, NekDouble>(nqtot);
228  }
229  for (int i = 0; i < ndim * nConvectiveFields; i++)
230  {
231  gradScaled[i] = Array<OneD, NekDouble>(nPointsTot);
232  }
233  for (int i = 0; i < nConvectiveFields; i++)
234  {
235  Outarray[i] = Array<OneD, NekDouble>(nPointsTot);
236  }
237 
238  for (int n = 0; n < nConvectiveFields; n++)
239  {
240  fields[0]->PhysDeriv(inarray[n], grad[0], grad[1], grad[2]);
241  for (int i = 0; i < ndim; i++)
242  {
243  timer.Start();
244  fields[0]->PhysInterp1DScaled(OneDptscale, grad[i],
245  gradScaled[n * ndim + i]);
246  timer.Stop();
247  timer.AccumulateRegion("Interp1DScaled");
248  }
249  }
250 
251  fields[0]->DealiasedDotProd(nPointsTot, AdvVel, gradScaled,
252  Outarray);
253 
254  timer.Start();
255  for (int n = 0; n < nConvectiveFields; n++)
256  {
257  fields[0]->PhysGalerkinProjection1DScaled(
258  OneDptscale, Outarray[n], outarray[n]);
259  }
260  timer.Stop();
261  timer.AccumulateRegion("GalerinProject");
262  }
263  else if (m_homogen_dealiasing == true &&
264  m_specHP_dealiasing == false)
265  {
266  Array<OneD, Array<OneD, NekDouble>> grad(ndim *
267  nConvectiveFields);
268  Array<OneD, Array<OneD, NekDouble>> Outarray(nConvectiveFields);
269  for (int i = 0; i < ndim * nConvectiveFields; i++)
270  {
271  grad[i] = Array<OneD, NekDouble>(nPointsTot);
272  }
273  for (int i = 0; i < nConvectiveFields; i++)
274  {
275  Outarray[i] = Array<OneD, NekDouble>(nPointsTot);
276  }
277 
278  for (int n = 0; n < nConvectiveFields; n++)
279  {
280  fields[0]->PhysDeriv(inarray[n], grad[n * ndim + 0],
281  grad[n * ndim + 1],
282  grad[n * ndim + 2]);
283  }
284 
285  fields[0]->DealiasedDotProd(nPointsTot, AdvVel, grad, outarray);
286  }
287  else
288  {
289  grad0 = Array<OneD, NekDouble>(nqtot);
290  grad1 = Array<OneD, NekDouble>(nqtot);
291  grad2 = Array<OneD, NekDouble>(nqtot);
292  Array<OneD, NekDouble> tmp = grad2;
293  for (int n = 0; n < nConvectiveFields; ++n)
294  {
295  if (fields[0]->GetWaveSpace() == true &&
296  fields[0]->GetExpType() == MultiRegions::e3DH1D)
297  {
298  fields[0]->HomogeneousBwdTrans(nqtot, inarray[n], tmp);
299  fields[0]->PhysDeriv(tmp, grad0, grad1);
300  // Take d/dz derivative using wave space field
301  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],
302  inarray[n], outarray[n]);
303  fields[0]->HomogeneousBwdTrans(nqtot, outarray[n],
304  grad2);
305  }
306  else if (fields[0]->GetWaveSpace() == true &&
307  fields[0]->GetExpType() == MultiRegions::e3DH2D)
308  {
309  fields[0]->HomogeneousBwdTrans(nqtot, inarray[n], tmp);
310  fields[0]->PhysDeriv(tmp, grad0);
311  // Take d/dy derivative using wave space field
312  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],
313  inarray[n], outarray[n]);
314  fields[0]->HomogeneousBwdTrans(nqtot, outarray[n],
315  grad1);
316  // Take d/dz derivative using wave space field
317  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],
318  inarray[n], outarray[n]);
319  fields[0]->HomogeneousBwdTrans(nqtot, outarray[n],
320  grad2);
321  }
322  else
323  {
324  fields[0]->PhysDeriv(inarray[n], grad0, grad1, grad2);
325  }
326  if (m_specHP_dealiasing) // interpolate spectral/hp gradient
327  // field
328  {
329  Array<OneD, NekDouble> Outarray(nPointsTot);
330  timer.Start();
331  fields[0]->PhysInterp1DScaled(OneDptscale, grad0, wkSp);
332  timer.Stop();
333  timer.AccumulateRegion("Interp1DScaled");
334  Vmath::Vmul(nPointsTot, wkSp, 1, AdvVel[0], 1, Outarray,
335  1);
336 
337  timer.Start();
338  fields[0]->PhysInterp1DScaled(OneDptscale, grad1, wkSp);
339  timer.Stop();
340  timer.AccumulateRegion("Interp1DScaled");
341  Vmath::Vvtvp(nPointsTot, wkSp, 1, AdvVel[1], 1,
342  Outarray, 1, Outarray, 1);
343 
344  timer.Start();
345  fields[0]->PhysInterp1DScaled(OneDptscale, grad2, wkSp);
346  timer.Stop();
347  timer.AccumulateRegion("Interp1DScaled");
348  Vmath::Vvtvp(nPointsTot, wkSp, 1, AdvVel[2], 1,
349  Outarray, 1, Outarray, 1);
350  timer.Start();
351  fields[0]->PhysGalerkinProjection1DScaled(
352  OneDptscale, Outarray, outarray[n]);
353  timer.Stop();
354  timer.AccumulateRegion("GalerinProject");
355  }
356  else
357  {
358  Vmath::Vmul(nPointsTot, grad0, 1, AdvVel[0], 1,
359  outarray[n], 1);
360  Vmath::Vvtvp(nPointsTot, grad1, 1, AdvVel[1], 1,
361  outarray[n], 1, outarray[n], 1);
362  Vmath::Vvtvp(nPointsTot, grad2, 1, AdvVel[2], 1,
363  outarray[n], 1, outarray[n], 1);
364  }
365 
366  if (fields[0]->GetWaveSpace() == true)
367  {
368  fields[0]->HomogeneousFwdTrans(nqtot, outarray[n],
369  outarray[n]);
370  }
371  }
372  }
373  break;
374  default:
375  ASSERTL0(false, "dimension unknown");
376  }
377 
378  for (int n = 0; n < nConvectiveFields; ++n)
379  {
380  Vmath::Neg(nqtot, outarray[n], 1);
381  }
382 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:215
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:249
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:91
double NekDouble
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:209
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.cpp:518
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:574

References Nektar::LibUtilities::Timer::AccumulateRegion(), ASSERTL0, ASSERTL1, Nektar::MultiRegions::DirCartesianMap, Nektar::MultiRegions::e3DH1D, Nektar::MultiRegions::e3DH2D, m_HalfMode, m_homogen_dealiasing, m_SingleMode, m_specHP_dealiasing, Vmath::Neg(), Nektar::LibUtilities::Timer::Start(), Nektar::LibUtilities::Timer::Stop(), Vmath::Vmul(), and Vmath::Vvtvp().

◆ v_InitObject()

void Nektar::NavierStokesAdvection::v_InitObject ( LibUtilities::SessionReaderSharedPtr  pSession,
Array< OneD, MultiRegions::ExpListSharedPtr pFields 
)
overrideprotectedvirtual

Initialises the advection object.

This function should be overridden in derived classes to initialise the specific advection data members. However, this base class function should be called as the first statement of the overridden function to ensure the base class is correctly initialised in order.

Parameters
pSessionSession information.
pFieldsExpansion lists for scalar fields.

Reimplemented from Nektar::SolverUtils::Advection.

Definition at line 72 of file NavierStokesAdvection.cpp.

75 {
76  m_homogen_dealiasing = pSession->DefinesSolverInfo("dealiasing");
77 
78  pSession->MatchSolverInfo("SPECTRALHPDEALIASING", "True",
79  m_specHP_dealiasing, false);
80  pSession->MatchSolverInfo("ModeType", "SingleMode", m_SingleMode, false);
81  pSession->MatchSolverInfo("ModeType", "HalfMode", m_HalfMode, false);
82 
83  Advection::v_InitObject(pSession, pFields);
84 }
virtual SOLVER_UTILS_EXPORT void v_InitObject(LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)
Initialises the advection object.
Definition: Advection.cpp:353

References m_HalfMode, m_homogen_dealiasing, m_SingleMode, m_specHP_dealiasing, and Nektar::SolverUtils::Advection::v_InitObject().

Friends And Related Function Documentation

◆ MemoryManager< NavierStokesAdvection >

friend class MemoryManager< NavierStokesAdvection >
friend

Definition at line 1 of file NavierStokesAdvection.h.

Member Data Documentation

◆ className

string Nektar::NavierStokesAdvection::className
static
Initial value:
=
"Convective", NavierStokesAdvection::create, "Convective")
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:198
static SolverUtils::AdvectionSharedPtr create(std::string)
Creates an instance of this class.
AdvectionFactory & GetAdvectionFactory()
Gets the factory for initialising advection objects.
Definition: Advection.cpp:47

Name of class.

Definition at line 56 of file NavierStokesAdvection.h.

◆ className2

string Nektar::NavierStokesAdvection::className2
static
Initial value:

Definition at line 57 of file NavierStokesAdvection.h.

◆ m_HalfMode

bool Nektar::NavierStokesAdvection::m_HalfMode
private

Definition at line 90 of file NavierStokesAdvection.h.

Referenced by v_Advect(), and v_InitObject().

◆ m_homogen_dealiasing

bool Nektar::NavierStokesAdvection::m_homogen_dealiasing
private

Definition at line 88 of file NavierStokesAdvection.h.

Referenced by v_Advect(), and v_InitObject().

◆ m_SingleMode

bool Nektar::NavierStokesAdvection::m_SingleMode
private

Definition at line 89 of file NavierStokesAdvection.h.

Referenced by v_Advect(), and v_InitObject().

◆ m_specHP_dealiasing

bool Nektar::NavierStokesAdvection::m_specHP_dealiasing
private

Definition at line 87 of file NavierStokesAdvection.h.

Referenced by SetSpecHPDealiasing(), v_Advect(), and v_InitObject().

◆ navierStokesAdvectionTypeLookupIds

std::string Nektar::NavierStokesAdvection::navierStokesAdvectionTypeLookupIds
staticprotected
Initial value:
= {
"True", 0),
"False", 1)}
static std::string RegisterEnumValue(std::string pEnum, std::string pString, int pEnumValue)
Registers an enumeration value.

Definition at line 84 of file NavierStokesAdvection.h.