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

Constructor & Destructor Documentation

◆ ~Advection()

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

Definition at line 83 of file Advection.h.

83{};

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

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

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

75{
76 v_Advect(nConvectiveFields, pFields, pAdvVel, pInarray, pOutarray, pTime,
77 pFwd, pBwd);
78}
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 241 of file Advection.cpp.

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

55{
56 v_InitObject(pSession, pFields);
57}
virtual SOLVER_UTILS_EXPORT void v_InitObject(LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)
Initialises the advection object.
Definition: Advection.cpp:295

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

144 {
145 v_SetBaseFlow(inarray, fields);
146 }
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:322

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

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

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

111 {
113 std::bind(func, obj, std::placeholders::_1, std::placeholders::_2);
114 }

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

122 {
123 m_riemann = riemann;
124 }
RiemannSolverSharedPtr m_riemann
Riemann solver for DG-type schemes.
Definition: Advection.h:170

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

298{
299 m_spaceDim = pFields[0]->GetCoordim(0);
300
301 if (pSession->DefinesSolverInfo("HOMOGENEOUS"))
302 {
303 std::string HomoStr = pSession->GetSolverInfo("HOMOGENEOUS");
304 if (HomoStr == "HOMOGENEOUS1D" || HomoStr == "Homogeneous1D" ||
305 HomoStr == "1D" || HomoStr == "Homo1D" ||
306 HomoStr == "HOMOGENEOUS2D" || HomoStr == "Homogeneous2D" ||
307 HomoStr == "2D" || HomoStr == "Homo2D")
308 {
309 m_spaceDim++;
310 }
311 else
312 {
314 "Only 1D homogeneous dimension supported.");
315 }
316 }
317}
#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:172

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

325{
327 "A baseflow is not appropriate for this advection type.");
328}

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