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::AlternateSkewAdvection Class Reference

#include <AlternateSkewAdvection.h>

Inheritance diagram for Nektar::AlternateSkewAdvection:
Inheritance graph
[legend]
Collaboration diagram for Nektar::AlternateSkewAdvection:
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

 AlternateSkewAdvection ()
 
virtual ~AlternateSkewAdvection ()
 
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

int m_advectioncalls
 
bool m_SingleMode
 
bool m_HalfMode
 

Friends

class MemoryManager< AlternateSkewAdvection >
 

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 44 of file AlternateSkewAdvection.h.

Constructor & Destructor Documentation

Nektar::AlternateSkewAdvection::AlternateSkewAdvection ( )
protected

Constructor. Creates ...

Parameters
param

Definition at line 53 of file AlternateSkewAdvection.cpp.

54  : Advection()
55 {
56 }
Nektar::AlternateSkewAdvection::~AlternateSkewAdvection ( )
protectedvirtual

Definition at line 58 of file AlternateSkewAdvection.cpp.

59 {
60 }

Member Function Documentation

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

Creates an instance of this class.

Definition at line 51 of file AlternateSkewAdvection.h.

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

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

References ASSERTL0, Nektar::MultiRegions::DirCartesianMap, m_advectioncalls, m_HalfMode, m_SingleMode, Vmath::Smul(), Vmath::Vadd(), Vmath::Vmul(), and Vmath::Vvtvp().

79 {
80  // use dimension of Velocity vector to dictate dimension of operation
81  int ndim = advVel.num_elements();
82  int nPointsTot = fields[0]->GetNpoints();
83  Array<OneD, Array<OneD, NekDouble> > velocity(ndim);
84  for(int i = 0; i < ndim; ++i)
85  {
86  if(fields[i]->GetWaveSpace() && !m_SingleMode && !m_HalfMode)
87  {
88  velocity[i] = Array<OneD, NekDouble>(nPointsTot,0.0);
89  fields[i]->HomogeneousBwdTrans(advVel[i],velocity[i]);
90  }
91  else
92  {
93  velocity[i] = advVel[i];
94  }
95  }
96  for(int n = 0; n < nConvectiveFields; ++n)
97  {
98  // ToDo: here we should add a check that V has right dimension
99  Array<OneD, NekDouble> gradV0,gradV1,gradV2, tmp, Up;
100 
101  gradV0 = Array<OneD, NekDouble> (nPointsTot);
102  tmp = Array<OneD, NekDouble> (nPointsTot);
103 
104  // Evaluate V\cdot Grad(u)
105  switch(ndim)
106  {
107  case 1:
108  if(m_advectioncalls % 2 == 0)
109  {
110  fields[0]->PhysDeriv(inarray[n],gradV0);
111  Vmath::Vmul(nPointsTot,gradV0,1,velocity[0],1,outarray[n],1);
112  }
113  else
114  {
115  Vmath::Vmul(nPointsTot,inarray[n],1,velocity[0],1,gradV0,1);
116  fields[0]->PhysDeriv(gradV0,outarray[n]);
117  }
118  Vmath::Smul(nPointsTot,0.5,outarray[n],1,outarray[n],1); //must be mult by 0.5????
119  break;
120  case 2:
121  gradV1 = Array<OneD, NekDouble> (nPointsTot);
122  if(m_advectioncalls % 2 == 0)
123  {
124  fields[0]->PhysDeriv(inarray[n],gradV0,gradV1);
125  Vmath::Vmul (nPointsTot,gradV0,1,velocity[0],1,outarray[n],1);
126  Vmath::Vvtvp(nPointsTot,gradV1,1,velocity[1],1,outarray[n],1,outarray[n],1);
127  }
128  else
129  {
130  Vmath::Vmul(nPointsTot,inarray[n],1,velocity[0],1,gradV0,1);
131  Vmath::Vmul(nPointsTot,inarray[n],1,velocity[1],1,gradV1,1);
132  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],gradV0,outarray[n]);
133  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],gradV1,tmp);
134  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
135  }
136  Vmath::Smul(nPointsTot,1.0,outarray[n],1,outarray[n],1); //must be mult by 0.5????
137  break;
138  case 3:
139  gradV1 = Array<OneD, NekDouble> (nPointsTot);
140  gradV2 = Array<OneD, NekDouble> (nPointsTot);
141 
142  //outarray[n] = 1/2(u*du/dx + v*du/dy + w*du/dz + duu/dx + duv/dy + duw/dz)
143 
144  if(fields[0]->GetWaveSpace() == true)
145  {
146  if(m_advectioncalls % 2 == 0)
147  {
148  //vector reused to avoid even more memory requirements
149  //names may be misleading
150  fields[0]->PhysDeriv(inarray[n],gradV0,gradV1,gradV2);
151  fields[0]->HomogeneousBwdTrans(gradV0,tmp);
152  Vmath::Vmul(nPointsTot,tmp,1,velocity[0],1,outarray[n],1); // + u*du/dx
153  fields[0]->HomogeneousBwdTrans(gradV1,tmp);
154  Vmath::Vvtvp(nPointsTot,tmp,1,velocity[1],1,outarray[n],1,outarray[n],1);// + v*du/dy
155  fields[0]->HomogeneousBwdTrans(gradV2,tmp);
156  Vmath::Vvtvp(nPointsTot,tmp,1,velocity[2],1,outarray[n],1,outarray[n],1);// + w*du/dz
157  }
158  else
159  {
160  Up = Array<OneD, NekDouble> (nPointsTot);
161  fields[0]->HomogeneousBwdTrans(inarray[n],Up);
162  Vmath::Vmul(nPointsTot,Up,1,velocity[0],1,gradV0,1);
163  Vmath::Vmul(nPointsTot,Up,1,velocity[1],1,gradV1,1);
164  Vmath::Vmul(nPointsTot,Up,1,velocity[2],1,gradV2,1);
165 
166  fields[0]->SetWaveSpace(false);
167  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],gradV0,outarray[n]);//duu/dx
168  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],gradV1,tmp);//duv/dy
169  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
170  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],gradV2,tmp);//duw/dz
171  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
172  fields[0]->SetWaveSpace(true);
173  }
174 
175  Vmath::Smul(nPointsTot,1.0,outarray[n],1,tmp,1); //must be mult by 0.5????
176  fields[0]->HomogeneousFwdTrans(tmp,outarray[n]);
177  }
178  else
179  {
180  if(m_advectioncalls % 2 == 0)
181  {
182  fields[0]->PhysDeriv(inarray[n],gradV0,gradV1,gradV2);
183  Vmath::Vmul(nPointsTot,gradV0,1,velocity[0],1,outarray[n],1);
184  Vmath::Vvtvp(nPointsTot,gradV1,1,velocity[1],1,outarray[n],1,outarray[n],1);
185  Vmath::Vvtvp(nPointsTot,gradV2,1,velocity[2],1,outarray[n],1,outarray[n],1);
186  }
187  else
188  {
189  Vmath::Vmul(nPointsTot,inarray[n],1,velocity[0],1,gradV0,1);
190  Vmath::Vmul(nPointsTot,inarray[n],1,velocity[1],1,gradV1,1);
191  Vmath::Vmul(nPointsTot,inarray[n],1,velocity[2],1,gradV2,1);
192  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],gradV0,outarray[n]);
193  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],gradV1,tmp);
194  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
195  fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],gradV2,tmp);
196  Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
197  }
198  Vmath::Smul(nPointsTot,1.0,outarray[n],1,outarray[n],1); //must be mult by 0.5????
199  }
200  break;
201  default:
202  ASSERTL0(false,"dimension unknown");
203  }
204  }
205 }
#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 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
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:86
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::AlternateSkewAdvection::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 62 of file AlternateSkewAdvection.cpp.

References m_HalfMode, and m_SingleMode.

65 {
66  pSession->MatchSolverInfo("ModeType","SingleMode",m_SingleMode,false);
67  pSession->MatchSolverInfo("ModeType","HalfMode",m_HalfMode,false);
68 }

Friends And Related Function Documentation

friend class MemoryManager< AlternateSkewAdvection >
friend

Definition at line 48 of file AlternateSkewAdvection.h.

Member Data Documentation

string Nektar::AlternateSkewAdvection::className
static
Initial value:

Name of class.

Definition at line 56 of file AlternateSkewAdvection.h.

std::string Nektar::AlternateSkewAdvection::className2
static

Definition at line 57 of file AlternateSkewAdvection.h.

int Nektar::AlternateSkewAdvection::m_advectioncalls
private

Definition at line 80 of file AlternateSkewAdvection.h.

Referenced by v_Advect().

bool Nektar::AlternateSkewAdvection::m_HalfMode
private

Definition at line 82 of file AlternateSkewAdvection.h.

Referenced by v_Advect(), and v_InitObject().

bool Nektar::AlternateSkewAdvection::m_SingleMode
private

Definition at line 81 of file AlternateSkewAdvection.h.

Referenced by v_Advect(), and v_InitObject().