Nektar++
Static Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | Friends | List of all members
Nektar::SubSteppingExtrapolate Class Reference

#include <SubSteppingExtrapolate.h>

Inheritance diagram for Nektar::SubSteppingExtrapolate:
[legend]

Static Public Member Functions

static ExtrapolateSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession, Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, MultiRegions::ExpListSharedPtr &pPressure, const Array< OneD, int > &pVel, const SolverUtils::AdvectionSharedPtr &advObject)
 Creates an instance of this class. More...
 

Static Public Attributes

static std::string className
 Name of class. More...
 

Protected Member Functions

 SubSteppingExtrapolate (const LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields, MultiRegions::ExpListSharedPtr pPressure, const Array< OneD, int > pVel, const SolverUtils::AdvectionSharedPtr advObject)
 
 ~SubSteppingExtrapolate () override=default
 
void v_EvaluatePressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis) override
 
void v_SubSteppingTimeIntegration (const LibUtilities::TimeIntegrationSchemeSharedPtr &IntegrationScheme) override
 
void v_SubStepSaveFields (int nstep) override
 
void v_SubStepSetPressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &inarray, NekDouble Aii_Dt, NekDouble kinvis) override
 
void v_AccelerationBDF (Array< OneD, Array< OneD, NekDouble > > &array) override
 
void v_SubStepAdvance (int nstep, NekDouble time) override
 
void v_MountHOPBCs (int HBCdata, NekDouble kinvis, Array< OneD, NekDouble > &Q, Array< OneD, const NekDouble > &Advection) override
 
std::string v_GetSubStepName (void) override
 
void SubStepAdvection (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
 
void SubStepProjection (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
 
void SubStepExtrapolateField (NekDouble toff, Array< OneD, Array< OneD, NekDouble > > &ExtVel, Array< OneD, NekDouble > &ExtVn)
 
void AddAdvectionPenaltyFlux (const Array< OneD, NekDouble > &Vn, const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
NekDouble GetSubstepTimeStep ()
 
- Protected Member Functions inherited from Nektar::Extrapolate
virtual void v_EvaluatePressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &inarray, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis)=0
 
virtual void v_SubSteppingTimeIntegration (const LibUtilities::TimeIntegrationSchemeSharedPtr &IntegrationScheme)=0
 
virtual void v_SubStepSaveFields (int nstep)=0
 
virtual void v_SubStepSetPressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &inarray, NekDouble Aii_DT, NekDouble kinvis)=0
 
virtual void v_SubStepAdvance (int nstep, NekDouble time)=0
 
virtual void v_MountHOPBCs (int HBCdata, NekDouble kinvis, Array< OneD, NekDouble > &Q, Array< OneD, const NekDouble > &Advection)=0
 
virtual std::string v_GetSubStepName (void)
 
virtual void v_AccelerationBDF (Array< OneD, Array< OneD, NekDouble > > &array)
 
virtual void v_CalcNeumannPressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis)
 
virtual void v_CorrectPressureBCs (const Array< OneD, NekDouble > &pressure)
 
virtual void v_AddNormVelOnOBC (const int nbcoeffs, const int nreg, Array< OneD, Array< OneD, NekDouble > > &u)
 
void CalcOutflowBCs (const Array< OneD, const Array< OneD, NekDouble > > &fields, NekDouble kinvis)
 
void RollOver (Array< OneD, Array< OneD, NekDouble > > &input)
 

Protected Attributes

LibUtilities::TimeIntegrationSchemeSharedPtr m_intScheme
 
LibUtilities::TimeIntegrationSchemeSharedPtr m_subStepIntegrationScheme
 
LibUtilities::TimeIntegrationSchemeOperators m_subStepIntegrationOps
 
Array< OneD, Array< OneD, NekDouble > > m_previousVelFields
 
Array< OneD, Array< OneD, NekDouble > > m_previousVnFields
 
NekDouble m_cflSafetyFactor
 
int m_infosteps
 
int m_minsubsteps
 
int m_maxsubsteps
 
- Protected Attributes inherited from Nektar::Extrapolate
LibUtilities::SessionReaderSharedPtr m_session
 
LibUtilities::CommSharedPtr m_comm
 
Array< OneD, HBCTypem_hbcType
 Array of type of high order BCs for splitting shemes. More...
 
Array< OneD, MultiRegions::ExpListSharedPtrm_fields
 Velocity fields. More...
 
MultiRegions::ExpListSharedPtr m_pressure
 Pointer to field holding pressure field. More...
 
Array< OneD, int > m_velocity
 int which identifies which components of m_fields contains the velocity (u,v,w); More...
 
SolverUtils::AdvectionSharedPtr m_advObject
 
std::vector< SolverUtils::ForcingSharedPtrm_forcing
 
Array< OneD, Array< OneD, NekDouble > > m_previousVelFields
 
int m_curl_dim
 Curl-curl dimensionality. More...
 
int m_bnd_dim
 bounday dimensionality More...
 
Array< OneD, const SpatialDomains::BoundaryConditionShPtrm_PBndConds
 pressure boundary conditions container More...
 
Array< OneD, MultiRegions::ExpListSharedPtrm_PBndExp
 pressure boundary conditions expansion container More...
 
Array< OneD, MultiRegions::ExpListSharedPtrm_bndElmtExps
 Boundary expansions on each domain boundary. More...
 
int m_pressureCalls
 number of times the high-order pressure BCs have been called More...
 
int m_numHBCDof
 
int m_HBCnumber
 
int m_intSteps
 Maximum points used in pressure BC evaluation. More...
 
NekDouble m_timestep
 
Array< OneD, Array< OneD, NekDouble > > m_pressureHBCs
 Storage for current and previous levels of high order pressure boundary conditions. More...
 
Array< OneD, Array< OneD, NekDouble > > m_iprodnormvel
 Storage for current and previous levels of the inner product of normal velocity. More...
 
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
 
HighOrderOutflowSharedPtr m_houtflow
 

Friends

class MemoryManager< SubSteppingExtrapolate >
 

Additional Inherited Members

- Public Member Functions inherited from Nektar::Extrapolate
 Extrapolate (const LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields, MultiRegions::ExpListSharedPtr pPressure, const Array< OneD, int > pVel, const SolverUtils::AdvectionSharedPtr advObject)
 
virtual ~Extrapolate ()=default
 
void SubSteppingTimeIntegration (const LibUtilities::TimeIntegrationSchemeSharedPtr &IntegrationScheme)
 
void SubStepSaveFields (const int nstep)
 
void SubStepSetPressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &inarray, const NekDouble Aii_DT, NekDouble kinvis)
 
void SubStepAdvance (const int nstep, NekDouble time)
 
void MountHOPBCs (int HBCdata, NekDouble kinvis, Array< OneD, NekDouble > &Q, Array< OneD, const NekDouble > &Advection)
 
void EvaluatePressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis)
 
void AddNormVelOnOBC (const int nbcoeffs, const int nreg, Array< OneD, Array< OneD, NekDouble > > &u)
 
void CorrectPressureBCs (const Array< OneD, NekDouble > &pressure)
 
void CalcNeumannPressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis)
 
std::string GetSubStepName (void)
 
void SetForcing (const std::vector< SolverUtils::ForcingSharedPtr > &forcing)
 
void GenerateHOPBCMap (const LibUtilities::SessionReaderSharedPtr &pSsession)
 
void UpdateRobinPrimCoeff (void)
 
void AddDuDt (void)
 
void AddVelBC (void)
 
void ExtrapolatePressureHBCs (void)
 
void CopyPressureHBCsToPbndExp (void)
 
Array< OneD, NekDoubleGetMaxStdVelocity (const Array< OneD, Array< OneD, NekDouble > > inarray)
 
void IProductNormVelocityOnHBC (const Array< OneD, const Array< OneD, NekDouble > > &Vel, Array< OneD, NekDouble > &IprodVn)
 
void IProductNormVelocityBCOnHBC (Array< OneD, NekDouble > &IprodVn)
 
void ExtrapolateArray (Array< OneD, Array< OneD, NekDouble > > &array)
 
void EvaluateBDFArray (Array< OneD, Array< OneD, NekDouble > > &array)
 
void ExtrapolateArray (Array< OneD, Array< OneD, NekDouble > > &oldarrays, Array< OneD, NekDouble > &newarray, Array< OneD, NekDouble > &outarray)
 
void AddPressureToOutflowBCs (NekDouble kinvis)
 
void GenerateBndElmtExpansion (void)
 
- Static Protected Attributes inherited from Nektar::Extrapolate
static NekDouble StifflyStable_Betaq_Coeffs [3][3]
 
static NekDouble StifflyStable_Alpha_Coeffs [3][3]
 
static NekDouble StifflyStable_Gamma0_Coeffs [3] = {1.0, 1.5, 11.0 / 6.0}
 

Detailed Description

Definition at line 57 of file SubSteppingExtrapolate.h.

Constructor & Destructor Documentation

◆ SubSteppingExtrapolate()

Nektar::SubSteppingExtrapolate::SubSteppingExtrapolate ( const LibUtilities::SessionReaderSharedPtr  pSession,
Array< OneD, MultiRegions::ExpListSharedPtr pFields,
MultiRegions::ExpListSharedPtr  pPressure,
const Array< OneD, int >  pVel,
const SolverUtils::AdvectionSharedPtr  advObject 
)
protected

Definition at line 51 of file SubSteppingExtrapolate.cpp.

56 : Extrapolate(pSession, pFields, pPressure, pVel, advObject)
57{
58 m_session->LoadParameter("IO_InfoSteps", m_infosteps, 0);
59 m_session->LoadParameter("SubStepCFL", m_cflSafetyFactor, 0.5);
60 m_session->LoadParameter("MinSubSteps", m_minsubsteps, 1);
61 m_session->LoadParameter("MaxSubSteps", m_maxsubsteps, 100);
62
63 size_t dim = m_fields[0]->GetCoordim(0);
64 m_traceNormals = Array<OneD, Array<OneD, NekDouble>>(dim);
65
66 size_t nTracePts = m_fields[0]->GetTrace()->GetNpoints();
67 for (size_t i = 0; i < dim; ++i)
68 {
69 m_traceNormals[i] = Array<OneD, NekDouble>(nTracePts);
70 }
71 m_fields[0]->GetTrace()->GetNormals(m_traceNormals);
72}
Extrapolate(const LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields, MultiRegions::ExpListSharedPtr pPressure, const Array< OneD, int > pVel, const SolverUtils::AdvectionSharedPtr advObject)
Definition: Extrapolate.cpp:53
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Velocity fields.
Definition: Extrapolate.h:231
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Definition: Extrapolate.h:283
LibUtilities::SessionReaderSharedPtr m_session
Definition: Extrapolate.h:223

References m_cflSafetyFactor, Nektar::Extrapolate::m_fields, m_infosteps, m_maxsubsteps, m_minsubsteps, Nektar::Extrapolate::m_session, and Nektar::Extrapolate::m_traceNormals.

◆ ~SubSteppingExtrapolate()

Nektar::SubSteppingExtrapolate::~SubSteppingExtrapolate ( )
overrideprotecteddefault

Member Function Documentation

◆ AddAdvectionPenaltyFlux()

void Nektar::SubSteppingExtrapolate::AddAdvectionPenaltyFlux ( const Array< OneD, NekDouble > &  Vn,
const Array< OneD, const Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, NekDouble > > &  Outarray 
)
protected

Add the advection penalty term \( (\hat{u} - u^e)V_n \) given the normal velocity Vn at this time level and the physfield values containing the velocity field at this time level

Number of trace points

Forward state array

Backward state array

upwind numerical flux state array

Extract forwards/backwards trace spaces Note it is important to use the zeroth field but with the specialised method to use boudnary conditions from other fields since trace spaces may not be the same if there are mixed boundary conditions

Upwind between elements

Construct difference between numflux and Fwd,Bwd

Calculate the numerical fluxes multipling Fwd, Bwd and numflux by the normal advection velocity

Definition at line 519 of file SubSteppingExtrapolate.cpp.

523{
524 ASSERTL1(physfield.size() == Outarray.size(),
525 "Physfield and outarray are of different dimensions");
526
527 size_t i;
529 timer.Start();
530
531 /// Number of trace points
532 size_t nTracePts = m_fields[0]->GetTrace()->GetNpoints();
533
534 /// Forward state array
535 Array<OneD, NekDouble> Fwd(3 * nTracePts);
536
537 /// Backward state array
538 Array<OneD, NekDouble> Bwd = Fwd + nTracePts;
539
540 /// upwind numerical flux state array
541 Array<OneD, NekDouble> numflux = Bwd + nTracePts;
542
543 for (i = 0; i < physfield.size(); ++i)
544 {
545 /// Extract forwards/backwards trace spaces
546 /// Note it is important to use the zeroth field but with the
547 /// specialised method to use boudnary conditions from other
548 /// fields since trace spaces may not be the same if there are
549 /// mixed boundary conditions
550 std::dynamic_pointer_cast<MultiRegions::DisContField>(m_fields[0])
551 ->GetFwdBwdTracePhys(physfield[i], Fwd, Bwd,
552 m_fields[i]->GetBndConditions(),
553 m_fields[i]->GetBndCondExpansions());
554
555 /// Upwind between elements
556 m_fields[0]->GetTrace()->Upwind(Vn, Fwd, Bwd, numflux);
557
558 /// Construct difference between numflux and Fwd,Bwd
559 Vmath::Vsub(nTracePts, numflux, 1, Fwd, 1, Fwd, 1);
560 Vmath::Vsub(nTracePts, numflux, 1, Bwd, 1, Bwd, 1);
561
562 /// Calculate the numerical fluxes multipling Fwd, Bwd and
563 /// numflux by the normal advection velocity
564 Vmath::Vmul(nTracePts, Fwd, 1, Vn, 1, Fwd, 1);
565 Vmath::Vmul(nTracePts, Bwd, 1, Vn, 1, Bwd, 1);
566
567 m_fields[0]->AddFwdBwdTraceIntegral(Fwd, Bwd, Outarray[i]);
568 }
569 timer.Stop();
570 timer.AccumulateRegion("SubSteppingExtrapolate:AddAdvectionPenaltyFlux",
571 10);
572}
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:242
void AccumulateRegion(std::string, int iolevel=0)
Accumulate elapsed time for a region.
Definition: Timer.cpp:70
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 Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.
Definition: Vmath.hpp:220

References Nektar::LibUtilities::Timer::AccumulateRegion(), ASSERTL1, Nektar::Extrapolate::m_fields, Nektar::LibUtilities::Timer::Start(), Nektar::LibUtilities::Timer::Stop(), Vmath::Vmul(), and Vmath::Vsub().

Referenced by SubStepAdvection().

◆ create()

static ExtrapolateSharedPtr Nektar::SubSteppingExtrapolate::create ( const LibUtilities::SessionReaderSharedPtr pSession,
Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
MultiRegions::ExpListSharedPtr pPressure,
const Array< OneD, int > &  pVel,
const SolverUtils::AdvectionSharedPtr advObject 
)
inlinestatic

Creates an instance of this class.

Definition at line 63 of file SubSteppingExtrapolate.h.

68 {
71 pSession, pFields, pPressure, pVel, advObject);
72 return p;
73 }
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
std::shared_ptr< Extrapolate > ExtrapolateSharedPtr
Definition: Extrapolate.h:60

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and CellMLToNektar.cellml_metadata::p.

◆ GetSubstepTimeStep()

NekDouble Nektar::SubSteppingExtrapolate::GetSubstepTimeStep ( )
protected

Definition at line 478 of file SubSteppingExtrapolate.cpp.

479{
481 timer.Start();
482 size_t n_element = m_fields[0]->GetExpSize();
483
484 const Array<OneD, int> ExpOrder = m_fields[0]->EvalBasisNumModesMaxPerExp();
485
486 const NekDouble cLambda = 0.2; // Spencer book pag. 317
487
488 Array<OneD, NekDouble> tstep(n_element, 0.0);
489 Array<OneD, NekDouble> stdVelocity(n_element, 0.0);
490 Array<OneD, Array<OneD, NekDouble>> velfields(m_velocity.size());
491
492 for (size_t i = 0; i < m_velocity.size(); ++i)
493 {
494 velfields[i] = m_fields[m_velocity[i]]->UpdatePhys();
495 }
496 stdVelocity = GetMaxStdVelocity(velfields);
497
498 for (size_t el = 0; el < n_element; ++el)
499 {
500 tstep[el] =
501 m_cflSafetyFactor / (stdVelocity[el] * cLambda *
502 (ExpOrder[el] - 1) * (ExpOrder[el] - 1));
503 }
504
505 NekDouble TimeStep = Vmath::Vmin(n_element, tstep, 1);
506 m_comm->AllReduce(TimeStep, LibUtilities::ReduceMin);
507 timer.Stop();
508 timer.AccumulateRegion("SubSteppingExtrapolate:GetSubStepTimeStep", 10);
509
510 return TimeStep;
511}
Array< OneD, NekDouble > GetMaxStdVelocity(const Array< OneD, Array< OneD, NekDouble > > inarray)
Array< OneD, int > m_velocity
int which identifies which components of m_fields contains the velocity (u,v,w);
Definition: Extrapolate.h:238
LibUtilities::CommSharedPtr m_comm
Definition: Extrapolate.h:225
double NekDouble
T Vmin(int n, const T *x, const int incx)
Return the minimum element in x - called vmin to avoid conflict with min.
Definition: Vmath.hpp:725

References Nektar::LibUtilities::Timer::AccumulateRegion(), Nektar::Extrapolate::GetMaxStdVelocity(), m_cflSafetyFactor, Nektar::Extrapolate::m_comm, Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_velocity, Nektar::LibUtilities::ReduceMin, Nektar::LibUtilities::Timer::Start(), Nektar::LibUtilities::Timer::Stop(), and Vmath::Vmin().

Referenced by v_SubStepAdvance().

◆ SubStepAdvection()

void Nektar::SubSteppingExtrapolate::SubStepAdvection ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray,
const NekDouble  time 
)
protected

Explicit Advection terms used by SubStepAdvance time integration

Get the number of coefficients

Define an auxiliary variable to compute the RHS

Operations to compute the RHS

Multiply the flux by the inverse of the mass matrix

Store in outarray the physical values of the RHS

Definition at line 158 of file SubSteppingExtrapolate.cpp.

161{
162 Nektar::LibUtilities::Timer timer, timer1;
163 size_t i;
164 size_t nVariables = inarray.size();
165 size_t nQuadraturePts = inarray[0].size();
166
167 timer.Start();
168 /// Get the number of coefficients
169 size_t ncoeffs = m_fields[0]->GetNcoeffs();
170
171 /// Define an auxiliary variable to compute the RHS
172 Array<OneD, Array<OneD, NekDouble>> WeakAdv(nVariables);
173 WeakAdv[0] = Array<OneD, NekDouble>(ncoeffs * nVariables);
174 for (i = 1; i < nVariables; ++i)
175 {
176 WeakAdv[i] = WeakAdv[i - 1] + ncoeffs;
177 }
178
179 Array<OneD, Array<OneD, NekDouble>> Velfields(m_velocity.size());
180
181 Velfields[0] = Array<OneD, NekDouble>(nQuadraturePts * m_velocity.size());
182
183 for (i = 1; i < m_velocity.size(); ++i)
184 {
185 Velfields[i] = Velfields[i - 1] + nQuadraturePts;
186 }
187
188 Array<OneD, NekDouble> Vn(m_fields[0]->GetTrace()->GetTotPoints());
189
190 SubStepExtrapolateField(fmod(time, m_timestep), Velfields, Vn);
191
192 for (auto &x : m_forcing)
193 {
194 x->PreApply(m_fields, Velfields, Velfields, time);
195 }
196 timer1.Start();
197 m_advObject->Advect(m_velocity.size(), m_fields, Velfields, inarray,
198 outarray, time);
199 timer1.Stop();
200 timer1.AccumulateRegion("SubStepAdvection:Advect", 10);
201
202 for (auto &x : m_forcing)
203 {
204 x->Apply(m_fields, outarray, outarray, time);
205 }
206
207 for (i = 0; i < nVariables; ++i)
208 {
209 m_fields[i]->IProductWRTBase(outarray[i], WeakAdv[i]);
210 // negation requried due to sign of DoAdvection term to be consistent
211 Vmath::Neg(ncoeffs, WeakAdv[i], 1);
212 }
213
214 AddAdvectionPenaltyFlux(Vn, inarray, WeakAdv);
215
216 /// Operations to compute the RHS
217 for (i = 0; i < nVariables; ++i)
218 {
219 // Negate the RHS
220 Vmath::Neg(ncoeffs, WeakAdv[i], 1);
221
222 /// Multiply the flux by the inverse of the mass matrix
223 m_fields[i]->MultiplyByElmtInvMass(WeakAdv[i], WeakAdv[i]);
224
225 /// Store in outarray the physical values of the RHS
226 m_fields[i]->BwdTrans(WeakAdv[i], outarray[i]);
227 }
228 timer.Stop();
229 timer.AccumulateRegion("SubSteppingExtrapolate:SubStepAdvection", 10);
230}
NekDouble m_timestep
Definition: Extrapolate.h:273
std::vector< SolverUtils::ForcingSharedPtr > m_forcing
Definition: Extrapolate.h:242
SolverUtils::AdvectionSharedPtr m_advObject
Definition: Extrapolate.h:240
void AddAdvectionPenaltyFlux(const Array< OneD, NekDouble > &Vn, const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &outarray)
void SubStepExtrapolateField(NekDouble toff, Array< OneD, Array< OneD, NekDouble > > &ExtVel, Array< OneD, NekDouble > &ExtVn)
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.hpp:292

References Nektar::LibUtilities::Timer::AccumulateRegion(), AddAdvectionPenaltyFlux(), Nektar::Extrapolate::m_advObject, Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_forcing, Nektar::Extrapolate::m_timestep, Nektar::Extrapolate::m_velocity, Vmath::Neg(), Nektar::LibUtilities::Timer::Start(), Nektar::LibUtilities::Timer::Stop(), and SubStepExtrapolateField().

Referenced by v_SubSteppingTimeIntegration().

◆ SubStepExtrapolateField()

void Nektar::SubSteppingExtrapolate::SubStepExtrapolateField ( NekDouble  toff,
Array< OneD, Array< OneD, NekDouble > > &  ExtVel,
Array< OneD, NekDouble > &  ExtVn 
)
protected

Extrapolate field using equally time spaced field un,un-1,un-2, (at dt intervals) to time for substep n+t at order m_intSteps. Also extrapolate the normal velocity along the trace at the same order

Definition at line 579 of file SubSteppingExtrapolate.cpp.

582{
583 size_t npts = m_fields[0]->GetTotPoints();
584 size_t ntracepts = m_fields[0]->GetTrace()->GetTotPoints();
585 size_t nvel = m_velocity.size();
586 size_t i, j;
587 Array<OneD, NekDouble> l(4);
589 timer.Start();
590
591 size_t ord = m_intSteps;
592
593 // calculate Lagrange interpolants
594 Vmath::Fill(4, 1.0, l, 1);
595
596 for (i = 0; i <= ord; ++i)
597 {
598 for (j = 0; j <= ord; ++j)
599 {
600 if (i != j)
601 {
602 l[i] *= (j * m_timestep + toff);
603 l[i] /= (j * m_timestep - i * m_timestep);
604 }
605 }
606 }
607
608 for (i = 0; i < nvel; ++i)
609 {
610 Vmath::Smul(npts, l[0], m_previousVelFields[i], 1, ExtVel[i], 1);
611
612 for (j = 1; j <= ord; ++j)
613 {
614 Blas::Daxpy(npts, l[j], m_previousVelFields[j * nvel + i], 1,
615 ExtVel[i], 1);
616 }
617 }
618
619 Vmath::Smul(ntracepts, l[0], m_previousVnFields[0], 1, ExtVn, 1);
620 for (j = 1; j <= ord; ++j)
621 {
622 Blas::Daxpy(ntracepts, l[j], m_previousVnFields[j], 1, ExtVn, 1);
623 }
624 timer.Stop();
625 timer.AccumulateRegion("SubSteppingExtrapolate:SubStepExtrapolateFields",
626 10);
627}
int m_intSteps
Maximum points used in pressure BC evaluation.
Definition: Extrapolate.h:271
Array< OneD, Array< OneD, NekDouble > > m_previousVnFields
Array< OneD, Array< OneD, NekDouble > > m_previousVelFields
static void Daxpy(const int &n, const double &alpha, const double *x, const int &incx, const double *y, const int &incy)
BLAS level 1: y = alpha x plus y.
Definition: Blas.hpp:135
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
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
Definition: Vmath.hpp:54

References Nektar::LibUtilities::Timer::AccumulateRegion(), Blas::Daxpy(), Vmath::Fill(), Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_intSteps, m_previousVelFields, m_previousVnFields, Nektar::Extrapolate::m_timestep, Nektar::Extrapolate::m_velocity, Vmath::Smul(), Nektar::LibUtilities::Timer::Start(), and Nektar::LibUtilities::Timer::Stop().

Referenced by SubStepAdvection().

◆ SubStepProjection()

void Nektar::SubSteppingExtrapolate::SubStepProjection ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray,
const NekDouble  time 
)
protected

Projection used by SubStepAdvance time integration

Definition at line 235 of file SubSteppingExtrapolate.cpp.

239{
241 ASSERTL1(inarray.size() == outarray.size(),
242 "Inarray and outarray of different sizes ");
243
244 timer.Start();
245 for (size_t i = 0; i < inarray.size(); ++i)
246 {
247 Vmath::Vcopy(inarray[i].size(), inarray[i], 1, outarray[i], 1);
248 }
249 timer.Stop();
250 timer.AccumulateRegion("SubSteppingExtrapolate:SubStepProjection", 10);
251}
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825

References Nektar::LibUtilities::Timer::AccumulateRegion(), ASSERTL1, Nektar::LibUtilities::Timer::Start(), Nektar::LibUtilities::Timer::Stop(), and Vmath::Vcopy().

Referenced by v_SubSteppingTimeIntegration().

◆ v_AccelerationBDF()

void Nektar::SubSteppingExtrapolate::v_AccelerationBDF ( Array< OneD, Array< OneD, NekDouble > > &  array)
overrideprotectedvirtual

At the start, the newest value is stored in array[nlevels-1] and the previous values in the first positions At the end, the acceleration from BDF is stored in array[nlevels-1] and the storage has been updated to included the new value

Reimplemented from Nektar::Extrapolate.

Definition at line 292 of file SubSteppingExtrapolate.cpp.

294{
296 int nlevels = array.size();
297 int nPts = array[0].size();
298
299 timer.Start();
300 if (nPts)
301 {
302 // Update array
303 RollOver(array);
304
305 // Calculate acceleration using Backward Differentiation Formula
306 Array<OneD, NekDouble> accelerationTerm(nPts, 0.0);
307
308 int acc_order = std::min(m_pressureCalls, m_intSteps);
309 Vmath::Smul(nPts, StifflyStable_Gamma0_Coeffs[acc_order - 1], array[0],
310 1, accelerationTerm, 1);
311
312 for (int i = 0; i < acc_order; i++)
313 {
315 nPts, -1 * StifflyStable_Alpha_Coeffs[acc_order - 1][i],
316 array[i + 1], 1, accelerationTerm, 1, accelerationTerm, 1);
317 }
318 array[nlevels - 1] = accelerationTerm;
319 }
320 timer.Stop();
321 timer.AccumulateRegion("SubSteppingExtrapolate:v_AccelerationBDF", 10);
322}
static NekDouble StifflyStable_Alpha_Coeffs[3][3]
Definition: Extrapolate.h:287
void RollOver(Array< OneD, Array< OneD, NekDouble > > &input)
int m_pressureCalls
number of times the high-order pressure BCs have been called
Definition: Extrapolate.h:262
static NekDouble StifflyStable_Gamma0_Coeffs[3]
Definition: Extrapolate.h:288
void Svtvp(int n, const T alpha, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Svtvp (scalar times vector plus vector): z = alpha*x + y.
Definition: Vmath.hpp:396

References Nektar::LibUtilities::Timer::AccumulateRegion(), Nektar::Extrapolate::m_intSteps, Nektar::Extrapolate::m_pressureCalls, Nektar::Extrapolate::RollOver(), Vmath::Smul(), Nektar::LibUtilities::Timer::Start(), Nektar::Extrapolate::StifflyStable_Alpha_Coeffs, Nektar::Extrapolate::StifflyStable_Gamma0_Coeffs, Nektar::LibUtilities::Timer::Stop(), and Vmath::Svtvp().

◆ v_EvaluatePressureBCs()

void Nektar::SubSteppingExtrapolate::v_EvaluatePressureBCs ( const Array< OneD, const Array< OneD, NekDouble > > &  fields,
const Array< OneD, const Array< OneD, NekDouble > > &  N,
NekDouble  kinvis 
)
overrideprotectedvirtual

Implements Nektar::Extrapolate.

Definition at line 74 of file SubSteppingExtrapolate.cpp.

78{
79 ASSERTL0(false, "This method should not be called by Substepping routine");
80}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208

References ASSERTL0.

◆ v_GetSubStepName()

std::string Nektar::SubSteppingExtrapolate::v_GetSubStepName ( void  )
overrideprotectedvirtual

Reimplemented from Nektar::Extrapolate.

Definition at line 639 of file SubSteppingExtrapolate.cpp.

640{
641 return m_subStepIntegrationScheme->GetName();
642}
LibUtilities::TimeIntegrationSchemeSharedPtr m_subStepIntegrationScheme

References m_subStepIntegrationScheme.

◆ v_MountHOPBCs()

void Nektar::SubSteppingExtrapolate::v_MountHOPBCs ( int  HBCdata,
NekDouble  kinvis,
Array< OneD, NekDouble > &  Q,
Array< OneD, const NekDouble > &  Advection 
)
overrideprotectedvirtual

Implements Nektar::Extrapolate.

Definition at line 632 of file SubSteppingExtrapolate.cpp.

635{
636 Vmath::Smul(HBCdata, -kinvis, Q, 1, Q, 1);
637}

References Vmath::Smul().

◆ v_SubStepAdvance()

void Nektar::SubSteppingExtrapolate::v_SubStepAdvance ( int  nstep,
NekDouble  time 
)
overrideprotectedvirtual

Implements Nektar::Extrapolate.

Definition at line 412 of file SubSteppingExtrapolate.cpp.

413{
414 int n;
415 int nsubsteps;
417
418 NekDouble dt;
419
420 Array<OneD, Array<OneD, NekDouble>> fields;
421
422 timer.Start();
423 static int ncalls = 1;
424 size_t nint = std::min(ncalls++, m_intSteps);
425
426 // this needs to change
427 m_comm = m_fields[0]->GetComm()->GetRowComm();
428
429 // Get the proper time step with CFL control
430 dt = GetSubstepTimeStep();
431
432 nsubsteps = (m_timestep > dt) ? ((int)(m_timestep / dt) + 1) : 1;
433 nsubsteps = std::max(m_minsubsteps, nsubsteps);
434
435 ASSERTL0(nsubsteps < m_maxsubsteps,
436 "Number of substeps has exceeded maximum");
437
438 dt = m_timestep / nsubsteps;
439
440 if (m_infosteps && !((nstep + 1) % m_infosteps) && m_comm->GetRank() == 0)
441 {
442 std::cout << "Sub-integrating using " << nsubsteps
443 << " steps over Dt = " << m_timestep
444 << " (SubStep CFL=" << m_cflSafetyFactor << ")" << std::endl;
445 }
446
447 const TripleArray &solutionVector = m_intScheme->GetSolutionVector();
448
449 for (size_t m = 0; m < nint; ++m)
450 {
451 // We need to update the fields held by the m_intScheme
452 fields = solutionVector[m];
453
454 // Initialise NS solver which is set up to use a GLM method
455 // with calls to EvaluateAdvection_SetPressureBCs and
456 // SolveUnsteadyStokesSystem
457 m_subStepIntegrationScheme->InitializeScheme(dt, fields, time,
459
460 for (n = 0; n < nsubsteps; ++n)
461 {
462 fields = m_subStepIntegrationScheme->TimeIntegrate(n, dt);
463 }
464
465 // set up HBC m_acceleration field for Pressure BCs
467
468 // Reset time integrated solution in m_intScheme
469 m_intScheme->SetSolutionVector(m, fields);
470 }
471 timer.Stop();
472 timer.AccumulateRegion("SubSteppingExtrapolate:v_SubStepAdvance", 10);
473}
Array< OneD, Array< OneD, NekDouble > > m_iprodnormvel
Storage for current and previous levels of the inner product of normal velocity.
Definition: Extrapolate.h:281
void IProductNormVelocityOnHBC(const Array< OneD, const Array< OneD, NekDouble > > &Vel, Array< OneD, NekDouble > &IprodVn)
LibUtilities::TimeIntegrationSchemeOperators m_subStepIntegrationOps
LibUtilities::TimeIntegrationSchemeSharedPtr m_intScheme
AT< AT< AT< NekDouble > > > TripleArray

References Nektar::LibUtilities::Timer::AccumulateRegion(), ASSERTL0, GetSubstepTimeStep(), Nektar::Extrapolate::IProductNormVelocityOnHBC(), m_cflSafetyFactor, Nektar::Extrapolate::m_comm, Nektar::Extrapolate::m_fields, m_infosteps, m_intScheme, Nektar::Extrapolate::m_intSteps, Nektar::Extrapolate::m_iprodnormvel, m_maxsubsteps, m_minsubsteps, m_subStepIntegrationOps, m_subStepIntegrationScheme, Nektar::Extrapolate::m_timestep, Nektar::LibUtilities::Timer::Start(), and Nektar::LibUtilities::Timer::Stop().

◆ v_SubSteppingTimeIntegration()

void Nektar::SubSteppingExtrapolate::v_SubSteppingTimeIntegration ( const LibUtilities::TimeIntegrationSchemeSharedPtr IntegrationScheme)
overrideprotectedvirtual

Implements Nektar::Extrapolate.

Definition at line 82 of file SubSteppingExtrapolate.cpp.

84{
86 m_intScheme = IntegrationScheme;
87
88 timer.Start();
89 size_t order = IntegrationScheme->GetOrder();
90
91 // Set to 1 for first step and it will then be increased in
92 // time advance routines
93 if ((IntegrationScheme->GetName() == "Euler" &&
94 IntegrationScheme->GetVariant() == "Backward") ||
95
96 (IntegrationScheme->GetName() == "BDFImplicit" &&
97 (order == 1 || order == 2)))
98 {
99 // Note RK first order SSP is just Forward Euler.
100 std::string vSubStepIntScheme = "RungeKutta";
101 std::string vSubStepIntSchemeVariant = "SSP";
102 int vSubStepIntSchemeOrder = order;
103
104 if (m_session->DefinesSolverInfo("SubStepIntScheme"))
105 {
106 vSubStepIntScheme = m_session->GetSolverInfo("SubStepIntScheme");
107 vSubStepIntSchemeVariant = "";
108 vSubStepIntSchemeOrder = order;
109 }
110
113 vSubStepIntScheme, vSubStepIntSchemeVariant,
114 vSubStepIntSchemeOrder, std::vector<NekDouble>());
115
116 size_t nvel = m_velocity.size();
117 size_t ndim = order + 1;
118
119 // Fields for linear/quadratic interpolation
120 m_previousVelFields = Array<OneD, Array<OneD, NekDouble>>(ndim * nvel);
121 int ntotpts = m_fields[0]->GetTotPoints();
122 m_previousVelFields[0] = Array<OneD, NekDouble>(ndim * nvel * ntotpts);
123
124 for (size_t i = 1; i < ndim * nvel; ++i)
125 {
126 m_previousVelFields[i] = m_previousVelFields[i - 1] + ntotpts;
127 }
128 // Vn fields
129 m_previousVnFields = Array<OneD, Array<OneD, NekDouble>>(ndim);
130 ntotpts = m_fields[0]->GetTrace()->GetTotPoints();
131 m_previousVnFields[0] = Array<OneD, NekDouble>(ndim * ntotpts);
132 for (size_t i = 1; i < ndim; ++i)
133 {
134 m_previousVnFields[i] = m_previousVnFields[i - 1] + ntotpts;
135 }
136 }
137 else
138 {
139 ASSERTL0(0, "Integration method not suitable: Options include "
140 "BackwardEuler or BDFImplicitOrder{1,2}");
141 }
142
143 m_intSteps = IntegrationScheme->GetNumIntegrationPhases();
144
145 // set explicit time-integration class operators
150 timer.Stop();
151 timer.AccumulateRegion(
152 "SubSteppingExtrapolate:v_SubSteppingTimeIntegration", 10);
153}
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
void DefineProjection(FuncPointerT func, ObjectPointerT obj)
void DefineOdeRhs(FuncPointerT func, ObjectPointerT obj)
void SubStepAdvection(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
void SubStepProjection(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
TimeIntegrationSchemeFactory & GetTimeIntegrationSchemeFactory()

References Nektar::LibUtilities::Timer::AccumulateRegion(), ASSERTL0, Nektar::LibUtilities::NekFactory< tKey, tBase, tParam >::CreateInstance(), Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineOdeRhs(), Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineProjection(), Nektar::LibUtilities::GetTimeIntegrationSchemeFactory(), Nektar::Extrapolate::m_fields, m_intScheme, Nektar::Extrapolate::m_intSteps, m_previousVelFields, m_previousVnFields, Nektar::Extrapolate::m_session, m_subStepIntegrationOps, m_subStepIntegrationScheme, Nektar::Extrapolate::m_velocity, Nektar::LibUtilities::Timer::Start(), Nektar::LibUtilities::Timer::Stop(), SubStepAdvection(), and SubStepProjection().

◆ v_SubStepSaveFields()

void Nektar::SubSteppingExtrapolate::v_SubStepSaveFields ( int  nstep)
overrideprotectedvirtual

Save the current fields in m_previousVelFields and cycle out older previous fields so that it can be extrapolated to the new substeps of the next time step. Also extract the normal trace of the velocity field into m_previousVnFields along the trace so that this can also be extrapolated.

Implements Nektar::Extrapolate.

Definition at line 331 of file SubSteppingExtrapolate.cpp.

332{
334 size_t i, n;
335 timer.Start();
336 size_t nvel = m_velocity.size();
337 size_t npts = m_fields[0]->GetTotPoints();
338 size_t ntracepts = m_fields[0]->GetTrace()->GetTotPoints();
339
340 // rotate fields
341 size_t nblocks = m_previousVelFields.size() / nvel;
342 Array<OneD, NekDouble> save;
343
344 // rotate storage space
345 for (n = 0; n < nvel; ++n)
346 {
347 save = m_previousVelFields[(nblocks - 1) * nvel + n];
348
349 for (i = nblocks - 1; i > 0; --i)
350 {
351 m_previousVelFields[i * nvel + n] =
352 m_previousVelFields[(i - 1) * nvel + n];
353 }
354
355 m_previousVelFields[n] = save;
356 }
357
358 save = m_previousVnFields[nblocks - 1];
359 for (i = nblocks - 1; i > 0; --i)
360 {
362 }
363 m_previousVnFields[0] = save;
364
365 // Put previous field
366 for (i = 0; i < nvel; ++i)
367 {
368 m_fields[m_velocity[i]]->BwdTrans(
369 m_fields[m_velocity[i]]->GetCoeffs(),
370 m_fields[m_velocity[i]]->UpdatePhys());
371 Vmath::Vcopy(npts, m_fields[m_velocity[i]]->GetPhys(), 1,
372 m_previousVelFields[i], 1);
373 }
374
375 Array<OneD, NekDouble> Fwd(ntracepts);
376
377 // calculate Vn
378 m_fields[0]->ExtractTracePhys(m_previousVelFields[0], Fwd);
379 Vmath::Vmul(ntracepts, m_traceNormals[0], 1, Fwd, 1, m_previousVnFields[0],
380 1);
381 for (i = 1; i < m_bnd_dim; ++i)
382 {
383 m_fields[0]->ExtractTracePhys(m_previousVelFields[i], Fwd);
384 Vmath::Vvtvp(ntracepts, m_traceNormals[i], 1, Fwd, 1,
386 }
387
388 if (nstep == 0) // initialise all levels with first field
389 {
390 for (n = 0; n < nvel; ++n)
391 {
392 for (i = 1; i < nblocks; ++i)
393 {
394 Vmath::Vcopy(npts, m_fields[m_velocity[n]]->GetPhys(), 1,
395 m_previousVelFields[i * nvel + n], 1);
396 }
397 }
398
399 for (i = 1; i < nblocks; ++i)
400 {
401 Vmath::Vcopy(ntracepts, m_previousVnFields[0], 1,
402 m_previousVnFields[i], 1);
403 }
404 }
405 timer.Stop();
406 timer.AccumulateRegion("SubSteppingExtrapolate:v_SubStepSaveFields", 10);
407}
int m_bnd_dim
bounday dimensionality
Definition: Extrapolate.h:250
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

References Nektar::LibUtilities::Timer::AccumulateRegion(), Nektar::Extrapolate::m_bnd_dim, Nektar::Extrapolate::m_fields, m_previousVelFields, m_previousVnFields, Nektar::Extrapolate::m_traceNormals, Nektar::Extrapolate::m_velocity, Nektar::LibUtilities::Timer::Start(), Nektar::LibUtilities::Timer::Stop(), Vmath::Vcopy(), Vmath::Vmul(), and Vmath::Vvtvp().

◆ v_SubStepSetPressureBCs()

void Nektar::SubSteppingExtrapolate::v_SubStepSetPressureBCs ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
NekDouble  Aii_Dt,
NekDouble  kinvis 
)
overrideprotectedvirtual

Implements Nektar::Extrapolate.

Reimplemented in Nektar::SubSteppingExtrapolateWeakPressure.

Definition at line 256 of file SubSteppingExtrapolate.cpp.

259{
261 // int nConvectiveFields =m_fields.size()-1;
262 Array<OneD, Array<OneD, NekDouble>> nullvelfields;
263
264 timer.Start();
266
267 // Calculate viscous BCs at current level and
268 // put in m_pressureHBCs[0]
269 CalcNeumannPressureBCs(m_previousVelFields, nullvelfields, kinvis);
270
271 // Extrapolate to m_pressureHBCs to n+1
273
274 // Add (phi,Du/Dt) term to m_presureHBC
275 AddDuDt();
276
277 // Copy m_pressureHBCs to m_PbndExp
279
280 // Evaluate High order outflow conditiosn if required.
281 CalcOutflowBCs(inarray, kinvis);
282 timer.Stop();
283 timer.AccumulateRegion("SubSteppingExtrapolate:SubStepSetPressureBCs", 10);
284}
Array< OneD, Array< OneD, NekDouble > > m_pressureHBCs
Storage for current and previous levels of high order pressure boundary conditions.
Definition: Extrapolate.h:277
void CopyPressureHBCsToPbndExp(void)
void CalcNeumannPressureBCs(const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis)
Definition: Extrapolate.h:132
void ExtrapolateArray(Array< OneD, Array< OneD, NekDouble > > &array)
void CalcOutflowBCs(const Array< OneD, const Array< OneD, NekDouble > > &fields, NekDouble kinvis)

References Nektar::LibUtilities::Timer::AccumulateRegion(), Nektar::Extrapolate::AddDuDt(), Nektar::Extrapolate::CalcNeumannPressureBCs(), Nektar::Extrapolate::CalcOutflowBCs(), Nektar::Extrapolate::CopyPressureHBCsToPbndExp(), Nektar::Extrapolate::ExtrapolateArray(), Nektar::Extrapolate::m_pressureCalls, Nektar::Extrapolate::m_pressureHBCs, m_previousVelFields, Nektar::LibUtilities::Timer::Start(), and Nektar::LibUtilities::Timer::Stop().

Friends And Related Function Documentation

◆ MemoryManager< SubSteppingExtrapolate >

friend class MemoryManager< SubSteppingExtrapolate >
friend

Definition at line 55 of file SubSteppingExtrapolate.h.

Member Data Documentation

◆ className

std::string Nektar::SubSteppingExtrapolate::className
static
Initial value:
=
"SubStepping", SubSteppingExtrapolate::create, "SubStepping")
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
static ExtrapolateSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, MultiRegions::ExpListSharedPtr &pPressure, const Array< OneD, int > &pVel, const SolverUtils::AdvectionSharedPtr &advObject)
Creates an instance of this class.
ExtrapolateFactory & GetExtrapolateFactory()
Definition: Extrapolate.cpp:47

Name of class.

Registers the class with the Factory.

Definition at line 76 of file SubSteppingExtrapolate.h.

◆ m_cflSafetyFactor

NekDouble Nektar::SubSteppingExtrapolate::m_cflSafetyFactor
protected

◆ m_infosteps

int Nektar::SubSteppingExtrapolate::m_infosteps
protected

Definition at line 87 of file SubSteppingExtrapolate.h.

Referenced by SubSteppingExtrapolate(), and v_SubStepAdvance().

◆ m_intScheme

LibUtilities::TimeIntegrationSchemeSharedPtr Nektar::SubSteppingExtrapolate::m_intScheme
protected

Definition at line 79 of file SubSteppingExtrapolate.h.

Referenced by v_SubStepAdvance(), and v_SubSteppingTimeIntegration().

◆ m_maxsubsteps

int Nektar::SubSteppingExtrapolate::m_maxsubsteps
protected

Definition at line 89 of file SubSteppingExtrapolate.h.

Referenced by SubSteppingExtrapolate(), and v_SubStepAdvance().

◆ m_minsubsteps

int Nektar::SubSteppingExtrapolate::m_minsubsteps
protected

Definition at line 88 of file SubSteppingExtrapolate.h.

Referenced by SubSteppingExtrapolate(), and v_SubStepAdvance().

◆ m_previousVelFields

Array<OneD, Array<OneD, NekDouble> > Nektar::SubSteppingExtrapolate::m_previousVelFields
protected

◆ m_previousVnFields

Array<OneD, Array<OneD, NekDouble> > Nektar::SubSteppingExtrapolate::m_previousVnFields
protected

◆ m_subStepIntegrationOps

LibUtilities::TimeIntegrationSchemeOperators Nektar::SubSteppingExtrapolate::m_subStepIntegrationOps
protected

Definition at line 81 of file SubSteppingExtrapolate.h.

Referenced by v_SubStepAdvance(), and v_SubSteppingTimeIntegration().

◆ m_subStepIntegrationScheme

LibUtilities::TimeIntegrationSchemeSharedPtr Nektar::SubSteppingExtrapolate::m_subStepIntegrationScheme
protected