Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Public Member Functions | Static Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | List of all members
Nektar::SubSteppingExtrapolate Class Reference

#include <SubSteppingExtrapolate.h>

Inheritance diagram for Nektar::SubSteppingExtrapolate:
Inheritance graph
[legend]
Collaboration diagram for Nektar::SubSteppingExtrapolate:
Collaboration graph
[legend]

Public Member Functions

 SubSteppingExtrapolate (const LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields, MultiRegions::ExpListSharedPtr pPressure, const Array< OneD, int > pVel, const SolverUtils::AdvectionSharedPtr advObject)
 
virtual ~SubSteppingExtrapolate ()
 
- 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 ()
 
void GenerateHOPBCMap (const LibUtilities::SessionReaderSharedPtr &pSsession)
 
void UpdateRobinPrimCoeff (void)
 
void SubSteppingTimeIntegration (const int intMethod, const LibUtilities::TimeIntegrationWrapperSharedPtr &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 LibUtilities::TimeIntegrationSolutionSharedPtr &integrationSoln, 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 AddDuDt (void)
 
void AddVelBC (void)
 
void ExtrapolatePressureHBCs (void)
 
void CopyPressureHBCsToPbndExp (void)
 
Array< OneD, NekDoubleGetMaxStdVelocity (const Array< OneD, Array< OneD, NekDouble > > inarray)
 
void CorrectPressureBCs (const Array< OneD, NekDouble > &pressure)
 
void IProductNormVelocityOnHBC (const Array< OneD, const Array< OneD, NekDouble > > &Vel, Array< OneD, NekDouble > &IprodVn)
 
void IProductNormVelocityBCOnHBC (Array< OneD, NekDouble > &IprodVn)
 
LibUtilities::TimeIntegrationMethod GetSubStepIntegrationMethod (void)
 
void ExtrapolateArray (Array< OneD, Array< OneD, NekDouble > > &array)
 
void EvaluateBDFArray (Array< OneD, Array< OneD, NekDouble > > &array)
 
void AccelerationBDF (Array< OneD, Array< OneD, NekDouble > > &array)
 
void ExtrapolateArray (Array< OneD, Array< OneD, NekDouble > > &oldarrays, Array< OneD, NekDouble > &newarray, Array< OneD, NekDouble > &outarray)
 
void AddNormVelOnOBC (const int nbcoeffs, const int nreg, Array< OneD, Array< OneD, NekDouble > > &u)
 
void AddPressureToOutflowBCs (NekDouble kinvis)
 

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

virtual void v_EvaluatePressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis)
 
virtual void v_SubSteppingTimeIntegration (int intMethod, const LibUtilities::TimeIntegrationWrapperSharedPtr &IntegrationScheme)
 
virtual void v_SubStepSaveFields (int nstep)
 
virtual void v_SubStepSetPressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &inarray, NekDouble Aii_Dt, NekDouble kinvis)
 
virtual void v_SubStepAdvance (const LibUtilities::TimeIntegrationSolutionSharedPtr &integrationSoln, int nstep, NekDouble time)
 
virtual void v_MountHOPBCs (int HBCdata, NekDouble kinvis, Array< OneD, NekDouble > &Q, Array< OneD, const NekDouble > &Advection)
 
virtual
LibUtilities::TimeIntegrationMethod 
v_GetSubStepIntegrationMethod (void)
 
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)
 
void AddAdvectionPenaltyFlux (const Array< OneD, const Array< OneD, NekDouble > > &velfield, const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
NekDouble GetSubstepTimeStep ()
 
- Protected Member Functions inherited from Nektar::Extrapolate
void CalcNeumannPressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis)
 
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::TimeIntegrationWrapperSharedPtr m_subStepIntegrationScheme
 
LibUtilities::TimeIntegrationSchemeOperators m_subStepIntegrationOps
 
Array< OneD, Array< OneD,
NekDouble > > 
m_previousVelFields
 
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::ExpListSharedPtr
m_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
 
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::BoundaryConditionShPtr
m_PBndConds
 pressure boundary conditions container More...
 
Array< OneD,
MultiRegions::ExpListSharedPtr
m_PBndExp
 pressure boundary conditions expansion container 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
 

Additional Inherited Members

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

Detailed Description

Definition at line 58 of file SubSteppingExtrapolate.h.

Constructor & Destructor Documentation

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

Definition at line 51 of file SubSteppingExtrapolate.cpp.

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

57  : Extrapolate(pSession,pFields,pPressure,pVel,advObject)
58  {
59  m_session->LoadParameter("IO_InfoSteps", m_infosteps, 0);
60  m_session->LoadParameter("SubStepCFL", m_cflSafetyFactor, 0.5);
61  m_session->LoadParameter("MinSubSteps", m_minsubsteps,1);
62  m_session->LoadParameter("MaxSubSteps", m_maxsubsteps,100);
63 
64  int dim = m_fields[0]->GetCoordim(0);
65  m_traceNormals = Array<OneD, Array<OneD, NekDouble> >(dim);
66 
67  int nTracePts = m_fields[0]->GetTrace()->GetNpoints();
68  for(int i = 0; i < dim; ++i)
69  {
70  m_traceNormals[i] = Array<OneD, NekDouble> (nTracePts);
71  }
72  m_fields[0]->GetTrace()->GetNormals(m_traceNormals);
73 
74  }
LibUtilities::SessionReaderSharedPtr m_session
Definition: Extrapolate.h:211
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Velocity fields.
Definition: Extrapolate.h:219
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:59
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Definition: Extrapolate.h:264
Nektar::SubSteppingExtrapolate::~SubSteppingExtrapolate ( )
virtual

Definition at line 76 of file SubSteppingExtrapolate.cpp.

77  {
78  }

Member Function Documentation

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

Number of trace points

Number of spatial dimensions

Forward state array

Backward state array

upwind numerical flux state array

Normal velocity array

Extract forwards/backwards trace spaces Note: Needs to have correct i value to get 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 416 of file SubSteppingExtrapolate.cpp.

References ASSERTL1, Nektar::Extrapolate::m_bnd_dim, Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_traceNormals, Vmath::Vmul(), Vmath::Vsub(), and Vmath::Vvtvp().

Referenced by SubStepAdvection().

420  {
421  ASSERTL1(
422  physfield.num_elements() == Outarray.num_elements(),
423  "Physfield and outarray are of different dimensions");
424 
425  int i;
426 
427  /// Number of trace points
428  int nTracePts = m_fields[0]->GetTrace()->GetNpoints();
429 
430  /// Number of spatial dimensions
431  int nDimensions = m_bnd_dim;
432 
433  /// Forward state array
434  Array<OneD, NekDouble> Fwd(3*nTracePts);
435 
436  /// Backward state array
437  Array<OneD, NekDouble> Bwd = Fwd + nTracePts;
438 
439  /// upwind numerical flux state array
440  Array<OneD, NekDouble> numflux = Bwd + nTracePts;
441 
442  /// Normal velocity array
443  Array<OneD, NekDouble> Vn (nTracePts, 0.0);
444 
445  // Extract velocity field along the trace space and multiply by trace normals
446  for(i = 0; i < nDimensions; ++i)
447  {
448  m_fields[0]->ExtractTracePhys(velfield[i], Fwd);
449  Vmath::Vvtvp(nTracePts, m_traceNormals[i], 1, Fwd, 1, Vn, 1, Vn, 1);
450  }
451 
452  for(i = 0; i < physfield.num_elements(); ++i)
453  {
454  /// Extract forwards/backwards trace spaces
455  /// Note: Needs to have correct i value to get boundary conditions
456  m_fields[i]->GetFwdBwdTracePhys(physfield[i], Fwd, Bwd);
457 
458  /// Upwind between elements
459  m_fields[0]->GetTrace()->Upwind(Vn, Fwd, Bwd, numflux);
460 
461  /// Construct difference between numflux and Fwd,Bwd
462  Vmath::Vsub(nTracePts, numflux, 1, Fwd, 1, Fwd, 1);
463  Vmath::Vsub(nTracePts, numflux, 1, Bwd, 1, Bwd, 1);
464 
465  /// Calculate the numerical fluxes multipling Fwd, Bwd and
466  /// numflux by the normal advection velocity
467  Vmath::Vmul(nTracePts, Fwd, 1, Vn, 1, Fwd, 1);
468  Vmath::Vmul(nTracePts, Bwd, 1, Vn, 1, Bwd, 1);
469 
470  m_fields[0]->AddFwdBwdTraceIntegral(Fwd,Bwd,Outarray[i]);
471  }
472  }
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.cpp:442
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Velocity fields.
Definition: Extrapolate.h:219
int m_bnd_dim
bounday dimensionality
Definition: Extrapolate.h:236
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.cpp:343
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Definition: Extrapolate.h:264
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
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.cpp:183
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.

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

69  {
71  ::AllocateSharedPtr(pSession,pFields,pPressure,pVel,advObject);
72  return p;
73  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
boost::shared_ptr< Extrapolate > ExtrapolateSharedPtr
Definition: Extrapolate.h:60
NekDouble Nektar::SubSteppingExtrapolate::GetSubstepTimeStep ( )
protected

Definition at line 382 of file SubSteppingExtrapolate.cpp.

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

Referenced by v_SubStepAdvance().

383  {
384  int n_element = m_fields[0]->GetExpSize();
385 
386  const Array<OneD, int> ExpOrder=m_fields[0]->EvalBasisNumModesMaxPerExp();
387 
388  const NekDouble cLambda = 0.2; // Spencer book pag. 317
389 
390  Array<OneD, NekDouble> tstep (n_element, 0.0);
391  Array<OneD, NekDouble> stdVelocity(n_element, 0.0);
392  Array<OneD, Array<OneD, NekDouble> > velfields(m_velocity.num_elements());
393 
394  for(int i = 0; i < m_velocity.num_elements(); ++i)
395  {
396  velfields[i] = m_fields[m_velocity[i]]->UpdatePhys();
397  }
398  stdVelocity = GetMaxStdVelocity(velfields);
399 
400  for(int el = 0; el < n_element; ++el)
401  {
402  tstep[el] = m_cflSafetyFactor /
403  (stdVelocity[el] * cLambda *
404  (ExpOrder[el]-1) * (ExpOrder[el]-1));
405  }
406 
407  NekDouble TimeStep = Vmath::Vmin(n_element, tstep, 1);
408  m_comm->AllReduce(TimeStep,LibUtilities::ReduceMin);
409 
410  return TimeStep;
411  }
Array< OneD, int > m_velocity
int which identifies which components of m_fields contains the velocity (u,v,w);
Definition: Extrapolate.h:226
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.cpp:871
LibUtilities::CommSharedPtr m_comm
Definition: Extrapolate.h:213
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Velocity fields.
Definition: Extrapolate.h:219
Array< OneD, NekDouble > GetMaxStdVelocity(const Array< OneD, Array< OneD, NekDouble > > inarray)
double NekDouble
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 161 of file SubSteppingExtrapolate.cpp.

References AddAdvectionPenaltyFlux(), Nektar::Extrapolate::m_advObject, Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_timestep, Nektar::Extrapolate::m_velocity, Vmath::Neg(), and SubStepExtrapolateField().

Referenced by v_SubSteppingTimeIntegration().

165  {
166  int i;
167  int nVariables = inarray.num_elements();
168  int nQuadraturePts = inarray[0].num_elements();
169 
170  /// Get the number of coefficients
171  int ncoeffs = m_fields[0]->GetNcoeffs();
172 
173  /// Define an auxiliary variable to compute the RHS
174  Array<OneD, Array<OneD, NekDouble> > WeakAdv(nVariables);
175  WeakAdv[0] = Array<OneD, NekDouble> (ncoeffs*nVariables);
176  for(i = 1; i < nVariables; ++i)
177  {
178  WeakAdv[i] = WeakAdv[i-1] + ncoeffs;
179  }
180 
181  Array<OneD, Array<OneD, NekDouble> > Velfields(m_velocity.num_elements());
182 
183  Velfields[0] = Array<OneD, NekDouble> (nQuadraturePts*m_velocity.num_elements());
184 
185  for(i = 1; i < m_velocity.num_elements(); ++i)
186  {
187  Velfields[i] = Velfields[i-1] + nQuadraturePts;
188  }
189 
190  SubStepExtrapolateField(fmod(time,m_timestep), Velfields);
191 
192  m_advObject->Advect(m_velocity.num_elements(), m_fields, Velfields, inarray, outarray, time);
193 
194  for(i = 0; i < nVariables; ++i)
195  {
196  m_fields[i]->IProductWRTBase(outarray[i],WeakAdv[i]);
197  // negation requried due to sign of DoAdvection term to be consistent
198  Vmath::Neg(ncoeffs, WeakAdv[i], 1);
199  }
200 
201  AddAdvectionPenaltyFlux(Velfields, inarray, WeakAdv);
202 
203  /// Operations to compute the RHS
204  for(i = 0; i < nVariables; ++i)
205  {
206  // Negate the RHS
207  Vmath::Neg(ncoeffs, WeakAdv[i], 1);
208 
209  /// Multiply the flux by the inverse of the mass matrix
210  m_fields[i]->MultiplyByElmtInvMass(WeakAdv[i], WeakAdv[i]);
211 
212  /// Store in outarray the physical values of the RHS
213  m_fields[i]->BwdTrans(WeakAdv[i], outarray[i]);
214  }
215  }
void SubStepExtrapolateField(NekDouble toff, Array< OneD, Array< OneD, NekDouble > > &ExtVel)
Array< OneD, int > m_velocity
int which identifies which components of m_fields contains the velocity (u,v,w);
Definition: Extrapolate.h:226
void AddAdvectionPenaltyFlux(const Array< OneD, const Array< OneD, NekDouble > > &velfield, const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &outarray)
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Velocity fields.
Definition: Extrapolate.h:219
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.cpp:396
SolverUtils::AdvectionSharedPtr m_advObject
Definition: Extrapolate.h:228
NekDouble m_timestep
Definition: Extrapolate.h:256
void Nektar::SubSteppingExtrapolate::SubStepExtrapolateField ( NekDouble  toff,
Array< OneD, Array< OneD, NekDouble > > &  ExtVel 
)
protected

Extrapolate field using equally time spaced field un,un-1,un-2, (at dt intervals) to time n+t at order Ord

Definition at line 478 of file SubSteppingExtrapolate.cpp.

References Vmath::Fill(), Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_intSteps, m_previousVelFields, Nektar::Extrapolate::m_timestep, Nektar::Extrapolate::m_velocity, npts, and Vmath::Smul().

Referenced by SubStepAdvection().

481  {
482  int npts = m_fields[0]->GetTotPoints();
483  int nvel = m_velocity.num_elements();
484  int i,j;
485  Array<OneD, NekDouble> l(4);
486 
487  int ord = m_intSteps;
488 
489  // calculate Lagrange interpolants
490  Vmath::Fill(4,1.0,l,1);
491 
492  for(i = 0; i <= ord; ++i)
493  {
494  for(j = 0; j <= ord; ++j)
495  {
496  if(i != j)
497  {
498  l[i] *= (j*m_timestep+toff);
499  l[i] /= (j*m_timestep-i*m_timestep);
500  }
501  }
502  }
503 
504  for(i = 0; i < nvel; ++i)
505  {
506  Vmath::Smul(npts,l[0],m_previousVelFields[i],1,ExtVel[i],1);
507 
508  for(j = 1; j <= ord; ++j)
509  {
510  Blas::Daxpy(npts,l[j],m_previousVelFields[j*nvel+i],1,
511  ExtVel[i],1);
512  }
513  }
514  }
Array< OneD, int > m_velocity
int which identifies which components of m_fields contains the velocity (u,v,w);
Definition: Extrapolate.h:226
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
Definition: Vmath.cpp:46
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Velocity fields.
Definition: Extrapolate.h:219
Array< OneD, Array< OneD, NekDouble > > m_previousVelFields
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
Definition: Vmath.cpp:213
int m_intSteps
Maximum points used in pressure BC evaluation.
Definition: Extrapolate.h:254
static std::string npts
Definition: InputFld.cpp:43
NekDouble m_timestep
Definition: Extrapolate.h:256
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 220 of file SubSteppingExtrapolate.cpp.

References ASSERTL1, and Vmath::Vcopy().

Referenced by v_SubSteppingTimeIntegration().

224  {
225  ASSERTL1(inarray.num_elements() == outarray.num_elements(),"Inarray and outarray of different sizes ");
226 
227  for(int i = 0; i < inarray.num_elements(); ++i)
228  {
229  Vmath::Vcopy(inarray[i].num_elements(),inarray[i],1,outarray[i],1);
230  }
231  }
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
void Nektar::SubSteppingExtrapolate::v_EvaluatePressureBCs ( const Array< OneD, const Array< OneD, NekDouble > > &  fields,
const Array< OneD, const Array< OneD, NekDouble > > &  N,
NekDouble  kinvis 
)
protectedvirtual

Implements Nektar::Extrapolate.

Definition at line 80 of file SubSteppingExtrapolate.cpp.

References ASSERTL0.

81  {
82  ASSERTL0(false,"This method should not be called by Substepping routine");
83  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
LibUtilities::TimeIntegrationMethod Nektar::SubSteppingExtrapolate::v_GetSubStepIntegrationMethod ( void  )
protectedvirtual

Reimplemented from Nektar::Extrapolate.

Definition at line 528 of file SubSteppingExtrapolate.cpp.

References m_subStepIntegrationScheme.

529  {
530  return m_subStepIntegrationScheme->GetIntegrationMethod();
531  }
LibUtilities::TimeIntegrationWrapperSharedPtr m_subStepIntegrationScheme
void Nektar::SubSteppingExtrapolate::v_MountHOPBCs ( int  HBCdata,
NekDouble  kinvis,
Array< OneD, NekDouble > &  Q,
Array< OneD, const NekDouble > &  Advection 
)
protectedvirtual

Implements Nektar::Extrapolate.

Definition at line 519 of file SubSteppingExtrapolate.cpp.

References Vmath::Smul().

524  {
525  Vmath::Smul(HBCdata,-kinvis,Q,1,Q,1);
526  }
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
Definition: Vmath.cpp:213
void Nektar::SubSteppingExtrapolate::v_SubStepAdvance ( const LibUtilities::TimeIntegrationSolutionSharedPtr integrationSoln,
int  nstep,
NekDouble  time 
)
protectedvirtual

Implements Nektar::Extrapolate.

Definition at line 317 of file SubSteppingExtrapolate.cpp.

References ASSERTL0, GetSubstepTimeStep(), Nektar::Extrapolate::IProductNormVelocityOnHBC(), m_cflSafetyFactor, Nektar::Extrapolate::m_comm, Nektar::Extrapolate::m_fields, m_infosteps, Nektar::Extrapolate::m_intSteps, Nektar::Extrapolate::m_iprodnormvel, m_maxsubsteps, m_minsubsteps, m_subStepIntegrationOps, m_subStepIntegrationScheme, and Nektar::Extrapolate::m_timestep.

321  {
322  int n;
323  int nsubsteps;
324 
325  NekDouble dt;
326 
327  Array<OneD, Array<OneD, NekDouble> > fields;
328 
329  static int ncalls = 1;
330  int nint = min(ncalls++, m_intSteps);
331 
332  //this needs to change
333  m_comm = m_fields[0]->GetComm()->GetRowComm();
334 
335  // Get the proper time step with CFL control
336  dt = GetSubstepTimeStep();
337 
338  nsubsteps = (m_timestep > dt)? ((int)(m_timestep/dt)+1):1;
339  nsubsteps = max(m_minsubsteps, nsubsteps);
340 
341  ASSERTL0(nsubsteps < m_maxsubsteps,"Number of substeps has exceeded maximum");
342 
343  dt = m_timestep/nsubsteps;
344 
345  if (m_infosteps && !((nstep+1)%m_infosteps) && m_comm->GetRank() == 0)
346  {
347  cout << "Sub-integrating using "<< nsubsteps
348  << " steps over Dt = " << m_timestep
349  << " (SubStep CFL=" << m_cflSafetyFactor << ")"<< endl;
350  }
351 
352  for (int m = 0; m < nint; ++m)
353  {
354  // We need to update the fields held by the m_integrationSoln
355  fields = integrationSoln->UpdateSolutionVector()[m];
356 
357  // Initialise NS solver which is set up to use a GLM method
358  // with calls to EvaluateAdvection_SetPressureBCs and
359  // SolveUnsteadyStokesSystem
361  SubIntegrationSoln = m_subStepIntegrationScheme->
362  InitializeScheme(dt, fields, time, m_subStepIntegrationOps);
363 
364  for(n = 0; n < nsubsteps; ++n)
365  {
366  fields = m_subStepIntegrationScheme->TimeIntegrate(n, dt, SubIntegrationSoln,
368  }
369 
370  // set up HBC m_acceleration field for Pressure BCs
372 
373  // Reset time integrated solution in m_integrationSoln
374  integrationSoln->SetSolVector(m,fields);
375  }
376  }
LibUtilities::TimeIntegrationWrapperSharedPtr m_subStepIntegrationScheme
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
void IProductNormVelocityOnHBC(const Array< OneD, const Array< OneD, NekDouble > > &Vel, Array< OneD, NekDouble > &IprodVn)
LibUtilities::CommSharedPtr m_comm
Definition: Extrapolate.h:213
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Velocity fields.
Definition: Extrapolate.h:219
int m_intSteps
Maximum points used in pressure BC evaluation.
Definition: Extrapolate.h:254
double NekDouble
boost::shared_ptr< TimeIntegrationSolution > TimeIntegrationSolutionSharedPtr
Array< OneD, Array< OneD, NekDouble > > m_iprodnormvel
Storage for current and previous levels of the inner product of normal velocity.
Definition: Extrapolate.h:262
NekDouble m_timestep
Definition: Extrapolate.h:256
LibUtilities::TimeIntegrationSchemeOperators m_subStepIntegrationOps
void Nektar::SubSteppingExtrapolate::v_SubSteppingTimeIntegration ( int  intMethod,
const LibUtilities::TimeIntegrationWrapperSharedPtr IntegrationScheme 
)
protectedvirtual

Implements Nektar::Extrapolate.

Definition at line 86 of file SubSteppingExtrapolate.cpp.

References ASSERTL0, Nektar::LibUtilities::NekFactory< tKey, tBase, >::CreateInstance(), Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineOdeRhs(), Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineProjection(), Nektar::LibUtilities::eBackwardEuler, Nektar::LibUtilities::eBDFImplicitOrder1, Nektar::LibUtilities::eBDFImplicitOrder2, Nektar::LibUtilities::GetTimeIntegrationWrapperFactory(), Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_intSteps, m_previousVelFields, Nektar::Extrapolate::m_session, m_subStepIntegrationOps, m_subStepIntegrationScheme, Nektar::Extrapolate::m_velocity, SubStepAdvection(), and SubStepProjection().

89  {
90  int i;
91 
92  // Set to 1 for first step and it will then be increased in
93  // time advance routines
94  switch(intMethod)
95  {
98  {
99  std::string vSubStepIntScheme = "ForwardEuler";
100 
101  if(m_session->DefinesSolverInfo("SubStepIntScheme"))
102  {
103  vSubStepIntScheme = m_session->GetSolverInfo("SubStepIntScheme");
104  }
105 
107 
108  int nvel = m_velocity.num_elements();
109 
110  // Fields for linear interpolation
111  m_previousVelFields = Array<OneD, Array<OneD, NekDouble> >(2*nvel);
112  int ntotpts = m_fields[0]->GetTotPoints();
113  m_previousVelFields[0] = Array<OneD, NekDouble>(2*nvel*ntotpts);
114  for(i = 1; i < 2*nvel; ++i)
115  {
116  m_previousVelFields[i] = m_previousVelFields[i-1] + ntotpts;
117  }
118 
119  }
120  break;
122  {
123  std::string vSubStepIntScheme = "RungeKutta2_ImprovedEuler";
124 
125  if(m_session->DefinesSolverInfo("SubStepIntScheme"))
126  {
127  vSubStepIntScheme = m_session->GetSolverInfo("SubStepIntScheme");
128  }
129 
131 
132  int nvel = m_velocity.num_elements();
133 
134  // Fields for quadratic interpolation
135  m_previousVelFields = Array<OneD, Array<OneD, NekDouble> >(3*nvel);
136 
137  int ntotpts = m_fields[0]->GetTotPoints();
138  m_previousVelFields[0] = Array<OneD, NekDouble>(3*nvel*ntotpts);
139  for(i = 1; i < 3*nvel; ++i)
140  {
141  m_previousVelFields[i] = m_previousVelFields[i-1] + ntotpts;
142  }
143 
144  }
145  break;
146  default:
147  ASSERTL0(0,"Integration method not suitable: Options include BackwardEuler or BDFImplicitOrder1");
148  break;
149  }
150 
151  m_intSteps = IntegrationScheme->GetIntegrationSteps();
152 
153  // set explicit time-integration class operators
156  }
void SubStepProjection(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
LibUtilities::TimeIntegrationWrapperSharedPtr m_subStepIntegrationScheme
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
LibUtilities::SessionReaderSharedPtr m_session
Definition: Extrapolate.h:211
Array< OneD, int > m_velocity
int which identifies which components of m_fields contains the velocity (u,v,w);
Definition: Extrapolate.h:226
BDF multi-step scheme of order 1 (implicit)
tBaseSharedPtr CreateInstance(tKey idKey BOOST_PP_COMMA_IF(MAX_PARAM) BOOST_PP_ENUM_BINARY_PARAMS(MAX_PARAM, tParam, x))
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:162
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Velocity fields.
Definition: Extrapolate.h:219
Array< OneD, Array< OneD, NekDouble > > m_previousVelFields
void DefineProjection(FuncPointerT func, ObjectPointerT obj)
void DefineOdeRhs(FuncPointerT func, ObjectPointerT obj)
int m_intSteps
Maximum points used in pressure BC evaluation.
Definition: Extrapolate.h:254
TimeIntegrationWrapperFactory & GetTimeIntegrationWrapperFactory()
void SubStepAdvection(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
BDF multi-step scheme of order 2 (implicit)
LibUtilities::TimeIntegrationSchemeOperators m_subStepIntegrationOps
void Nektar::SubSteppingExtrapolate::v_SubStepSaveFields ( int  nstep)
protectedvirtual

Implements Nektar::Extrapolate.

Definition at line 268 of file SubSteppingExtrapolate.cpp.

References Nektar::Extrapolate::m_fields, m_previousVelFields, Nektar::Extrapolate::m_velocity, npts, and Vmath::Vcopy().

269  {
270  int i,n;
271  int nvel = m_velocity.num_elements();
272  int npts = m_fields[0]->GetTotPoints();
273 
274  // rotate fields
275  int nblocks = m_previousVelFields.num_elements()/nvel;
276  Array<OneD, NekDouble> save;
277 
278  // rotate storage space
279  for(n = 0; n < nvel; ++n)
280  {
281  save = m_previousVelFields[(nblocks-1)*nvel+n];
282 
283  for(i = nblocks-1; i > 0; --i)
284  {
285  m_previousVelFields[i*nvel+n] = m_previousVelFields[(i-1)*nvel+n];
286  }
287 
288  m_previousVelFields[n] = save;
289  }
290 
291  // Put previous field
292  for(i = 0; i < nvel; ++i)
293  {
294  m_fields[m_velocity[i]]->BwdTrans(m_fields[m_velocity[i]]->GetCoeffs(),
295  m_fields[m_velocity[i]]->UpdatePhys());
296  Vmath::Vcopy(npts,m_fields[m_velocity[i]]->GetPhys(),1,
297  m_previousVelFields[i],1);
298  }
299 
300  if(nstep == 0)// initialise all levels with first field
301  {
302  for(n = 0; n < nvel; ++n)
303  {
304  for(i = 1; i < nblocks; ++i)
305  {
306  Vmath::Vcopy(npts,m_fields[m_velocity[n]]->GetPhys(),1,
307  m_previousVelFields[i*nvel+n],1);
308 
309  }
310  }
311  }
312  }
Array< OneD, int > m_velocity
int which identifies which components of m_fields contains the velocity (u,v,w);
Definition: Extrapolate.h:226
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Velocity fields.
Definition: Extrapolate.h:219
Array< OneD, Array< OneD, NekDouble > > m_previousVelFields
static std::string npts
Definition: InputFld.cpp:43
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
void Nektar::SubSteppingExtrapolate::v_SubStepSetPressureBCs ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
NekDouble  Aii_Dt,
NekDouble  kinvis 
)
protectedvirtual

Implements Nektar::Extrapolate.

Reimplemented in Nektar::SubSteppingExtrapolateWeakPressure.

Definition at line 237 of file SubSteppingExtrapolate.cpp.

References Nektar::Extrapolate::AddDuDt(), Nektar::Extrapolate::CalcNeumannPressureBCs(), Nektar::Extrapolate::CalcOutflowBCs(), Nektar::Extrapolate::CopyPressureHBCsToPbndExp(), Nektar::Extrapolate::ExtrapolateArray(), Nektar::Extrapolate::m_pressureCalls, and Nektar::Extrapolate::m_pressureHBCs.

241  {
242  //int nConvectiveFields =m_fields.num_elements()-1;
243  Array<OneD, Array<OneD, NekDouble> > nullvelfields;
244 
245  m_pressureCalls++;
246 
247  // Calculate non-linear and viscous BCs at current level and
248  // put in m_pressureHBCs[0]
249  CalcNeumannPressureBCs(inarray,nullvelfields,kinvis);
250 
251  // Extrapolate to m_pressureHBCs to n+1
253 
254  // Add (phi,Du/Dt) term to m_presureHBC
255  AddDuDt();
256 
257  // Copy m_pressureHBCs to m_PbndExp
259 
260  // Evaluate High order outflow conditiosn if required.
261  CalcOutflowBCs(inarray, kinvis);
262  }
void ExtrapolateArray(Array< OneD, Array< OneD, NekDouble > > &array)
void CopyPressureHBCsToPbndExp(void)
int m_pressureCalls
number of times the high-order pressure BCs have been called
Definition: Extrapolate.h:245
Array< OneD, Array< OneD, NekDouble > > m_pressureHBCs
Storage for current and previous levels of high order pressure boundary conditions.
Definition: Extrapolate.h:259
void CalcOutflowBCs(const Array< OneD, const Array< OneD, NekDouble > > &fields, NekDouble kinvis)
void CalcNeumannPressureBCs(const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis)
Definition: Extrapolate.h:188

Member Data Documentation

std::string Nektar::SubSteppingExtrapolate::className
static
Initial value:

Name of class.

Registers the class with the Factory.

Definition at line 76 of file SubSteppingExtrapolate.h.

NekDouble Nektar::SubSteppingExtrapolate::m_cflSafetyFactor
protected
int Nektar::SubSteppingExtrapolate::m_infosteps
protected

Definition at line 145 of file SubSteppingExtrapolate.h.

Referenced by SubSteppingExtrapolate(), and v_SubStepAdvance().

int Nektar::SubSteppingExtrapolate::m_maxsubsteps
protected

Definition at line 147 of file SubSteppingExtrapolate.h.

Referenced by SubSteppingExtrapolate(), and v_SubStepAdvance().

int Nektar::SubSteppingExtrapolate::m_minsubsteps
protected

Definition at line 146 of file SubSteppingExtrapolate.h.

Referenced by SubSteppingExtrapolate(), and v_SubStepAdvance().

Array<OneD, Array<OneD, NekDouble> > Nektar::SubSteppingExtrapolate::m_previousVelFields
protected
LibUtilities::TimeIntegrationSchemeOperators Nektar::SubSteppingExtrapolate::m_subStepIntegrationOps
protected

Definition at line 140 of file SubSteppingExtrapolate.h.

Referenced by v_SubStepAdvance(), and v_SubSteppingTimeIntegration().

LibUtilities::TimeIntegrationWrapperSharedPtr Nektar::SubSteppingExtrapolate::m_subStepIntegrationScheme
protected