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;
185 for(i = 1; i <
m_velocity.num_elements(); ++i)
187 Velfields[i] = Velfields[i-1] + nQuadraturePts;
194 for(i = 0; i < nVariables; ++i)
196 m_fields[i]->IProductWRTBase(outarray[i],WeakAdv[i]);
204 for(i = 0; i < nVariables; ++i)
210 m_fields[i]->MultiplyByElmtInvMass(WeakAdv[i], WeakAdv[i]);
213 m_fields[i]->BwdTrans(WeakAdv[i], outarray[i]);
225 ASSERTL1(inarray.num_elements() == outarray.num_elements(),
"Inarray and outarray of different sizes ");
227 for(
int i = 0; i < inarray.num_elements(); ++i)
229 Vmath::Vcopy(inarray[i].num_elements(),inarray[i],1,outarray[i],1);
279 for(n = 0; n < nvel; ++n)
283 for(i = nblocks-1; i > 0; --i)
292 for(i = 0; i < nvel; ++i)
302 for(n = 0; n < nvel; ++n)
304 for(i = 1; i < nblocks; ++i)
329 static int ncalls = 1;
347 cout <<
"Sub-integrating using "<< nsubsteps
352 for (
int m = 0; m < nint; ++m)
355 fields = integrationSoln->UpdateSolutionVector()[m];
364 for(n = 0; n < nsubsteps; ++n)
374 integrationSoln->SetSolVector(m,fields);
384 int n_element =
m_fields[0]->GetExpSize();
394 for(
int i = 0; i <
m_velocity.num_elements(); ++i)
400 for(
int el = 0; el < n_element; ++el)
403 (stdVelocity[el] * cLambda *
404 (ExpOrder[el]-1) * (ExpOrder[el]-1));
422 physfield.num_elements() == Outarray.num_elements(),
423 "Physfield and outarray are of different dimensions");
428 int nTracePts =
m_fields[0]->GetTrace()->GetNpoints();
446 for(i = 0; i < nDimensions; ++i)
448 m_fields[0]->ExtractTracePhys(velfield[i], Fwd);
452 for(i = 0; i < physfield.num_elements(); ++i)
456 m_fields[i]->GetFwdBwdTracePhys(physfield[i], Fwd, Bwd);
459 m_fields[0]->GetTrace()->Upwind(Vn, Fwd, Bwd, numflux);
462 Vmath::Vsub(nTracePts, numflux, 1, Fwd, 1, Fwd, 1);
463 Vmath::Vsub(nTracePts, numflux, 1, Bwd, 1, Bwd, 1);
470 m_fields[0]->AddFwdBwdTraceIntegral(Fwd,Bwd,Outarray[i]);
492 for(i = 0; i <= ord; ++i)
494 for(j = 0; j <= ord; ++j)
504 for(i = 0; i < nvel; ++i)
508 for(j = 1; j <= ord; ++j)
#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
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.