Nektar++
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
Nektar::SolverUtils::Advection Class Referenceabstract

An abstract base class encapsulating the concept of advection of a vector field. More...

#include <Advection.h>

Inheritance diagram for Nektar::SolverUtils::Advection:
[legend]

Public Member Functions

virtual SOLVER_UTILS_EXPORT ~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 Member Functions

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

Protected Attributes

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

An abstract base class encapsulating the concept of advection of a vector field.

Subclasses override the Advection::v_InitObject function to initialise the object and the Advection::v_Advect function to evaluate the advection of the vector field.

Definition at line 82 of file Advection.h.

Constructor & Destructor Documentation

◆ ~Advection()

virtual SOLVER_UTILS_EXPORT Nektar::SolverUtils::Advection::~Advection ( )
inlinevirtual

Definition at line 85 of file Advection.h.

85{};

Member Function Documentation

◆ AddTraceJacToMat() [1/2]

template<typename DataType , typename TypeNekBlkMatSharedPtr >
void Nektar::SolverUtils::Advection::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 
)

Definition at line 85 of file Advection.cpp.

92{
93 MultiRegions::ExpListSharedPtr tracelist = pFields[0]->GetTrace();
94 std::shared_ptr<LocalRegions::ExpansionVector> traceExp =
95 tracelist->GetExp();
96 int ntotTrac = (*traceExp).size();
97 int nTracePnts = tracelist->GetTotPoints();
98 int nTracPnt, nTracCoef;
99
100 DNekMatSharedPtr tmp2Add;
101 Array<OneD, DataType> MatData0;
102 Array<OneD, NekDouble> MatData1;
103 Array<OneD, NekDouble> MatData2;
104
105 Array<OneD, DNekMatSharedPtr> TraceJacFwd(ntotTrac);
106 Array<OneD, DNekMatSharedPtr> TraceJacBwd(ntotTrac);
107
108 for (int nelmt = 0; nelmt < ntotTrac; nelmt++)
109 {
110 nTracCoef = (*traceExp)[nelmt]->GetNcoeffs();
111 nTracPnt = (*traceExp)[nelmt]->GetTotPoints();
112 TraceJacFwd[nelmt] =
113 MemoryManager<DNekMat>::AllocateSharedPtr(nTracCoef, nTracPnt, 0.0);
114 TraceJacBwd[nelmt] =
115 MemoryManager<DNekMat>::AllocateSharedPtr(nTracCoef, nTracPnt, 0.0);
116 }
117
118 std::shared_ptr<LocalRegions::ExpansionVector> fieldExp =
119 pFields[0]->GetExp();
120 int nTotElmt = (*fieldExp).size();
121 int nElmtPnt, nElmtCoef;
122
123 Array<OneD, DNekMatSharedPtr> ElmtJacQuad(nTotElmt);
124 Array<OneD, DNekMatSharedPtr> ElmtJacCoef(nTotElmt);
125
126 Array<OneD, NekDouble> SymmMatData;
127
128 for (int nelmt = 0; nelmt < nTotElmt; nelmt++)
129 {
130 nElmtCoef = (*fieldExp)[nelmt]->GetNcoeffs();
131 nElmtPnt = (*fieldExp)[nelmt]->GetTotPoints();
132
133 ElmtJacQuad[nelmt] =
134 MemoryManager<DNekMat>::AllocateSharedPtr(nElmtCoef, nElmtPnt, 0.0);
136 nElmtCoef, nElmtCoef, 0.0);
137 }
138
139 bool TracePntJacGradflag = true;
140
141 Array<OneD, Array<OneD, DataType>> TraceJacConsSign(2);
142 for (int i = 0; i < 2; i++)
143 {
144 TraceJacConsSign[i] = Array<OneD, DataType>(nTracePnts, 1.0);
145 }
146
147 if (0 == TracePntJacGrad.size())
148 {
149 TracePntJacGradflag = false;
150 }
151
152 for (int m = 0; m < nConvectiveFields; m++)
153 {
154 for (int n = 0; n < nConvectiveFields; n++)
155 {
156 // ElmtJacCons to set 0
157 for (int nelmt = 0; nelmt < nTotElmt; nelmt++)
158 {
159 (*ElmtJacCoef[nelmt]) = 0.0;
160 (*ElmtJacQuad[nelmt]) = 0.0;
161 }
162
163 if (TracePntJacGradflag)
164 {
165 // only BJac is stored here because BJac = - FJac
166 for (int ndir = 0; ndir < nSpaceDim; ndir++)
167 {
168 // ElmtJacGrad to set 0
169 for (int nelmt = 0; nelmt < nTotElmt; nelmt++)
170 {
171 (*ElmtJacQuad[nelmt]) = 0.0;
172 }
173 int ngrad = n * nSpaceDim + ndir;
174
175 CalcJacobTraceInteg(pFields, m, ngrad, TracePntJacGrad,
176 TracePntJacGradSign, TraceJacFwd,
177 TraceJacBwd);
178 pFields[0]->AddTraceJacToElmtJac(TraceJacFwd, TraceJacBwd,
179 ElmtJacQuad);
180 // cost can be lower down by doing 3 directions together
181 pFields[0]->AddRightIPTPhysDerivBase(ndir, ElmtJacQuad,
182 ElmtJacCoef);
183 }
184 for (int nelmt = 0; nelmt < nTotElmt; nelmt++)
185 {
186 nElmtCoef = (*fieldExp)[nelmt]->GetNcoeffs();
187 if (SymmMatData.size() < nElmtCoef)
188 {
189 SymmMatData =
190 Array<OneD, NekDouble>(nElmtCoef * nElmtCoef);
191 }
192 MatData1 = ElmtJacCoef[nelmt]->GetPtr();
193 for (int i = 0; i < nElmtCoef; i++)
194 {
195 Vmath::Vcopy(nElmtCoef, &MatData1[i * nElmtCoef], 1,
196 &SymmMatData[i], nElmtCoef);
197 }
198 Vmath::Vadd(nElmtCoef * nElmtCoef, SymmMatData, 1, MatData1,
199 1, MatData1, 1);
200 }
201 }
202
203 CalcJacobTraceInteg(pFields, m, n, TracePntJacCons,
204 TraceJacConsSign, TraceJacFwd, TraceJacBwd);
205
206 pFields[0]->AddTraceJacToElmtJac(TraceJacFwd, TraceJacBwd,
207 ElmtJacQuad);
208 // Memory can be save to explore reusing ElmtJacQuad as output
209 pFields[0]->AddRightIPTBaseMatrix(ElmtJacQuad, ElmtJacCoef);
210
211 // add ElmtJacCons to gmtxarray[m][n]
212 for (int nelmt = 0; nelmt < nTotElmt; nelmt++)
213 {
214 tmp2Add = ElmtJacCoef[nelmt];
215 MatData0 = gmtxarray[m][n]->GetBlock(nelmt, nelmt)->GetPtr();
216 MatData1 = tmp2Add->GetPtr();
217 for (int i = 0; i < MatData0.size(); i++)
218 {
219 MatData0[i] += DataType(MatData1[i]);
220 }
221 }
222 }
223 }
224}
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
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)
Definition: Advection.cpp:245
std::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
std::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:75
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:354
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1191

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), CalcJacobTraceInteg(), Nektar::Array< OneD, const DataType >::size(), Vmath::Vadd(), and Vmath::Vcopy().

◆ AddTraceJacToMat() [2/2]

template<typename DataType , typename TypeNekBlkMatSharedPtr >
SOLVER_UTILS_EXPORT void Nektar::SolverUtils::Advection::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 
)

◆ Advect()

void Nektar::SolverUtils::Advection::Advect ( const int  nConvectiveFields,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const Array< OneD, Array< OneD, NekDouble > > &  pAdvVel,
const Array< OneD, Array< OneD, NekDouble > > &  pInarray,
Array< OneD, Array< OneD, NekDouble > > &  pOutarray,
const NekDouble pTime,
const Array< OneD, Array< OneD, NekDouble > > &  pFwd = NullNekDoubleArrayOfArray,
const Array< OneD, Array< OneD, NekDouble > > &  pBwd = NullNekDoubleArrayOfArray 
)

Interface function to advect the vector field.

Parameters
nConvectiveFieldsNumber of velocity components.
pFieldsExpansion lists for scalar fields.
pAdvVelAdvection velocity.
pInarrayScalar data to advect.
pOutarrayAdvected scalar data.
pTimeSimulation time.

Definition at line 71 of file Advection.cpp.

79{
80 v_Advect(nConvectiveFields, pFields, pAdvVel, pInarray, pOutarray, pTime,
81 pFwd, pBwd);
82}
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.

References v_Advect().

◆ CalcJacobTraceInteg()

template<typename DataType , typename TypeNekBlkMatSharedPtr >
void Nektar::SolverUtils::Advection::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 
)

Definition at line 245 of file Advection.cpp.

251{
252 MultiRegions::ExpListSharedPtr tracelist = pFields[0]->GetTrace();
253 std::shared_ptr<LocalRegions::ExpansionVector> traceExp =
254 tracelist->GetExp();
255 int ntotTrac = (*traceExp).size();
256 int nTracPnt, noffset, pntoffset;
257
258 Array<OneD, int> tracepnts(ntotTrac);
259 Array<OneD, Array<OneD, NekDouble>> JacFwd(ntotTrac);
260 Array<OneD, Array<OneD, NekDouble>> JacBwd(ntotTrac);
261
262 for (int nelmt = 0; nelmt < ntotTrac; nelmt++)
263 {
264 nTracPnt = (*traceExp)[nelmt]->GetTotPoints();
265 tracepnts[nelmt] = nTracPnt;
266
267 JacFwd[nelmt] = Array<OneD, NekDouble>(nTracPnt, 0.0);
268 JacBwd[nelmt] = Array<OneD, NekDouble>(nTracPnt, 0.0);
269 }
270
271 DataType ftmp, btmp;
272 for (int nelmt = 0; nelmt < ntotTrac; nelmt++)
273 {
274 nTracPnt = tracepnts[nelmt];
275 noffset = tracelist->GetPhys_Offset(nelmt);
276 for (int npnt = 0; npnt < nTracPnt; npnt++)
277 {
278 pntoffset = noffset + npnt;
279 ftmp = (*PntJac[0]->GetBlock(pntoffset, pntoffset))(m, n);
280 JacFwd[nelmt][npnt] = NekDouble(PntJacSign[0][pntoffset] * ftmp);
281
282 btmp = (*PntJac[1]->GetBlock(pntoffset, pntoffset))(m, n);
283 JacBwd[nelmt][npnt] = NekDouble(PntJacSign[1][pntoffset] * btmp);
284 }
285 }
286 tracelist->GetDiagMatIpwrtBase(JacFwd, TraceJacFwd);
287 tracelist->GetDiagMatIpwrtBase(JacBwd, TraceJacBwd);
288}
double NekDouble

Referenced by AddTraceJacToMat().

◆ InitObject()

void Nektar::SolverUtils::Advection::InitObject ( LibUtilities::SessionReaderSharedPtr  pSession,
Array< OneD, MultiRegions::ExpListSharedPtr pFields 
)

Interface function to initialise the advection object.

Parameters
pSessionSession configuration data.
pFieldsArray of ExpList objects.

Definition at line 57 of file Advection.cpp.

59{
60 v_InitObject(pSession, pFields);
61}
virtual SOLVER_UTILS_EXPORT void v_InitObject(LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)
Initialises the advection object.
Definition: Advection.cpp:299

References v_InitObject().

◆ SetBaseFlow()

void Nektar::SolverUtils::Advection::SetBaseFlow ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  fields 
)
inline

Set the base flow used for linearised advection objects.

Parameters
inarrayVector to use as baseflow

Definition at line 143 of file Advection.h.

146 {
147 v_SetBaseFlow(inarray, fields);
148 }
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.
Definition: Advection.cpp:326

References v_SetBaseFlow().

◆ SetFluxVector() [1/2]

void Nektar::SolverUtils::Advection::SetFluxVector ( AdvectionFluxVecCB  fluxVector)
inline

Set the flux vector callback function.

Parameters
fluxVectorThe callback function to override.

Definition at line 133 of file Advection.h.

134 {
135 m_fluxVector = fluxVector;
136 }
AdvectionFluxVecCB m_fluxVector
Callback function to the flux vector (set when advection is in conservative form).
Definition: Advection.h:170

References m_fluxVector.

◆ SetFluxVector() [2/2]

template<typename FuncPointerT , typename ObjectPointerT >
void Nektar::SolverUtils::Advection::SetFluxVector ( FuncPointerT  func,
ObjectPointerT  obj 
)
inline

Set the flux vector callback function.

This routine is a utility function to avoid the explicit use of std::bind. A function and object can be passed to this function instead.

Definition at line 112 of file Advection.h.

113 {
115 std::bind(func, obj, std::placeholders::_1, std::placeholders::_2);
116 }

References m_fluxVector.

◆ SetRiemannSolver()

void Nektar::SolverUtils::Advection::SetRiemannSolver ( RiemannSolverSharedPtr  riemann)
inline

Set a Riemann solver object for this advection object.

Parameters
riemannThe RiemannSolver object.

Definition at line 123 of file Advection.h.

124 {
125 m_riemann = riemann;
126 }
RiemannSolverSharedPtr m_riemann
Riemann solver for DG-type schemes.
Definition: Advection.h:172

References m_riemann.

◆ v_Advect()

virtual SOLVER_UTILS_EXPORT void Nektar::SolverUtils::Advection::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 
)
protectedpure virtual

◆ v_InitObject()

void Nektar::SolverUtils::Advection::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 in Nektar::SolverUtils::Advection3DHomogeneous1D, Nektar::SolverUtils::AdvectionFR, Nektar::SolverUtils::AdvectionNonConservative, Nektar::SolverUtils::AdvectionWeakDG, Nektar::AlternateSkewAdvection, Nektar::LinearisedAdvection, Nektar::NavierStokesAdvection, Nektar::NoAdvection, and Nektar::SkewSymmetricAdvection.

Definition at line 299 of file Advection.cpp.

302{
303 m_spaceDim = pFields[0]->GetCoordim(0);
304
305 if (pSession->DefinesSolverInfo("HOMOGENEOUS"))
306 {
307 std::string HomoStr = pSession->GetSolverInfo("HOMOGENEOUS");
308 if (HomoStr == "HOMOGENEOUS1D" || HomoStr == "Homogeneous1D" ||
309 HomoStr == "1D" || HomoStr == "Homo1D" ||
310 HomoStr == "HOMOGENEOUS2D" || HomoStr == "Homogeneous2D" ||
311 HomoStr == "2D" || HomoStr == "Homo2D")
312 {
313 m_spaceDim++;
314 }
315 else
316 {
318 "Only 1D homogeneous dimension supported.");
319 }
320 }
321}
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:209
int m_spaceDim
Storage for space dimension. Used for homogeneous extension.
Definition: Advection.h:174

References Nektar::ErrorUtil::efatal, m_spaceDim, and NEKERROR.

Referenced by InitObject(), Nektar::SolverUtils::AdvectionFR::v_InitObject(), Nektar::SolverUtils::AdvectionNonConservative::v_InitObject(), Nektar::SolverUtils::AdvectionWeakDG::v_InitObject(), Nektar::LinearisedAdvection::v_InitObject(), Nektar::NavierStokesAdvection::v_InitObject(), and Nektar::SkewSymmetricAdvection::v_InitObject().

◆ v_SetBaseFlow()

void Nektar::SolverUtils::Advection::v_SetBaseFlow ( const Array< OneD, Array< OneD, NekDouble > > &  inarray,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  fields 
)
protectedvirtual

Overrides the base flow used during linearised advection.

Reimplemented in Nektar::LinearisedAdvection.

Definition at line 326 of file Advection.cpp.

329{
330 boost::ignore_unused(inarray, fields);
332 "A baseflow is not appropriate for this advection type.");
333}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by SetBaseFlow().

Member Data Documentation

◆ m_fluxVector

AdvectionFluxVecCB Nektar::SolverUtils::Advection::m_fluxVector
protected

Callback function to the flux vector (set when advection is in conservative form).

Definition at line 170 of file Advection.h.

Referenced by Nektar::SolverUtils::AdvectionWeakDG::AdvectVolumeFlux(), SetFluxVector(), Nektar::SolverUtils::AdvectionFR::v_Advect(), and Nektar::SolverUtils::Advection3DHomogeneous1D::v_Advect().

◆ m_riemann

RiemannSolverSharedPtr Nektar::SolverUtils::Advection::m_riemann
protected

◆ m_spaceDim

int Nektar::SolverUtils::Advection::m_spaceDim
protected