Rythmos - Transient Integration for Differential Equations  Version of the Day
 All Classes Functions Variables Typedefs Pages
Public Types | Public Member Functions | Related Functions | List of all members
Rythmos::InterpolationBufferBase< Scalar > Class Template Referenceabstract

Base class for an interpolation buffer. More...

#include <Rythmos_InterpolationBufferBase.hpp>

Inheritance diagram for Rythmos::InterpolationBufferBase< Scalar >:
Inheritance graph
[legend]

Public Types

typedef Teuchos::ScalarTraits
< Scalar >::magnitudeType 
ScalarMag
 

Public Member Functions

virtual RCP< const
Thyra::VectorSpaceBase< Scalar > > 
get_x_space () const =0
 Return the space for x and x_dot. More...
 
virtual void addPoints (const Array< Scalar > &time_vec, const Array< RCP< const Thyra::VectorBase< Scalar > > > &x_vec, const Array< RCP< const Thyra::VectorBase< Scalar > > > &xdot_vec)=0
 Add points to the buffer. More...
 
virtual TimeRange< Scalar > getTimeRange () const =0
 Return the range of time values where interpolation calls can be performed. More...
 
virtual void getPoints (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) const =0
 Get values from the buffer at different time points. More...
 
virtual void getNodes (Array< Scalar > *time_vec) const =0
 Get interpolation nodes. More...
 
virtual void removeNodes (Array< Scalar > &time_vec)=0
 Remove nodes from the interpolation buffer. More...
 
virtual int getOrder () const =0
 Get order of interpolation. More...
 

Related Functions

(Note that these are not member functions.)

template<class Scalar >
RCP< const Thyra::VectorBase
< Scalar > > 
get_x (const InterpolationBufferBase< Scalar > &interpBuffer, const Scalar &t)
 Get a single point x(t) from an interpolation buffer. More...
 
template<class Scalar >
RCP< const Thyra::VectorBase
< Scalar > > 
get_xdot (const InterpolationBufferBase< Scalar > &interpBuffer, const Scalar &t)
 Get a single point xdot(t) from an interpolation buffer. More...
 
template<class Scalar >
void get_x_and_x_dot (const InterpolationBufferBase< Scalar > &interpBuffer, 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 x_dot at t. More...
 
template<class Scalar >
void assertTimePointsAreSorted (const Array< Scalar > &time_vec)
 Assert that a time point vector is sorted. More...
 
template<class Scalar >
void assertNoTimePointsBeforeCurrentTimeRange (const InterpolationBufferBase< Scalar > &interpBuffer, const Array< Scalar > &time_vec, const int &startingTimePointIndex=0)
 Assert that none of the time points fall before the current time range for an interpolation buffer object. More...
 
template<class Scalar >
void assertNoTimePointsInsideCurrentTimeRange (const InterpolationBufferBase< Scalar > &interpBuffer, const Array< Scalar > &time_vec)
 Assert that none of the time points fall inside the current time range for an interpolation buffer object. More...
 
template<class TimeType >
void selectPointsInTimeRange (const Array< TimeType > &points_in, const TimeRange< TimeType > &range, const Ptr< Array< TimeType > > &points_out)
 Select points from an Array that sit in a TimeRange. More...
 
template<class TimeType >
void removePointsInTimeRange (Array< TimeType > *points_in, const TimeRange< TimeType > &range)
 Remove points from an Array that sit in a TimeRange. More...
 
template<class Scalar >
bool getCurrentPoints (const InterpolationBufferBase< Scalar > &interpBuffer, const Array< Scalar > &time_vec, Array< RCP< const Thyra::VectorBase< Scalar > > > *x_vec, Array< RCP< const Thyra::VectorBase< Scalar > > > *xdot_vec, int *nextTimePointIndex)
 Get time points in the current range of an interpolation buffer object. More...
 

Detailed Description

template<class Scalar>
class Rythmos::InterpolationBufferBase< Scalar >

Base class for an interpolation buffer.

An interpolation buffer represents the data and the ability to represent and interpolate the values of some transient solution x and x_dot over a contiguous range of time.

Definitions

ToDo: Finish documentation!

Definition at line 68 of file Rythmos_InterpolationBufferBase.hpp.

Member Typedef Documentation

template<class Scalar>
typedef Teuchos::ScalarTraits<Scalar>::magnitudeType Rythmos::InterpolationBufferBase< Scalar >::ScalarMag

Definition at line 76 of file Rythmos_InterpolationBufferBase.hpp.

Member Function Documentation

template<class Scalar>
virtual RCP<const Thyra::VectorSpaceBase<Scalar> > Rythmos::InterpolationBufferBase< Scalar >::get_x_space ( ) const
pure virtual

Return the space for x and x_dot.

This space can be used to create vectors for calling addPoints() for instance and is also useful for writing unit testing software.

Also note that this space may not be the same as the space returned from StepperBase::getModel()->get_x_space() in some concrete StepperBase subclasses.

Implemented in Rythmos::ForwardSensitivityStepper< Scalar >, Rythmos::BackwardEulerStepper< Scalar >, Rythmos::ForwardEulerStepper< Scalar >, Rythmos::DefaultIntegrator< Scalar >, Rythmos::ThetaStepper< Scalar >, Rythmos::ImplicitBDFStepper< Scalar >, Rythmos::ExplicitTaylorPolynomialStepper< Scalar >, Rythmos::ImplicitRKStepper< Scalar >, Rythmos::ExplicitRKStepper< Scalar >, and Rythmos::InterpolationBuffer< Scalar >.

template<class Scalar>
virtual void Rythmos::InterpolationBufferBase< Scalar >::addPoints ( const Array< Scalar > &  time_vec,
const Array< RCP< const Thyra::VectorBase< Scalar > > > &  x_vec,
const Array< RCP< const Thyra::VectorBase< Scalar > > > &  xdot_vec 
)
pure virtual

Add points to the buffer.

Parameters
time_vec[in] Array (length n) of time points.
x_vec[in] Array (length n) of state vectors at each time point in time_vec. Specifically, x_vec[i] is the state vector at time time_vec[i], for i=0...n-1. The RCP for the vectors may or may not be copied by *this. An implementation is not required to copy the RCP's to the vector objects but instead might just use the vectors or do a deep copy.
xdot_vec[in] Array (length n) of state time differential vector at each time point in time_vec. Specifically, xdot_vec[i] is the state differential vector at time time_vec[i], for i=0...n-1. The RCP for the vectors may or may not be copied by *this. An implementation is not required to copy the RCP's to the vector objects. The implementation may copy the actually vectors or might just perform a shallow copy by copying the RCP objects.

Preconditions:

  • time_vec.size()!=0
  • time_vec must have unique and sorted values in ascending order
  • time_vec.size()==x_vec.size()
  • time_vec.size()==xdot_vec.size()
  • x_vec[i] != null for i = 0...n-1
  • xdot_vec[i] != null for i = 0...n-1
  • getTimeRange().isInRange(time_vec[i]) == false, for i=0...n-1

Postconditions:

Implemented in Rythmos::BackwardEulerStepper< Scalar >, Rythmos::ForwardEulerStepper< Scalar >, Rythmos::ExplicitTaylorPolynomialStepper< Scalar >, Rythmos::DefaultIntegrator< Scalar >, Rythmos::ThetaStepper< Scalar >, Rythmos::ImplicitBDFStepper< Scalar >, Rythmos::ImplicitRKStepper< Scalar >, and Rythmos::InterpolationBuffer< Scalar >.

template<class Scalar>
virtual TimeRange<Scalar> Rythmos::InterpolationBufferBase< Scalar >::getTimeRange ( ) const
pure virtual

Return the range of time values where interpolation calls can be performed.

A return value of returnVal.isValid()==false means that there is no time range for which interpolation can be performed. Otherwise, returnVal gives the time range for which state information can be returned.

Implemented in Rythmos::ForwardSensitivityStepper< Scalar >, Rythmos::BackwardEulerStepper< Scalar >, Rythmos::ForwardEulerStepper< Scalar >, Rythmos::ExplicitTaylorPolynomialStepper< Scalar >, Rythmos::DefaultIntegrator< Scalar >, Rythmos::ThetaStepper< Scalar >, Rythmos::ImplicitBDFStepper< Scalar >, Rythmos::ImplicitRKStepper< Scalar >, Rythmos::ExplicitRKStepper< Scalar >, and Rythmos::InterpolationBuffer< Scalar >.

template<class Scalar>
virtual void Rythmos::InterpolationBufferBase< Scalar >::getPoints ( 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 
) const
pure virtual

Get values from the buffer at different time points.

Parameters
time_vec[in] Array (length n) of time points to get.
x_vec[out] On output, if x_vec != 0, *x_vec will be resized to n = time_vec.size() and (*x_vec)[i] will be the state vector at time time_vec[i], for i=0...n-1. This argument can be left NULL in which case it will not be filled.
xdot_vec[out] On output, if xdot_vec != 0, *xdot_vec will be resized to n = time_vec.size() and (*xdot_vec)[i] will be the state derivative vector at time time_vec[i], for i=0...n-1. This argument can be left NULL in which case it will not be filled.
accuracy_vec[out] This contains an estimate of the accuracy of the interpolation. If you asked for a node, this should be zero.

Preconditions:

  • range.isInRange(time_vec[i]), for i=0...n-1, where range = this->getTimeRange().
  • time_vec must have unique and sorted values in ascending order

Implemented in Rythmos::ForwardSensitivityStepper< Scalar >, Rythmos::BackwardEulerStepper< Scalar >, Rythmos::ForwardEulerStepper< Scalar >, Rythmos::ExplicitTaylorPolynomialStepper< Scalar >, Rythmos::ThetaStepper< Scalar >, Rythmos::ImplicitBDFStepper< Scalar >, and Rythmos::InterpolationBuffer< Scalar >.

template<class Scalar>
virtual void Rythmos::InterpolationBufferBase< Scalar >::getNodes ( Array< Scalar > *  time_vec) const
pure virtual

Get interpolation nodes.

This function will return the time points where actual data is stored. This information can be used to get the actual nodal values themselves using the getPoints() function.

This function may return no nodes at all, and yet still have a valid timeRange.

Implemented in Rythmos::ForwardSensitivityStepper< Scalar >, Rythmos::BackwardEulerStepper< Scalar >, Rythmos::ForwardEulerStepper< Scalar >, Rythmos::ExplicitTaylorPolynomialStepper< Scalar >, Rythmos::DefaultIntegrator< Scalar >, Rythmos::ThetaStepper< Scalar >, Rythmos::ImplicitBDFStepper< Scalar >, Rythmos::ImplicitRKStepper< Scalar >, Rythmos::ExplicitRKStepper< Scalar >, and Rythmos::InterpolationBuffer< Scalar >.

template<class Scalar>
virtual void Rythmos::InterpolationBufferBase< Scalar >::removeNodes ( Array< Scalar > &  time_vec)
pure virtual

Remove nodes from the interpolation buffer.

Preconditions:

  • time_vec must have unique and sorted values in ascending order
  • range.isInRange(time_vec[i]), for i=0..n-1, where range = this->getTimeRange().
  • Each time_vec[i] must equal a node in the interpolation buffer.

Postconditions:

  • Each time_vec[i] is no longer a node in the interpolation buffer.

Implemented in Rythmos::ForwardSensitivityStepper< Scalar >, Rythmos::BackwardEulerStepper< Scalar >, Rythmos::ForwardEulerStepper< Scalar >, Rythmos::ExplicitTaylorPolynomialStepper< Scalar >, Rythmos::DefaultIntegrator< Scalar >, Rythmos::ThetaStepper< Scalar >, Rythmos::ImplicitBDFStepper< Scalar >, Rythmos::ImplicitRKStepper< Scalar >, Rythmos::ExplicitRKStepper< Scalar >, and Rythmos::InterpolationBuffer< Scalar >.

template<class Scalar>
virtual int Rythmos::InterpolationBufferBase< Scalar >::getOrder ( ) const
pure virtual

Friends And Related Function Documentation

template<class Scalar >
RCP< const Thyra::VectorBase< Scalar > > get_x ( const InterpolationBufferBase< Scalar > &  interpBuffer,
const Scalar &  t 
)
related

Get a single point x(t) from an interpolation buffer.

Definition at line 253 of file Rythmos_InterpolationBufferBase.hpp.

template<class Scalar >
RCP< const Thyra::VectorBase< Scalar > > get_xdot ( const InterpolationBufferBase< Scalar > &  interpBuffer,
const Scalar &  t 
)
related

Get a single point xdot(t) from an interpolation buffer.

Definition at line 277 of file Rythmos_InterpolationBufferBase.hpp.

template<class Scalar >
void get_x_and_x_dot ( const InterpolationBufferBase< Scalar > &  interpBuffer,
const Scalar  t,
const Ptr< RCP< const Thyra::VectorBase< Scalar > > > &  x,
const Ptr< RCP< const Thyra::VectorBase< Scalar > > > &  x_dot 
)
related

Nonmember helper function to get x and x_dot at t.

Definition at line 299 of file Rythmos_InterpolationBufferBase.hpp.

template<class Scalar >
void assertTimePointsAreSorted ( const Array< Scalar > &  time_vec)
related

Assert that a time point vector is sorted.

template<class Scalar >
void assertNoTimePointsBeforeCurrentTimeRange ( const InterpolationBufferBase< Scalar > &  interpBuffer,
const Array< Scalar > &  time_vec,
const int &  startingTimePointIndex = 0 
)
related

Assert that none of the time points fall before the current time range for an interpolation buffer object.

Parameters
interpBuffer[in] The interpolation buffer defining the time range. The reason that we accept the interpolation buffer and not just the time range is so we can create a better error message.
time_vec[in] The array of time points
startingTimePointIndex[in] The time index in time_vec to begin asserting time points. The default is 0.

This function with throw an exception if any of the time points are found to be before the current time range.

template<class Scalar >
void assertNoTimePointsInsideCurrentTimeRange ( const InterpolationBufferBase< Scalar > &  interpBuffer,
const Array< Scalar > &  time_vec 
)
related

Assert that none of the time points fall inside the current time range for an interpolation buffer object.

Parameters
interpBuffer[in] The interpolation buffer defining the time range. The reason that we accept the interpolation buffer and not just the time range is so we can create a better error message.
time_vec[in] The array of time points

This function with throw an exception if any of the time points are found to be inside the current time range.

template<class TimeType >
void selectPointsInTimeRange ( const Array< TimeType > &  points_in,
const TimeRange< TimeType > &  range,
const Ptr< Array< TimeType > > &  points_out 
)
related

Select points from an Array that sit in a TimeRange.

template<class TimeType >
void removePointsInTimeRange ( Array< TimeType > *  points_in,
const TimeRange< TimeType > &  range 
)
related

Remove points from an Array that sit in a TimeRange.

template<class Scalar >
bool getCurrentPoints ( const InterpolationBufferBase< Scalar > &  interpBuffer,
const Array< Scalar > &  time_vec,
Array< RCP< const Thyra::VectorBase< Scalar > > > *  x_vec,
Array< RCP< const Thyra::VectorBase< Scalar > > > *  xdot_vec,
int *  nextTimePointIndex 
)
related

Get time points in the current range of an interpolation buffer object.

Parameters
interpBuffer[in] The interpolation buffer object that will be used to get state values at different time points
time_vec[in] The whole time vector, including time points that fall outside of the current time range in interpBuffer.
x_vec[in/out] This argument is optional and it is allowed for x_vec==0. However, if x_vec!=0, then on input x_vec->size()==time_vec.size() must be true. On output, *x_vec will be filled with the state vectors for the current time points given in time_vec.
xdot_vec[in/out] This argument is optional and it is allowed for xdot_vec==0. However, if xdot_vec!=0, then on input xdot_vec->size()==time_vec.size() must be true. On output, *xdot_vec will be filled with the state derivative vectors for the current time points given in time_vec.
nextTimePointIndex[in/out] On input, *nextTimePointIndex gives first time point time_vec[*nextTimePointIndex] to extract state values for. On output, *nextTimePointIndex will be incremented such that there are no more time points left (i.e. nextTimePointIndex==time_vec.size() or such that time_vec[*nextTimePointIndex] > interpBuffer.getTimeRange().upper()

Preconditions:

nextTimePointIndex!=0 0 <= *nextTimePointIndex < time_vec.size() time_vec[*nextTimePointIndex] >= interpBuffer.getTimeRange().lower()

Preconditions:

  • [returnVal==true] *nextTimePointIndex is greater on output than on input
  • [returnVal==false] *nextTimePointIndex is unchanged on output

Returns
Returns true if one or more time points where extracted from interpBuffer and returns false otherwise.


The documentation for this class was generated from the following files: