Nektar++
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:
[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 ()
 
 ~SkewSymmetricAdvection () 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

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

Constructor & Destructor Documentation

◆ SkewSymmetricAdvection()

Nektar::SkewSymmetricAdvection::SkewSymmetricAdvection ( )
protected

Definition at line 46 of file SkewSymmetricAdvection.cpp.

46 : Advection()
47
48{
49}

◆ ~SkewSymmetricAdvection()

Nektar::SkewSymmetricAdvection::~SkewSymmetricAdvection ( )
overrideprotecteddefault

Member Function Documentation

◆ create()

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

Creates an instance of this class.

Definition at line 50 of file SkewSymmetricAdvection.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::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 
)
overrideprotectedvirtual

Advects a vector field.

Implements Nektar::SolverUtils::Advection.

Definition at line 68 of file SkewSymmetricAdvection.cpp.

77{
78 // use dimension of Velocity vector to dictate dimension of operation
79 int ndim = advVel.size();
80 int nqtot = fields[0]->GetTotPoints();
81 ASSERTL1(nConvectiveFields == inarray.size(),
82 "Number of convective fields and Inarray are not compatible");
83
84 Array<OneD, Array<OneD, NekDouble>> velocity(ndim);
85 for (int i = 0; i < ndim; ++i)
86 {
87 if (fields[i]->GetWaveSpace() && !m_SingleMode && !m_HalfMode)
88 {
89 velocity[i] = Array<OneD, NekDouble>(nqtot, 0.0);
90 fields[i]->HomogeneousBwdTrans(nqtot, advVel[i], velocity[i]);
91 }
92 else
93 {
94 velocity[i] = advVel[i];
95 }
96 }
97
98 for (int n = 0; n < nConvectiveFields; ++n)
99 {
100 // ToDo: here we should add a check that V has right dimension
101
102 int nPointsTot = fields[0]->GetNpoints();
103 Array<OneD, NekDouble> gradV0, gradV1, gradV2, tmp, Up;
104
105 gradV0 = Array<OneD, NekDouble>(nPointsTot);
106 tmp = Array<OneD, NekDouble>(nPointsTot);
107
108 // Evaluate V\cdot Grad(u)
109 switch (ndim)
110 {
111 case 1:
112 fields[0]->PhysDeriv(inarray[n], gradV0);
113 Vmath::Vmul(nPointsTot, gradV0, 1, velocity[0], 1, outarray[n],
114 1);
115 Vmath::Vmul(nPointsTot, inarray[n], 1, velocity[0], 1, gradV0,
116 1);
117 fields[0]->PhysDeriv(gradV0, tmp);
118 Vmath::Vadd(nPointsTot, tmp, 1, outarray[n], 1, outarray[n], 1);
119 Vmath::Smul(nPointsTot, 0.5, outarray[n], 1, outarray[n], 1);
120 break;
121 case 2:
122 gradV1 = Array<OneD, NekDouble>(nPointsTot);
123 fields[0]->PhysDeriv(inarray[n], gradV0, gradV1);
124 Vmath::Vmul(nPointsTot, gradV0, 1, velocity[0], 1, outarray[n],
125 1);
126 Vmath::Vvtvp(nPointsTot, gradV1, 1, velocity[1], 1, outarray[n],
127 1, outarray[n], 1);
128 Vmath::Vmul(nPointsTot, inarray[n], 1, velocity[0], 1, gradV0,
129 1);
130 Vmath::Vmul(nPointsTot, inarray[n], 1, velocity[1], 1, gradV1,
131 1);
132 fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0], gradV0,
133 tmp);
134 Vmath::Vadd(nPointsTot, tmp, 1, outarray[n], 1, outarray[n], 1);
135 fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1], gradV1,
136 tmp);
137 Vmath::Vadd(nPointsTot, tmp, 1, outarray[n], 1, outarray[n], 1);
138 Vmath::Smul(nPointsTot, 0.5, outarray[n], 1, outarray[n], 1);
139 break;
140 case 3:
141 gradV1 = Array<OneD, NekDouble>(nPointsTot);
142 gradV2 = Array<OneD, NekDouble>(nPointsTot);
143
144 fields[0]->PhysDeriv(inarray[n], gradV0, gradV1, gradV2);
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 (m_homogen_dealiasing == true &&
150 fields[0]->GetWaveSpace() == false)
151 {
152 fields[0]->DealiasedProd(nPointsTot, velocity[0], gradV0,
153 gradV0);
154 fields[0]->DealiasedProd(nPointsTot, velocity[1], gradV1,
155 gradV1);
156 fields[0]->DealiasedProd(nPointsTot, velocity[2], gradV2,
157 gradV2);
158 Vmath::Vadd(nPointsTot, gradV0, 1, gradV1, 1, outarray[n],
159 1);
160 Vmath::Vadd(nPointsTot, gradV2, 1, outarray[n], 1,
161 outarray[n], 1);
162 fields[0]->DealiasedProd(nPointsTot, inarray[n],
163 velocity[0], gradV0);
164 fields[0]->DealiasedProd(nPointsTot, inarray[n],
165 velocity[1], gradV1);
166 fields[0]->DealiasedProd(nPointsTot, inarray[n],
167 velocity[2], gradV2);
168 fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],
169 gradV0, tmp);
170 Vmath::Vadd(nPointsTot, tmp, 1, outarray[n], 1, outarray[n],
171 1);
172 fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],
173 gradV1, tmp);
174 Vmath::Vadd(nPointsTot, tmp, 1, outarray[n], 1, outarray[n],
175 1);
176 fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],
177 gradV2, tmp);
178 Vmath::Vadd(nPointsTot, tmp, 1, outarray[n], 1, outarray[n],
179 1);
180 Vmath::Smul(nPointsTot, 0.5, outarray[n], 1, outarray[n],
181 1);
182 }
183 else if (fields[0]->GetWaveSpace() == true &&
184 m_homogen_dealiasing == false)
185 {
186 Up = Array<OneD, NekDouble>(nPointsTot);
187 // vector reused to avoid even more memory requirements
188 // names may be misleading
189 fields[0]->HomogeneousBwdTrans(nPointsTot, gradV0, tmp);
190 Vmath::Vmul(nPointsTot, tmp, 1, velocity[0], 1, outarray[n],
191 1); // + u*du/dx
192 fields[0]->HomogeneousBwdTrans(nPointsTot, gradV1, tmp);
193 Vmath::Vvtvp(nPointsTot, tmp, 1, velocity[1], 1,
194 outarray[n], 1, outarray[n], 1); // + v*du/dy
195 fields[0]->HomogeneousBwdTrans(nPointsTot, gradV2, tmp);
196 Vmath::Vvtvp(nPointsTot, tmp, 1, velocity[2], 1,
197 outarray[n], 1, outarray[n], 1); // + w*du/dz
198
199 fields[0]->HomogeneousBwdTrans(nPointsTot, inarray[n], Up);
200 Vmath::Vmul(nPointsTot, Up, 1, velocity[0], 1, gradV0, 1);
201 Vmath::Vmul(nPointsTot, Up, 1, velocity[1], 1, gradV1, 1);
202 Vmath::Vmul(nPointsTot, Up, 1, velocity[2], 1, gradV2, 1);
203
204 fields[0]->SetWaveSpace(false);
205 fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],
206 gradV0, tmp); // duu/dx
207 Vmath::Vadd(nPointsTot, tmp, 1, outarray[n], 1, outarray[n],
208 1);
209 fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],
210 gradV1, tmp); // duv/dy
211 Vmath::Vadd(nPointsTot, tmp, 1, outarray[n], 1, outarray[n],
212 1);
213 fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],
214 gradV2, tmp); // duw/dz
215 Vmath::Vadd(nPointsTot, tmp, 1, outarray[n], 1, outarray[n],
216 1);
217 fields[0]->SetWaveSpace(true);
218
219 Vmath::Smul(nPointsTot, 0.5, outarray[n], 1, tmp, 1);
220 fields[0]->HomogeneousFwdTrans(nPointsTot, tmp,
221 outarray[n]);
222 }
223 else if (fields[0]->GetWaveSpace() == false &&
224 m_homogen_dealiasing == false)
225 {
226 Vmath::Vmul(nPointsTot, gradV0, 1, velocity[0], 1,
227 outarray[n], 1);
228 Vmath::Vvtvp(nPointsTot, gradV1, 1, velocity[1], 1,
229 outarray[n], 1, outarray[n], 1);
230 Vmath::Vvtvp(nPointsTot, gradV2, 1, velocity[2], 1,
231 outarray[n], 1, outarray[n], 1);
232 Vmath::Vmul(nPointsTot, inarray[n], 1, velocity[0], 1,
233 gradV0, 1);
234 Vmath::Vmul(nPointsTot, inarray[n], 1, velocity[1], 1,
235 gradV1, 1);
236 Vmath::Vmul(nPointsTot, inarray[n], 1, velocity[2], 1,
237 gradV2, 1);
238 fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],
239 gradV0, tmp);
240 Vmath::Vadd(nPointsTot, tmp, 1, outarray[n], 1, outarray[n],
241 1);
242 fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],
243 gradV1, tmp);
244 Vmath::Vadd(nPointsTot, tmp, 1, outarray[n], 1, outarray[n],
245 1);
246 fields[0]->PhysDeriv(MultiRegions::DirCartesianMap[2],
247 gradV2, tmp);
248 Vmath::Vadd(nPointsTot, tmp, 1, outarray[n], 1, outarray[n],
249 1);
250 Vmath::Smul(nPointsTot, 0.5, outarray[n], 1, outarray[n],
251 1);
252 }
253 else
254 {
255 ASSERTL0(false,
256 "Dealiasing is not allowed in combination "
257 "with the Skew-Symmetric advection form for "
258 "efficiency reasons.");
259 }
260 break;
261 default:
262 ASSERTL0(false, "dimension unknown");
263 }
264
265 Vmath::Neg(nqtot, outarray[n], 1);
266 }
267}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:242
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 Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.hpp:292
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, ASSERTL1, Nektar::MultiRegions::DirCartesianMap, m_HalfMode, m_homogen_dealiasing, m_SingleMode, Vmath::Neg(), Vmath::Smul(), Vmath::Vadd(), Vmath::Vmul(), and Vmath::Vvtvp().

◆ v_InitObject()

void Nektar::SkewSymmetricAdvection::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 54 of file SkewSymmetricAdvection.cpp.

57{
58 Advection::v_InitObject(pSession, pFields);
59
60 m_homogen_dealiasing = pSession->DefinesSolverInfo("dealiasing");
61 pSession->MatchSolverInfo("ModeType", "SingleMode", m_SingleMode, false);
62 pSession->MatchSolverInfo("ModeType", "HalfMode", m_HalfMode, false);
63}
virtual SOLVER_UTILS_EXPORT void v_InitObject(LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)
Initialises the advection object.
Definition: Advection.cpp:264

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

Friends And Related Function Documentation

◆ MemoryManager< SkewSymmetricAdvection >

friend class MemoryManager< SkewSymmetricAdvection >
friend

Definition at line 1 of file SkewSymmetricAdvection.h.

Member Data Documentation

◆ className

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

Name of class.

Definition at line 56 of file SkewSymmetricAdvection.h.

◆ className2

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

Definition at line 57 of file SkewSymmetricAdvection.h.

◆ m_HalfMode

bool Nektar::SkewSymmetricAdvection::m_HalfMode
private

Definition at line 82 of file SkewSymmetricAdvection.h.

Referenced by v_Advect(), and v_InitObject().

◆ m_homogen_dealiasing

bool Nektar::SkewSymmetricAdvection::m_homogen_dealiasing
private

Definition at line 80 of file SkewSymmetricAdvection.h.

Referenced by v_Advect(), and v_InitObject().

◆ m_SingleMode

bool Nektar::SkewSymmetricAdvection::m_SingleMode
private

Definition at line 81 of file SkewSymmetricAdvection.h.

Referenced by v_Advect(), and v_InitObject().