29 #ifndef Rythmos_INTEGRATOR_BASE_H 30 #define Rythmos_INTEGRATOR_BASE_H 32 #include "Rythmos_InterpolationBufferBase.hpp" 33 #include "Rythmos_StepperBase.hpp" 34 #include "Teuchos_as.hpp" 40 namespace Exceptions {
61 template<
class Scalar>
67 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType
ScalarMag;
71 {
return Teuchos::null; }
103 virtual void setStepper(
105 const Scalar &finalTime,
106 const bool landOnFinalTime =
true 114 virtual Teuchos::RCP<const StepperBase<Scalar> > getStepper()
const =0;
121 virtual Teuchos::RCP<StepperBase<Scalar> > getNonconstStepper()
const =0;
131 virtual RCP<StepperBase<Scalar> > unSetStepper() =0;
175 virtual void getFwdPoints(
176 const Array<Scalar>& time_vec,
177 Array<RCP<
const Thyra::VectorBase<Scalar> > >* x_vec,
178 Array<RCP<
const Thyra::VectorBase<Scalar> > >* xdot_vec,
179 Array<ScalarMag>* accuracy_vec
203 template<
class Scalar>
204 RCP<const Thyra::VectorBase<Scalar> >
207 Array<Scalar> time_vec;
208 time_vec.push_back(t);
209 Array<RCP<const Thyra::VectorBase<Scalar> > > x_vec;
239 template<
class Scalar>
243 const Ptr<RCP<
const Thyra::VectorBase<Scalar> > > &x,
244 const Ptr<RCP<
const Thyra::VectorBase<Scalar> > > &x_dot
247 Array<Scalar> time_vec;
248 time_vec.push_back(t);
249 Array<RCP<const Thyra::VectorBase<Scalar> > > x_vec;
250 Array<RCP<const Thyra::VectorBase<Scalar> > > x_dot_vec;
253 nonnull(x) ? &x_vec : 0,
254 nonnull(x_dot) ? &x_dot_vec : 0,
260 *x_dot = x_dot_vec[0];
265 template<
class Scalar>
266 void get_fwd_x_and_x_dot(
269 RCP<
const Thyra::VectorBase<Scalar> > *x,
270 RCP<
const Thyra::VectorBase<Scalar> > *x_dot
273 get_fwd_x_and_x_dot(integrator, t, ptr(x), ptr(x_dot));
280 #endif // Rythmos_INTEGRATOR_BASE_H
RCP< const Thyra::VectorBase< Scalar > > get_fwd_x(IntegratorBase< Scalar > &integrator, const Scalar t)
Nonmember helper function to get x at a (forward) time t.
Base class for defining stepper functionality.
Abstract interface for time integrators.
Teuchos::ScalarTraits< Scalar >::magnitudeType ScalarMag
virtual RCP< IntegratorBase< Scalar > > cloneIntegrator() const
Base class for an interpolation buffer.
void get_fwd_x_and_x_dot(IntegratorBase< Scalar > &integrator, const Scalar t, const Ptr< RCP< const Thyra::VectorBase< Scalar > > > &x, const Ptr< RCP< const Thyra::VectorBase< Scalar > > > &x_dot)
Nonmember helper function to get x and/or x_dot at s (forward) time t.
virtual void getFwdPoints(const Array< Scalar > &time_vec, Array< RCP< const Thyra::VectorBase< Scalar > > > *x_vec, Array< RCP< const Thyra::VectorBase< Scalar > > > *xdot_vec, Array< ScalarMag > *accuracy_vec)=0
Get values at time points both inside and outside (forward) of current TimeRange. ...