Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Public Member Functions | Static Public Member Functions | Static Public Attributes | Protected Member Functions | Private Attributes | Friends | List of all members
Nektar::NavierStokesAdvection Class Reference

#include <NavierStokesAdvection.h>

Inheritance diagram for Nektar::NavierStokesAdvection:
Inheritance graph
[legend]
Collaboration diagram for Nektar::NavierStokesAdvection:
Collaboration graph
[legend]

Public Member Functions

void SetSpecHPDealiasing (bool value)
 
- Public Member Functions inherited from Nektar::SolverUtils::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...
 
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...
 

Static Public Member Functions

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

Static Public Attributes

static std::string className = SolverUtils::GetAdvectionFactory().RegisterCreatorFunction("Convective", NavierStokesAdvection::create)
 Name of class. More...
 
static std::string className2 = SolverUtils::GetAdvectionFactory().RegisterCreatorFunction("NonConservative", NavierStokesAdvection::create)
 

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_SetBaseFlow (const Array< OneD, Array< OneD, NekDouble > > &inarray, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields)
 Overrides the base flow used during linearised advection. More...
 

Private Attributes

MultiRegions::CoeffState m_CoeffState
 
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 45 of file NavierStokesAdvection.h.

Constructor & Destructor Documentation

Nektar::NavierStokesAdvection::NavierStokesAdvection ( )
protected

Constructor. Creates ...

Parameters
param

Definition at line 52 of file NavierStokesAdvection.cpp.

52  :
53  Advection()
54 
55  {
56 
57  }
Nektar::NavierStokesAdvection::~NavierStokesAdvection ( )
protectedvirtual

Definition at line 59 of file NavierStokesAdvection.cpp.

60  {
61  }

Member Function Documentation

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

Creates an instance of this class.

Definition at line 52 of file NavierStokesAdvection.h.

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

52  {
54  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
void Nektar::NavierStokesAdvection::SetSpecHPDealiasing ( bool  value)
inline

Definition at line 60 of file NavierStokesAdvection.h.

References m_specHP_dealiasing.

61  {
62  m_specHP_dealiasing = value;
63  }
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 83 of file NavierStokesAdvection.cpp.

References ASSERTL0, ASSERTL1, Nektar::MultiRegions::DirCartesianMap, Nektar::MultiRegions::e3DH1D, Nektar::MultiRegions::e3DH2D, m_CoeffState, m_HalfMode, m_homogen_dealiasing, m_SingleMode, m_specHP_dealiasing, Vmath::Neg(), Vmath::Vmul(), and Vmath::Vvtvp().

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

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

67  {
69  m_homogen_dealiasing = pSession->DefinesSolverInfo("dealiasing");
70 
71  pSession->MatchSolverInfo("SPECTRALHPDEALIASING","True",m_specHP_dealiasing,false);
72  if(m_specHP_dealiasing == false)
73  {
74  pSession->MatchSolverInfo("SPECTRALHPDEALIASING","On",m_specHP_dealiasing,false);
75  }
76  pSession->MatchSolverInfo("ModeType","SingleMode",m_SingleMode,false);
77  pSession->MatchSolverInfo("ModeType","HalfMode",m_HalfMode,false);
78 
79  Advection::v_InitObject(pSession, pFields);
80  }
Local coefficients.
MultiRegions::CoeffState m_CoeffState
virtual SOLVER_UTILS_EXPORT void v_InitObject(LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)
Initialises the advection object.
Definition: Advection.cpp:100

Friends And Related Function Documentation

friend class MemoryManager< NavierStokesAdvection >
friend

Definition at line 49 of file NavierStokesAdvection.h.

Member Data Documentation

string Nektar::NavierStokesAdvection::className = SolverUtils::GetAdvectionFactory().RegisterCreatorFunction("Convective", NavierStokesAdvection::create)
static

Name of class.

Definition at line 57 of file NavierStokesAdvection.h.

string Nektar::NavierStokesAdvection::className2 = SolverUtils::GetAdvectionFactory().RegisterCreatorFunction("NonConservative", NavierStokesAdvection::create)
static

Definition at line 58 of file NavierStokesAdvection.h.

MultiRegions::CoeffState Nektar::NavierStokesAdvection::m_CoeffState
private

Definition at line 86 of file NavierStokesAdvection.h.

Referenced by v_Advect(), and v_InitObject().

bool Nektar::NavierStokesAdvection::m_HalfMode
private

Definition at line 90 of file NavierStokesAdvection.h.

Referenced by v_Advect(), and v_InitObject().

bool Nektar::NavierStokesAdvection::m_homogen_dealiasing
private

Definition at line 88 of file NavierStokesAdvection.h.

Referenced by v_Advect(), and v_InitObject().

bool Nektar::NavierStokesAdvection::m_SingleMode
private

Definition at line 89 of file NavierStokesAdvection.h.

Referenced by v_Advect(), and v_InitObject().

bool Nektar::NavierStokesAdvection::m_specHP_dealiasing
private

Definition at line 87 of file NavierStokesAdvection.h.

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