48 SubSteppingExtrapolate::create,
51 SubSteppingExtrapolate::SubSteppingExtrapolate(
57 :
Extrapolate(pSession,pFields,pPressure,pVel,advObject)
64 int dim =
m_fields[0]->GetCoordim(0);
67 int nTracePts =
m_fields[0]->GetTrace()->GetNpoints();
68 for(
int i = 0; i < dim; ++i)
82 ASSERTL0(
false,
"This method should not be called by Substepping routine");
99 std::string vSubStepIntScheme =
"ForwardEuler";
101 if(
m_session->DefinesSolverInfo(
"SubStepIntScheme"))
103 vSubStepIntScheme =
m_session->GetSolverInfo(
"SubStepIntScheme");
112 int ntotpts =
m_fields[0]->GetTotPoints();
114 for(i = 1; i < 2*nvel; ++i)
123 std::string vSubStepIntScheme =
"RungeKutta2_ImprovedEuler";
125 if(
m_session->DefinesSolverInfo(
"SubStepIntScheme"))
127 vSubStepIntScheme =
m_session->GetSolverInfo(
"SubStepIntScheme");
137 int ntotpts =
m_fields[0]->GetTotPoints();
139 for(i = 1; i < 3*nvel; ++i)
147 ASSERTL0(0,
"Integration method not suitable: Options include BackwardEuler or BDFImplicitOrder1");
151 m_intSteps = IntegrationScheme->GetIntegrationSteps();
167 int nVariables = inarray.num_elements();
168 int nQuadraturePts = inarray[0].num_elements();
171 int ncoeffs =
m_fields[0]->GetNcoeffs();
176 for(i = 1; i < nVariables; ++i)
178 WeakAdv[i] = WeakAdv[i-1] + ncoeffs;
186 for(i = 1; i <
m_velocity.num_elements(); ++i)
188 Velfields[i] = Velfields[i-1] + nQuadraturePts;
195 for(i = 0; i < nVariables; ++i)
197 m_fields[i]->IProductWRTBase(outarray[i],WeakAdv[i]);
205 for(i = 0; i < nVariables; ++i)
211 m_fields[i]->MultiplyByElmtInvMass(WeakAdv[i], WeakAdv[i]);
214 m_fields[i]->BwdTrans(WeakAdv[i], outarray[i]);
226 ASSERTL1(inarray.num_elements() == outarray.num_elements(),
"Inarray and outarray of different sizes ");
228 for(
int i = 0; i < inarray.num_elements(); ++i)
230 Vmath::Vcopy(inarray[i].num_elements(),inarray[i],1,outarray[i],1);
243 int nConvectiveFields =
m_fields.num_elements()-1;
246 for(
int i = 0; i < nConvectiveFields; ++i)
287 for(n = 0; n < nvel; ++n)
291 for(i = nblocks-1; i > 0; --i)
300 for(i = 0; i < nvel; ++i)
310 for(n = 0; n < nvel; ++n)
312 for(i = 1; i < nblocks; ++i)
337 static int ncalls = 1;
357 cout <<
"Sub-integrating using "<< nsubsteps
362 for (
int m = 0; m < nint; ++m)
365 fields = integrationSoln->UpdateSolutionVector()[m];
374 for(n = 0; n < nsubsteps; ++n)
384 integrationSoln->SetSolVector(m,fields);
394 int n_element =
m_fields[0]->GetExpSize();
405 for(
int i = 0; i <
m_velocity.num_elements(); ++i)
411 for(
int el = 0; el < n_element; ++el)
414 (stdVelocity[el] * cLambda *
415 (ExpOrder[el]-1) * (ExpOrder[el]-1));
433 physfield.num_elements() == Outarray.num_elements(),
434 "Physfield and outarray are of different dimensions");
439 int nTracePts =
m_fields[0]->GetTrace()->GetNpoints();
457 for(i = 0; i < nDimensions; ++i)
459 m_fields[0]->ExtractTracePhys(velfield[i], Fwd);
463 for(i = 0; i < physfield.num_elements(); ++i)
467 m_fields[i]->GetFwdBwdTracePhys(physfield[i], Fwd, Bwd);
470 m_fields[0]->GetTrace()->Upwind(Vn, Fwd, Bwd, numflux);
473 Vmath::Vsub(nTracePts, numflux, 1, Fwd, 1, Fwd, 1);
474 Vmath::Vsub(nTracePts, numflux, 1, Bwd, 1, Bwd, 1);
481 m_fields[0]->AddFwdBwdTraceIntegral(Fwd,Bwd,Outarray[i]);
503 for(i = 0; i <= ord; ++i)
505 for(j = 0; j <= ord; ++j)
515 for(i = 0; i < nvel; ++i)
519 for(j = 1; j <= ord; ++j)
556 for(
int i = 0; i < acc_order; i++)
562 accelerationTerm, 1);
#define ASSERTL0(condition, msg)
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.
ExtrapolateFactory & GetExtrapolateFactory()
T Vmin(int n, const T *x, const int incx)
Return the minimum element in x - called vmin to avoid conflict with min.
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
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
boost::shared_ptr< Advection > AdvectionSharedPtr
A shared pointer to an Advection object.
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
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
void DefineProjection(FuncPointerT func, ObjectPointerT obj)
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
void DefineOdeRhs(FuncPointerT func, ObjectPointerT obj)
boost::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
void Neg(int n, T *x, const int incx)
Negate x = -x.
TimeIntegrationWrapperFactory & GetTimeIntegrationWrapperFactory()
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.
BDF multi-step scheme of order 2 (implicit)
boost::shared_ptr< TimeIntegrationSolution > TimeIntegrationSolutionSharedPtr
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
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.
Defines a callback function which evaluates the flux vector.
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, tDescription pDesc="")
Register a class with the factory.