Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
SubSteppingExtrapolate.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File: SubSteppingExtrapolate.cpp
4 //a
5 // For more information, please see: http://www.nektar.info
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description: Abstract base class for SubSteppingExtrapolate.
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
38 
39 using namespace std;
40 
41 namespace Nektar
42 {
43  /**
44  * Registers the class with the Factory.
45  */
46  std::string SubSteppingExtrapolate::className = GetExtrapolateFactory().RegisterCreatorFunction(
47  "SubStepping",
48  SubSteppingExtrapolate::create,
49  "SubStepping");
50 
51  SubSteppingExtrapolate::SubSteppingExtrapolate(
55  const Array<OneD, int> pVel,
56  const SolverUtils::AdvectionSharedPtr advObject)
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);
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  }
75 
77  {
78  }
79 
81  {
82  ASSERTL0(false,"This method should not be called by Substepping routine");
83  }
84 
85 
87  int intMethod,
88  const LibUtilities::TimeIntegrationWrapperSharedPtr &IntegrationScheme)
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
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
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  }
157 
158  /**
159  * Explicit Advection terms used by SubStepAdvance time integration
160  */
162  const Array<OneD, const Array<OneD, NekDouble> > &inarray,
163  Array<OneD, Array<OneD, NekDouble> > &outarray,
164  const NekDouble time)
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  Array<OneD, int> VelIds(m_velocity.num_elements());
183 
184  Velfields[0] = Array<OneD, NekDouble> (nQuadraturePts*m_velocity.num_elements());
185 
186  for(i = 1; i < m_velocity.num_elements(); ++i)
187  {
188  Velfields[i] = Velfields[i-1] + nQuadraturePts;
189  }
190 
191  SubStepExtrapolateField(fmod(time,m_timestep), Velfields);
192 
193  m_advObject->Advect(m_velocity.num_elements(), m_fields, Velfields, inarray, outarray, time);
194 
195  for(i = 0; i < nVariables; ++i)
196  {
197  m_fields[i]->IProductWRTBase(outarray[i],WeakAdv[i]);
198  // negation requried due to sign of DoAdvection term to be consistent
199  Vmath::Neg(ncoeffs, WeakAdv[i], 1);
200  }
201 
202  AddAdvectionPenaltyFlux(Velfields, inarray, WeakAdv);
203 
204  /// Operations to compute the RHS
205  for(i = 0; i < nVariables; ++i)
206  {
207  // Negate the RHS
208  Vmath::Neg(ncoeffs, WeakAdv[i], 1);
209 
210  /// Multiply the flux by the inverse of the mass matrix
211  m_fields[i]->MultiplyByElmtInvMass(WeakAdv[i], WeakAdv[i]);
212 
213  /// Store in outarray the physical values of the RHS
214  m_fields[i]->BwdTrans(WeakAdv[i], outarray[i]);
215  }
216  }
217 
218  /**
219  * Projection used by SubStepAdvance time integration
220  */
222  const Array<OneD, const Array<OneD, NekDouble> > &inarray,
223  Array<OneD, Array<OneD, NekDouble> > &outarray,
224  const NekDouble time)
225  {
226  ASSERTL1(inarray.num_elements() == outarray.num_elements(),"Inarray and outarray of different sizes ");
227 
228  for(int i = 0; i < inarray.num_elements(); ++i)
229  {
230  Vmath::Vcopy(inarray[i].num_elements(),inarray[i],1,outarray[i],1);
231  }
232  }
233 
234 
235  /**
236  *
237  */
239  const Array<OneD, const Array<OneD, NekDouble> > &inarray,
240  const NekDouble Aii_Dt,
241  NekDouble kinvis)
242  {
243  int nConvectiveFields =m_fields.num_elements()-1;
244  Array<OneD, Array<OneD, NekDouble> > velfields(nConvectiveFields);
245 
246  for(int i = 0; i < nConvectiveFields; ++i)
247  {
248  velfields[i] = m_fields[m_velocity[i]]->GetPhys();
249  }
250 
251  m_pressureCalls++;
252 
253  // Rotate HOPBCs storage
255 
256  // Calculate non-linear and viscous BCs at current level and
257  // put in m_pressureHBCs[0]
258  CalcNeumannPressureBCs(inarray,velfields,kinvis);
259 
260  // Extrapolate to m_pressureHBCs to n+1
262 
263  // Add (phi,Du/Dt) term to m_presureHBC
264  AddDuDt();
265 
266  // Copy m_pressureHBCs to m_PbndExp
268 
269  // Evaluate High order outflow conditiosn if required.
270  CalcOutflowBCs(inarray, kinvis);
271  }
272 
273  /**
274  *
275  */
277  {
278  int i,n;
279  int nvel = m_velocity.num_elements();
280  int npts = m_fields[0]->GetTotPoints();
281 
282  // rotate fields
283  int nblocks = m_previousVelFields.num_elements()/nvel;
285 
286  // rotate storage space
287  for(n = 0; n < nvel; ++n)
288  {
289  save = m_previousVelFields[(nblocks-1)*nvel+n];
290 
291  for(i = nblocks-1; i > 0; --i)
292  {
293  m_previousVelFields[i*nvel+n] = m_previousVelFields[(i-1)*nvel+n];
294  }
295 
296  m_previousVelFields[n] = save;
297  }
298 
299  // Put previous field
300  for(i = 0; i < nvel; ++i)
301  {
302  m_fields[m_velocity[i]]->BwdTrans(m_fields[m_velocity[i]]->GetCoeffs(),
303  m_fields[m_velocity[i]]->UpdatePhys());
304  Vmath::Vcopy(npts,m_fields[m_velocity[i]]->GetPhys(),1,
305  m_previousVelFields[i],1);
306  }
307 
308  if(nstep == 0)// initialise all levels with first field
309  {
310  for(n = 0; n < nvel; ++n)
311  {
312  for(i = 1; i < nblocks; ++i)
313  {
314  Vmath::Vcopy(npts,m_fields[m_velocity[n]]->GetPhys(),1,
315  m_previousVelFields[i*nvel+n],1);
316 
317  }
318  }
319  }
320  }
321 
322  /**
323  *
324  */
326  const LibUtilities::TimeIntegrationSolutionSharedPtr &integrationSoln,
327  int nstep,
328  NekDouble time)
329  {
330  int n;
331  int nsubsteps;
332 
333  NekDouble dt;
334 
335  Array<OneD, Array<OneD, NekDouble> > fields, velfields;
336 
337  static int ncalls = 1;
338  int nint = min(ncalls++, m_intSteps);
339 
340  Array<OneD, NekDouble> CFL(m_fields[0]->GetExpSize(),
342  //this needs to change
343  m_comm = m_fields[0]->GetComm()->GetRowComm();
344 
345  // Get the proper time step with CFL control
346  dt = GetSubstepTimeStep();
347 
348  nsubsteps = (m_timestep > dt)? ((int)(m_timestep/dt)+1):1;
349  nsubsteps = max(m_minsubsteps, nsubsteps);
350 
351  ASSERTL0(nsubsteps < m_maxsubsteps,"Number of substeps has exceeded maximum");
352 
353  dt = m_timestep/nsubsteps;
354 
355  if (m_infosteps && !((nstep+1)%m_infosteps) && m_comm->GetRank() == 0)
356  {
357  cout << "Sub-integrating using "<< nsubsteps
358  << " steps over Dt = " << m_timestep
359  << " (SubStep CFL=" << m_cflSafetyFactor << ")"<< endl;
360  }
361 
362  for (int m = 0; m < nint; ++m)
363  {
364  // We need to update the fields held by the m_integrationSoln
365  fields = integrationSoln->UpdateSolutionVector()[m];
366 
367  // Initialise NS solver which is set up to use a GLM method
368  // with calls to EvaluateAdvection_SetPressureBCs and
369  // SolveUnsteadyStokesSystem
371  SubIntegrationSoln = m_subStepIntegrationScheme->
372  InitializeScheme(dt, fields, time, m_subStepIntegrationOps);
373 
374  for(n = 0; n < nsubsteps; ++n)
375  {
376  fields = m_subStepIntegrationScheme->TimeIntegrate(n, dt, SubIntegrationSoln,
378  }
379 
380  // set up HBC m_acceleration field for Pressure BCs
382 
383  // Reset time integrated solution in m_integrationSoln
384  integrationSoln->SetSolVector(m,fields);
385  }
386  }
387 
388 
389  /**
390  *
391  */
393  {
394  int n_element = m_fields[0]->GetExpSize();
395 
396  const Array<OneD, int> ExpOrder=m_fields[0]->EvalBasisNumModesMaxPerExp();
397  Array<OneD, int> ExpOrderList (n_element, ExpOrder);
398 
399  const NekDouble cLambda = 0.2; // Spencer book pag. 317
400 
401  Array<OneD, NekDouble> tstep (n_element, 0.0);
402  Array<OneD, NekDouble> stdVelocity(n_element, 0.0);
403  Array<OneD, Array<OneD, NekDouble> > velfields(m_velocity.num_elements());
404 
405  for(int i = 0; i < m_velocity.num_elements(); ++i)
406  {
407  velfields[i] = m_fields[m_velocity[i]]->UpdatePhys();
408  }
409  stdVelocity = GetMaxStdVelocity(velfields);
410 
411  for(int el = 0; el < n_element; ++el)
412  {
413  tstep[el] = m_cflSafetyFactor /
414  (stdVelocity[el] * cLambda *
415  (ExpOrder[el]-1) * (ExpOrder[el]-1));
416  }
417 
418  NekDouble TimeStep = Vmath::Vmin(n_element, tstep, 1);
419  m_comm->AllReduce(TimeStep,LibUtilities::ReduceMin);
420 
421  return TimeStep;
422  }
423 
424 
425 
426 
428  const Array<OneD, const Array<OneD, NekDouble> > &velfield,
429  const Array<OneD, const Array<OneD, NekDouble> > &physfield,
430  Array<OneD, Array<OneD, NekDouble> > &Outarray)
431  {
432  ASSERTL1(
433  physfield.num_elements() == Outarray.num_elements(),
434  "Physfield and outarray are of different dimensions");
435 
436  int i;
437 
438  /// Number of trace points
439  int nTracePts = m_fields[0]->GetTrace()->GetNpoints();
440 
441  /// Number of spatial dimensions
442  int nDimensions = m_bnd_dim;
443 
444  /// Forward state array
445  Array<OneD, NekDouble> Fwd(3*nTracePts);
446 
447  /// Backward state array
448  Array<OneD, NekDouble> Bwd = Fwd + nTracePts;
449 
450  /// upwind numerical flux state array
451  Array<OneD, NekDouble> numflux = Bwd + nTracePts;
452 
453  /// Normal velocity array
454  Array<OneD, NekDouble> Vn (nTracePts, 0.0);
455 
456  // Extract velocity field along the trace space and multiply by trace normals
457  for(i = 0; i < nDimensions; ++i)
458  {
459  m_fields[0]->ExtractTracePhys(velfield[i], Fwd);
460  Vmath::Vvtvp(nTracePts, m_traceNormals[i], 1, Fwd, 1, Vn, 1, Vn, 1);
461  }
462 
463  for(i = 0; i < physfield.num_elements(); ++i)
464  {
465  /// Extract forwards/backwards trace spaces
466  /// Note: Needs to have correct i value to get boundary conditions
467  m_fields[i]->GetFwdBwdTracePhys(physfield[i], Fwd, Bwd);
468 
469  /// Upwind between elements
470  m_fields[0]->GetTrace()->Upwind(Vn, Fwd, Bwd, numflux);
471 
472  /// Construct difference between numflux and Fwd,Bwd
473  Vmath::Vsub(nTracePts, numflux, 1, Fwd, 1, Fwd, 1);
474  Vmath::Vsub(nTracePts, numflux, 1, Bwd, 1, Bwd, 1);
475 
476  /// Calculate the numerical fluxes multipling Fwd, Bwd and
477  /// numflux by the normal advection velocity
478  Vmath::Vmul(nTracePts, Fwd, 1, Vn, 1, Fwd, 1);
479  Vmath::Vmul(nTracePts, Bwd, 1, Vn, 1, Bwd, 1);
480 
481  m_fields[0]->AddFwdBwdTraceIntegral(Fwd,Bwd,Outarray[i]);
482  }
483  }
484 
485  /**
486  * Extrapolate field using equally time spaced field un,un-1,un-2, (at
487  * dt intervals) to time n+t at order Ord
488  */
490  NekDouble toff,
491  Array< OneD, Array<OneD, NekDouble> > &ExtVel)
492  {
493  int npts = m_fields[0]->GetTotPoints();
494  int nvel = m_velocity.num_elements();
495  int i,j;
497 
498  int ord = m_intSteps;
499 
500  // calculate Lagrange interpolants
501  Vmath::Fill(4,1.0,l,1);
502 
503  for(i = 0; i <= ord; ++i)
504  {
505  for(j = 0; j <= ord; ++j)
506  {
507  if(i != j)
508  {
509  l[i] *= (j*m_timestep+toff);
510  l[i] /= (j*m_timestep-i*m_timestep);
511  }
512  }
513  }
514 
515  for(i = 0; i < nvel; ++i)
516  {
517  Vmath::Smul(npts,l[0],m_previousVelFields[i],1,ExtVel[i],1);
518 
519  for(j = 1; j <= ord; ++j)
520  {
521  Blas::Daxpy(npts,l[j],m_previousVelFields[j*nvel+i],1,
522  ExtVel[i],1);
523  }
524  }
525  }
526 
527  /**
528  *
529  */
531  int HBCdata,
532  NekDouble kinvis,
535  {
536  Vmath::Smul(HBCdata,-kinvis,Q,1,Q,1);
537  }
538 
539  /**
540  *
541  */
543  {
544  int HBCPts = m_acceleration[0].num_elements();
545 
546  Array<OneD, NekDouble> accelerationTerm(HBCPts, 0.0);
547 
548  // Update velocity BF at n+1 (actually only needs doing if velocity is time dependent on HBCs)
550 
551  //Calculate acceleration term at level n based on previous steps
552  int acc_order = min(m_pressureCalls,m_intSteps);
554  m_acceleration[0], 1, accelerationTerm, 1);
555 
556  for(int i = 0; i < acc_order; i++)
557  {
558  Vmath::Svtvp(HBCPts,
559  -1*StifflyStable_Alpha_Coeffs[acc_order-1][i]/m_timestep,
560  m_acceleration[i+1], 1,
561  accelerationTerm, 1,
562  accelerationTerm, 1);
563  }
564 
565  // Subtract accleration term off m_pressureHBCs[nlevels-1]
566  int nlevels = m_pressureHBCs.num_elements();
567  Vmath::Vsub(HBCPts, m_pressureHBCs[nlevels-1],1,accelerationTerm, 1,
568  m_pressureHBCs[nlevels-1],1);
569  }
570 
572  {
573  return m_subStepIntegrationScheme->GetIntegrationMethod();
574  }
575 
576 }
577 
virtual void v_EvaluatePressureBCs(const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis)
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)
LibUtilities::TimeIntegrationWrapperSharedPtr m_subStepIntegrationScheme
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
void IProductNormVelocityOnHBC(const Array< OneD, const Array< OneD, NekDouble > > &Vel, Array< OneD, NekDouble > &IprodVn)
LibUtilities::SessionReaderSharedPtr m_session
Definition: Extrapolate.h:208
Array< OneD, int > m_velocity
int which identifies which components of m_fields contains the velocity (u,v,w);
Definition: Extrapolate.h:219
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
virtual void v_SubStepSetPressureBCs(const Array< OneD, const Array< OneD, NekDouble > > &inarray, NekDouble Aii_Dt, NekDouble kinvis)
ExtrapolateFactory & GetExtrapolateFactory()
Definition: Extrapolate.cpp:50
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:857
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
Definition: Vmath.cpp:46
boost::shared_ptr< TimeIntegrationWrapper > TimeIntegrationWrapperSharedPtr
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.cpp:471
boost::shared_ptr< Advection > AdvectionSharedPtr
A shared pointer to an Advection object.
Definition: Advection.h:158
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:428
void RollOver(Array< OneD, Array< OneD, NekDouble > > &input)
STL namespace.
void AddAdvectionPenaltyFlux(const Array< OneD, const Array< OneD, NekDouble > > &velfield, const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &outarray)
void CopyPressureHBCsToPbndExp(void)
void ExtrapolatePressureHBCs(void)
LibUtilities::CommSharedPtr m_comm
Definition: Extrapolate.h:210
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Definition: Extrapolate.h:212
static NekDouble StifflyStable_Gamma0_Coeffs[3]
Definition: Extrapolate.h:314
virtual void v_SubStepSaveFields(int nstep)
Array< OneD, Array< OneD, NekDouble > > m_previousVelFields
virtual void v_SubSteppingTimeIntegration(int intMethod, const LibUtilities::TimeIntegrationWrapperSharedPtr &IntegrationScheme)
void DefineProjection(FuncPointerT func, ObjectPointerT obj)
int m_bnd_dim
bounday dimensionality
Definition: Extrapolate.h:229
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:199
Array< OneD, Array< OneD, NekDouble > > m_acceleration
Storage for current and previous levels of the acceleration term.
Definition: Extrapolate.h:274
Array< OneD, NekDouble > GetMaxStdVelocity(const Array< OneD, Array< OneD, NekDouble > > inarray)
void DefineOdeRhs(FuncPointerT func, ObjectPointerT obj)
int m_intSteps
Maximum points used in pressure BC evaluation.
Definition: Extrapolate.h:247
boost::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
static std::string npts
Definition: InputFld.cpp:43
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.cpp:382
virtual void v_SubStepAdvance(const LibUtilities::TimeIntegrationSolutionSharedPtr &integrationSoln, int nstep, NekDouble time)
double NekDouble
TimeIntegrationWrapperFactory & GetTimeIntegrationWrapperFactory()
int m_pressureCalls
number of times the high-order pressure BCs have been called
Definition: Extrapolate.h:238
void SubStepAdvection(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
void IProductNormVelocityBCOnHBC(Array< OneD, NekDouble > &IprodVn)
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:329
virtual void v_MountHOPBCs(int HBCdata, NekDouble kinvis, Array< OneD, NekDouble > &Q, Array< OneD, const NekDouble > &Advection)
BDF multi-step scheme of order 2 (implicit)
Array< OneD, Array< OneD, NekDouble > > m_pressureHBCs
Storage for current and previous levels of high order pressure boundary conditions.
Definition: Extrapolate.h:271
boost::shared_ptr< TimeIntegrationSolution > TimeIntegrationSolutionSharedPtr
static NekDouble StifflyStable_Alpha_Coeffs[3][3]
Definition: Extrapolate.h:313
SolverUtils::AdvectionSharedPtr m_advObject
Definition: Extrapolate.h:221
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Definition: Extrapolate.h:288
virtual LibUtilities::TimeIntegrationMethod v_GetSubStepIntegrationMethod(void)
void CalcOutflowBCs(const Array< OneD, const Array< OneD, NekDouble > > &fields, NekDouble kinvis)
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
NekDouble m_timestep
Definition: Extrapolate.h:249
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:169
Defines a callback function which evaluates the flux vector.
Definition: Advection.h:69
void CalcNeumannPressureBCs(const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis)
Definition: Extrapolate.h:179
LibUtilities::TimeIntegrationSchemeOperators m_subStepIntegrationOps
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, tDescription pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:215