Nektar++
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:
[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 ()
 
 ~AlternateSkewAdvection () override=default
 
void v_InitObject (LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields) override
 Initialises the advection object. More...
 
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 ~Advection ()=default
 
virtual SOLVER_UTILS_EXPORT void v_InitObject (LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)
 Initialises the advection object. More...
 
virtual SOLVER_UTILS_EXPORT 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)=0
 Advects a vector field. More...
 
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...
 
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)
 
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)
 
- 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 AlternateSkewAdvection.h.

Constructor & Destructor Documentation

◆ AlternateSkewAdvection()

Nektar::AlternateSkewAdvection::AlternateSkewAdvection ( )
protected

Constructor. Creates ...

Parameters

param

Definition at line 51 of file AlternateSkewAdvection.cpp.

51 : Advection()
52{
53}

◆ ~AlternateSkewAdvection()

Nektar::AlternateSkewAdvection::~AlternateSkewAdvection ( )
overrideprotecteddefault

Member Function Documentation

◆ create()

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

Creates an instance of this class.

Definition at line 50 of file AlternateSkewAdvection.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().

◆ v_Advect()

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

Advects a vector field.

Implements Nektar::SolverUtils::Advection.

Definition at line 63 of file AlternateSkewAdvection.cpp.

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

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

◆ v_InitObject()

void Nektar::AlternateSkewAdvection::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 55 of file AlternateSkewAdvection.cpp.

58{
59 pSession->MatchSolverInfo("ModeType", "SingleMode", m_SingleMode, false);
60 pSession->MatchSolverInfo("ModeType", "HalfMode", m_HalfMode, false);
61}

References m_HalfMode, and m_SingleMode.

Friends And Related Function Documentation

◆ MemoryManager< AlternateSkewAdvection >

friend class MemoryManager< AlternateSkewAdvection >
friend

Definition at line 1 of file AlternateSkewAdvection.h.

Member Data Documentation

◆ className

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

Name of class.

Definition at line 56 of file AlternateSkewAdvection.h.

◆ className2

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

Definition at line 57 of file AlternateSkewAdvection.h.

◆ m_advectioncalls

int Nektar::AlternateSkewAdvection::m_advectioncalls
private

Definition at line 80 of file AlternateSkewAdvection.h.

Referenced by v_Advect().

◆ m_HalfMode

bool Nektar::AlternateSkewAdvection::m_HalfMode
private

Definition at line 82 of file AlternateSkewAdvection.h.

Referenced by v_Advect(), and v_InitObject().

◆ m_SingleMode

bool Nektar::AlternateSkewAdvection::m_SingleMode
private

Definition at line 81 of file AlternateSkewAdvection.h.

Referenced by v_Advect(), and v_InitObject().