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

Advects a vector field.

Implements Nektar::SolverUtils::Advection.

Definition at line 86 of file NavierStokesAdvection.cpp.

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

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.