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

#include <SkewSymmetricAdvection.h>

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

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

 SkewSymmetricAdvection ()
 
virtual ~SkewSymmetricAdvection ()
 
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_homogen_dealiasing
 
bool m_SingleMode
 
bool m_HalfMode
 

Friends

class MemoryManager< SkewSymmetricAdvection >
 

Additional Inherited Members

- 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...
 
- 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 SkewSymmetricAdvection.h.

Constructor & Destructor Documentation

Nektar::SkewSymmetricAdvection::SkewSymmetricAdvection ( )
protected

Definition at line 50 of file SkewSymmetricAdvection.cpp.

50  :
51  Advection()
52 
53 {
54 }
Nektar::SkewSymmetricAdvection::~SkewSymmetricAdvection ( )
protectedvirtual

Definition at line 60 of file SkewSymmetricAdvection.cpp.

61 {
62 }

Member Function Documentation

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

Creates an instance of this class.

Definition at line 52 of file SkewSymmetricAdvection.h.

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

52  {
54  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
void Nektar::SkewSymmetricAdvection::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 84 of file SkewSymmetricAdvection.cpp.

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

93 {
94  // use dimension of Velocity vector to dictate dimension of operation
95  int ndim = advVel.num_elements();
96  int nqtot = fields[0]->GetTotPoints();
97  ASSERTL1(nConvectiveFields == inarray.num_elements(),"Number of convective fields and Inarray are not compatible");
98 
99  Array<OneD, Array<OneD, NekDouble> > velocity(ndim);
100  for(int i = 0; i < ndim; ++i)
101  {
102  if(fields[i]->GetWaveSpace() && !m_SingleMode && !m_HalfMode)
103  {
104  velocity[i] = Array<OneD, NekDouble>(nqtot,0.0);
105  fields[i]->HomogeneousBwdTrans(advVel[i],velocity[i]);
106  }
107  else
108  {
109  velocity[i] = advVel[i];
110  }
111  }
112 
113  for(int n = 0; n < nConvectiveFields; ++n)
114  {
115  // ToDo: here we should add a check that V has right dimension
116 
117  int nPointsTot = fields[0]->GetNpoints();
118  Array<OneD, NekDouble> gradV0,gradV1,gradV2, tmp, Up;
119 
120  gradV0 = Array<OneD, NekDouble> (nPointsTot);
121  tmp = Array<OneD, NekDouble> (nPointsTot);
122 
123  // Evaluate V\cdot Grad(u)
124  switch(ndim)
125  {
126  case 1:
127  fields[0]->PhysDeriv(inarray[n],gradV0);
128  Vmath::Vmul(nPointsTot,gradV0,1,velocity[0],1,outarray[n],1);
129  Vmath::Vmul(nPointsTot,inarray[n],1,velocity[0],1,gradV0,1);
130  fields[0]->PhysDeriv(gradV0,tmp);
131  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
132  Vmath::Smul(nPointsTot,0.5,outarray[n],1,outarray[n],1);
133  break;
134  case 2:
135  gradV1 = Array<OneD, NekDouble> (nPointsTot);
136  fields[0]->PhysDeriv(inarray[n],gradV0,gradV1);
137  Vmath::Vmul (nPointsTot,gradV0,1,velocity[0],1,outarray[n],1);
138  Vmath::Vvtvp(nPointsTot,gradV1,1,velocity[1],1,outarray[n],1,outarray[n],1);
139  Vmath::Vmul(nPointsTot,inarray[n],1,velocity[0],1,gradV0,1);
140  Vmath::Vmul(nPointsTot,inarray[n],1,velocity[1],1,gradV1,1);
141  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],gradV0,tmp);
142  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
143  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],gradV1,tmp);
144  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
145  Vmath::Smul(nPointsTot,0.5,outarray[n],1,outarray[n],1);
146  break;
147  case 3:
148  gradV1 = Array<OneD, NekDouble> (nPointsTot);
149  gradV2 = Array<OneD, NekDouble> (nPointsTot);
150 
151  fields[0]->PhysDeriv(inarray[n],gradV0,gradV1,gradV2);
152 
153  //outarray[n] = 1/2(u*du/dx + v*du/dy + w*du/dz + duu/dx + duv/dy + duw/dz)
154 
155  if(m_homogen_dealiasing == true && fields[0]->GetWaveSpace() == false)
156  {
157  fields[0]->DealiasedProd(velocity[0],gradV0,gradV0,m_CoeffState);
158  fields[0]->DealiasedProd(velocity[1],gradV1,gradV1,m_CoeffState);
159  fields[0]->DealiasedProd(velocity[2],gradV2,gradV2,m_CoeffState);
160  Vmath::Vadd(nPointsTot,gradV0,1,gradV1,1,outarray[n],1);
161  Vmath::Vadd(nPointsTot,gradV2,1,outarray[n],1,outarray[n],1);
162  fields[0]->DealiasedProd(inarray[n],velocity[0],gradV0,m_CoeffState);
163  fields[0]->DealiasedProd(inarray[n],velocity[1],gradV1,m_CoeffState);
164  fields[0]->DealiasedProd(inarray[n],velocity[2],gradV2,m_CoeffState);
165  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],gradV0,tmp);
166  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
167  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],gradV1,tmp);
168  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
169  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],gradV2,tmp);
170  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
171  Vmath::Smul(nPointsTot,0.5,outarray[n],1,outarray[n],1);
172  }
173  else if(fields[0]->GetWaveSpace() == true && m_homogen_dealiasing == false)
174  {
175  Up = Array<OneD, NekDouble> (nPointsTot);
176  //vector reused to avoid even more memory requirements
177  //names may be misleading
178  fields[0]->HomogeneousBwdTrans(gradV0,tmp);
179  Vmath::Vmul(nPointsTot,tmp,1,velocity[0],1,outarray[n],1); // + u*du/dx
180  fields[0]->HomogeneousBwdTrans(gradV1,tmp);
181  Vmath::Vvtvp(nPointsTot,tmp,1,velocity[1],1,outarray[n],1,outarray[n],1);// + v*du/dy
182  fields[0]->HomogeneousBwdTrans(gradV2,tmp);
183  Vmath::Vvtvp(nPointsTot,tmp,1,velocity[2],1,outarray[n],1,outarray[n],1);// + w*du/dz
184 
185  fields[0]->HomogeneousBwdTrans(inarray[n],Up);
186  Vmath::Vmul(nPointsTot,Up,1,velocity[0],1,gradV0,1);
187  Vmath::Vmul(nPointsTot,Up,1,velocity[1],1,gradV1,1);
188  Vmath::Vmul(nPointsTot,Up,1,velocity[2],1,gradV2,1);
189 
190  fields[0]->SetWaveSpace(false);
191  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],gradV0,tmp);//duu/dx
192  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
193  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],gradV1,tmp);//duv/dy
194  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
195  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],gradV2,tmp);//duw/dz
196  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
197  fields[0]->SetWaveSpace(true);
198 
199  Vmath::Smul(nPointsTot,0.5,outarray[n],1,tmp,1);
200  fields[0]->HomogeneousFwdTrans(tmp,outarray[n]);
201  }
202  else if(fields[0]->GetWaveSpace() == false && m_homogen_dealiasing == false)
203  {
204  Vmath::Vmul(nPointsTot,gradV0,1,velocity[0],1,outarray[n],1);
205  Vmath::Vvtvp(nPointsTot,gradV1,1,velocity[1],1,outarray[n],1,outarray[n],1);
206  Vmath::Vvtvp(nPointsTot,gradV2,1,velocity[2],1,outarray[n],1,outarray[n],1);
207  Vmath::Vmul(nPointsTot,inarray[n],1,velocity[0],1,gradV0,1);
208  Vmath::Vmul(nPointsTot,inarray[n],1,velocity[1],1,gradV1,1);
209  Vmath::Vmul(nPointsTot,inarray[n],1,velocity[2],1,gradV2,1);
210  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],gradV0,tmp);
211  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
212  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],gradV1,tmp);
213  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
214  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],gradV2,tmp);
215  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
216  Vmath::Smul(nPointsTot,0.5,outarray[n],1,outarray[n],1);
217  }
218  else
219  {
220  ASSERTL0(false, "Dealiasing is not allowed in combination "
221  "with the Skew-Symmetric advection form for "
222  "efficiency reasons.");
223  }
224  break;
225  default:
226  ASSERTL0(false,"dimension unknown");
227  }
228 
229  Vmath::Neg(nqtot,outarray[n],1);
230  }
231 
232 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
MultiRegions::CoeffState m_CoeffState
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 Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
Definition: Vmath.cpp:213
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.cpp:396
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 Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.cpp:299
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::SkewSymmetricAdvection::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 68 of file SkewSymmetricAdvection.cpp.

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

71 {
72  Advection::v_InitObject(pSession, pFields);
73 
75  m_homogen_dealiasing = pSession->DefinesSolverInfo("dealiasing");
76  pSession->MatchSolverInfo("ModeType","SingleMode",m_SingleMode,false);
77  pSession->MatchSolverInfo("ModeType","HalfMode",m_HalfMode,false);
78 }
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< SkewSymmetricAdvection >
friend

Definition at line 49 of file SkewSymmetricAdvection.h.

Member Data Documentation

string Nektar::SkewSymmetricAdvection::className
static
Initial value:

Name of class.

Definition at line 57 of file SkewSymmetricAdvection.h.

std::string Nektar::SkewSymmetricAdvection::className2
static

Definition at line 58 of file SkewSymmetricAdvection.h.

MultiRegions::CoeffState Nektar::SkewSymmetricAdvection::m_CoeffState
private

Definition at line 81 of file SkewSymmetricAdvection.h.

Referenced by v_Advect(), and v_InitObject().

bool Nektar::SkewSymmetricAdvection::m_HalfMode
private

Definition at line 84 of file SkewSymmetricAdvection.h.

Referenced by v_Advect(), and v_InitObject().

bool Nektar::SkewSymmetricAdvection::m_homogen_dealiasing
private

Definition at line 82 of file SkewSymmetricAdvection.h.

Referenced by v_Advect(), and v_InitObject().

bool Nektar::SkewSymmetricAdvection::m_SingleMode
private

Definition at line 83 of file SkewSymmetricAdvection.h.

Referenced by v_Advect(), and v_InitObject().