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

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

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 80 of file Advection.h.

Constructor & Destructor Documentation

◆ ~Advection()

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

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 50 of file Advection.cpp.

57{
58 MultiRegions::ExpListSharedPtr tracelist = pFields[0]->GetTrace();
59 std::shared_ptr<LocalRegions::ExpansionVector> traceExp =
60 tracelist->GetExp();
61 int ntotTrac = (*traceExp).size();
62 int nTracePnts = tracelist->GetTotPoints();
63 int nTracPnt, nTracCoef;
64
65 DNekMatSharedPtr tmp2Add;
66 Array<OneD, DataType> MatData0;
67 Array<OneD, NekDouble> MatData1;
68 Array<OneD, NekDouble> MatData2;
69
70 Array<OneD, DNekMatSharedPtr> TraceJacFwd(ntotTrac);
71 Array<OneD, DNekMatSharedPtr> TraceJacBwd(ntotTrac);
72
73 for (int nelmt = 0; nelmt < ntotTrac; nelmt++)
74 {
75 nTracCoef = (*traceExp)[nelmt]->GetNcoeffs();
76 nTracPnt = (*traceExp)[nelmt]->GetTotPoints();
77 TraceJacFwd[nelmt] =
78 MemoryManager<DNekMat>::AllocateSharedPtr(nTracCoef, nTracPnt, 0.0);
79 TraceJacBwd[nelmt] =
80 MemoryManager<DNekMat>::AllocateSharedPtr(nTracCoef, nTracPnt, 0.0);
81 }
82
83 std::shared_ptr<LocalRegions::ExpansionVector> fieldExp =
84 pFields[0]->GetExp();
85 int nTotElmt = (*fieldExp).size();
86 int nElmtPnt, nElmtCoef;
87
88 Array<OneD, DNekMatSharedPtr> ElmtJacQuad(nTotElmt);
89 Array<OneD, DNekMatSharedPtr> ElmtJacCoef(nTotElmt);
90
91 Array<OneD, NekDouble> SymmMatData;
92
93 for (int nelmt = 0; nelmt < nTotElmt; nelmt++)
94 {
95 nElmtCoef = (*fieldExp)[nelmt]->GetNcoeffs();
96 nElmtPnt = (*fieldExp)[nelmt]->GetTotPoints();
97
98 ElmtJacQuad[nelmt] =
99 MemoryManager<DNekMat>::AllocateSharedPtr(nElmtCoef, nElmtPnt, 0.0);
101 nElmtCoef, nElmtCoef, 0.0);
102 }
103
104 bool TracePntJacGradflag = true;
105
106 Array<OneD, Array<OneD, DataType>> TraceJacConsSign(2);
107 for (int i = 0; i < 2; i++)
108 {
109 TraceJacConsSign[i] = Array<OneD, DataType>(nTracePnts, 1.0);
110 }
111
112 if (0 == TracePntJacGrad.size())
113 {
114 TracePntJacGradflag = false;
115 }
116
117 for (int m = 0; m < nConvectiveFields; m++)
118 {
119 for (int n = 0; n < nConvectiveFields; n++)
120 {
121 // ElmtJacCons to set 0
122 for (int nelmt = 0; nelmt < nTotElmt; nelmt++)
123 {
124 (*ElmtJacCoef[nelmt]) = 0.0;
125 (*ElmtJacQuad[nelmt]) = 0.0;
126 }
127
128 if (TracePntJacGradflag)
129 {
130 // only BJac is stored here because BJac = - FJac
131 for (int ndir = 0; ndir < nSpaceDim; ndir++)
132 {
133 // ElmtJacGrad to set 0
134 for (int nelmt = 0; nelmt < nTotElmt; nelmt++)
135 {
136 (*ElmtJacQuad[nelmt]) = 0.0;
137 }
138 int ngrad = n * nSpaceDim + ndir;
139
140 CalcJacobTraceInteg(pFields, m, ngrad, TracePntJacGrad,
141 TracePntJacGradSign, TraceJacFwd,
142 TraceJacBwd);
143 pFields[0]->AddTraceJacToElmtJac(TraceJacFwd, TraceJacBwd,
144 ElmtJacQuad);
145 // cost can be lower down by doing 3 directions together
146 pFields[0]->AddRightIPTPhysDerivBase(ndir, ElmtJacQuad,
147 ElmtJacCoef);
148 }
149 for (int nelmt = 0; nelmt < nTotElmt; nelmt++)
150 {
151 nElmtCoef = (*fieldExp)[nelmt]->GetNcoeffs();
152 if (SymmMatData.size() < nElmtCoef)
153 {
154 SymmMatData =
155 Array<OneD, NekDouble>(nElmtCoef * nElmtCoef);
156 }
157 MatData1 = ElmtJacCoef[nelmt]->GetPtr();
158 for (int i = 0; i < nElmtCoef; i++)
159 {
160 Vmath::Vcopy(nElmtCoef, &MatData1[i * nElmtCoef], 1,
161 &SymmMatData[i], nElmtCoef);
162 }
163 Vmath::Vadd(nElmtCoef * nElmtCoef, SymmMatData, 1, MatData1,
164 1, MatData1, 1);
165 }
166 }
167
168 CalcJacobTraceInteg(pFields, m, n, TracePntJacCons,
169 TraceJacConsSign, TraceJacFwd, TraceJacBwd);
170
171 pFields[0]->AddTraceJacToElmtJac(TraceJacFwd, TraceJacBwd,
172 ElmtJacQuad);
173 // Memory can be save to explore reusing ElmtJacQuad as output
174 pFields[0]->AddRightIPTBaseMatrix(ElmtJacQuad, ElmtJacCoef);
175
176 // add ElmtJacCons to gmtxarray[m][n]
177 for (int nelmt = 0; nelmt < nTotElmt; nelmt++)
178 {
179 tmp2Add = ElmtJacCoef[nelmt];
180 MatData0 = gmtxarray[m][n]->GetBlock(nelmt, nelmt)->GetPtr();
181 MatData1 = tmp2Add->GetPtr();
182 for (int i = 0; i < MatData0.size(); i++)
183 {
184 MatData0[i] += DataType(MatData1[i]);
185 }
186 }
187 }
188 }
189}
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:210
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.hpp:180
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825

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

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

Interface function to advect the vector field.

Definition at line 92 of file Advection.h.

102 {
103 v_Advect(nConvectiveFields, fields, advVel, inarray, outarray, time,
104 pFwd, pBwd);
105 }
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 210 of file Advection.cpp.

216{
217 MultiRegions::ExpListSharedPtr tracelist = pFields[0]->GetTrace();
218 std::shared_ptr<LocalRegions::ExpansionVector> traceExp =
219 tracelist->GetExp();
220 int ntotTrac = (*traceExp).size();
221 int nTracPnt, noffset, pntoffset;
222
223 Array<OneD, int> tracepnts(ntotTrac);
224 Array<OneD, Array<OneD, NekDouble>> JacFwd(ntotTrac);
225 Array<OneD, Array<OneD, NekDouble>> JacBwd(ntotTrac);
226
227 for (int nelmt = 0; nelmt < ntotTrac; nelmt++)
228 {
229 nTracPnt = (*traceExp)[nelmt]->GetTotPoints();
230 tracepnts[nelmt] = nTracPnt;
231
232 JacFwd[nelmt] = Array<OneD, NekDouble>(nTracPnt, 0.0);
233 JacBwd[nelmt] = Array<OneD, NekDouble>(nTracPnt, 0.0);
234 }
235
236 DataType ftmp, btmp;
237 for (int nelmt = 0; nelmt < ntotTrac; nelmt++)
238 {
239 nTracPnt = tracepnts[nelmt];
240 noffset = tracelist->GetPhys_Offset(nelmt);
241 for (int npnt = 0; npnt < nTracPnt; npnt++)
242 {
243 pntoffset = noffset + npnt;
244 ftmp = (*PntJac[0]->GetBlock(pntoffset, pntoffset))(m, n);
245 JacFwd[nelmt][npnt] = NekDouble(PntJacSign[0][pntoffset] * ftmp);
246
247 btmp = (*PntJac[1]->GetBlock(pntoffset, pntoffset))(m, n);
248 JacBwd[nelmt][npnt] = NekDouble(PntJacSign[1][pntoffset] * btmp);
249 }
250 }
251 tracelist->GetDiagMatIpwrtBase(JacFwd, TraceJacFwd);
252 tracelist->GetDiagMatIpwrtBase(JacBwd, TraceJacBwd);
253}
double NekDouble

Referenced by AddTraceJacToMat().

◆ InitObject()

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

Interface function to initialise the advection object.

Definition at line 84 of file Advection.h.

87 {
88 v_InitObject(pSession, pFields);
89 }
virtual SOLVER_UTILS_EXPORT void v_InitObject(LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)
Initialises the advection object.
Definition: Advection.cpp:264

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 146 of file Advection.h.

149 {
150 v_SetBaseFlow(inarray, fields);
151 }
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:291

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 136 of file Advection.h.

137 {
138 m_fluxVector = fluxVector;
139 }
AdvectionFluxVecCB m_fluxVector
Callback function to the flux vector (set when advection is in conservative form).
Definition: Advection.h:175

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 115 of file Advection.h.

116 {
118 std::bind(func, obj, std::placeholders::_1, std::placeholders::_2);
119 }

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 126 of file Advection.h.

127 {
128 m_riemann = riemann;
129 }
RiemannSolverSharedPtr m_riemann
Riemann solver for DG-type schemes.
Definition: Advection.h:177

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 264 of file Advection.cpp.

267{
268 m_spaceDim = pFields[0]->GetCoordim(0);
269
270 if (pSession->DefinesSolverInfo("HOMOGENEOUS"))
271 {
272 std::string HomoStr = pSession->GetSolverInfo("HOMOGENEOUS");
273 if (HomoStr == "HOMOGENEOUS1D" || HomoStr == "Homogeneous1D" ||
274 HomoStr == "1D" || HomoStr == "Homo1D" ||
275 HomoStr == "HOMOGENEOUS2D" || HomoStr == "Homogeneous2D" ||
276 HomoStr == "2D" || HomoStr == "Homo2D")
277 {
278 m_spaceDim++;
279 }
280 else
281 {
283 "Only 1D homogeneous dimension supported.");
284 }
285 }
286}
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:202
int m_spaceDim
Storage for space dimension. Used for homogeneous extension.
Definition: Advection.h:179

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 291 of file Advection.cpp.

294{
296 "A baseflow is not appropriate for this advection type.");
297}

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