Extrapolate a polynomial passing through up to three points. More...

Public Types
using	value_type = get_value_type<ContainerType>

using	container_type = ContainerType

Public Member Functions
	Extrapolation ()
	Leave values uninitialized.

	Extrapolation (unsigned max, const ContainerType &copyable)
	Set maximum extrapolation order and allocate memory.

void	set_max (unsigned max, const ContainerType &copyable)
	Set maximum extrapolation order and allocate memory.

unsigned	get_max () const
	Current extrapolation count.

bool	exists (value_type t) const
	Check if time exists in current points.

template<class ContainerType0 >
void	extrapolate (value_type t, ContainerType0 &new_x) const
	Extrapolate value to given time.

template<class ContainerType0 >
void	derive (value_type t, ContainerType0 &dot_x) const
	Evaluate first derivative of interpolating polynomial

template<class ContainerType0 >
void	extrapolate (ContainerType0 &new_x) const
	Extrapolate value (equidistant version)

template<class MatrixType0 , class ContainerType0 , class ContainerType1 >
void	extrapolate_least_squares (MatrixType0 &&A, const ContainerType0 &b, ContainerType0 &new_x, const ContainerType1 &weights)
	EXPERIMENTAL Perform a least squares extrapolation.

template<class ContainerType0 >
void	derive (ContainerType0 &dot_x) const
	Evaluate first derivative of interpolating polynomial (equidistant version)

template<class ContainerType0 >
void	update (value_type t_new, const ContainerType0 &new_entry)
	insert a new entry, deleting the oldest entry or update existing entry

template<class ContainerType0 >
void	update (const ContainerType0 &new_entry)
	insert a new entry

const ContainerType &	head () const
	return the current head (the one most recently inserted)

ContainerType &	tail ()
	DEPRECATED write access to tail value ( the one that will be deleted in the next update, undefined if max==0)

const ContainerType &	tail () const
	DEPRECATED read access to tail value ( the one that will be deleted in the next update, undefined if max==0)

Detailed Description

template<class ContainerType>
struct dg::Extrapolation< ContainerType >

Extrapolate a polynomial passing through up to three points.

This class constructs an interpolating polynomial through up to three given points and evaluates its value or its derivative at a new point. The points can be updated to get a new polynomial.

The intention of this class is to provide an initial guess for iterative solvers based on past solutions:

\[ x_{init} = \alpha_0 x_0 + \alpha_{-1}x_{-1} + \alpha_{-2} x_{-2}\]

where the indices indicate the current (0) and past (negative) solutions. Choose between 1 (constant), 2 (linear) or 3 (parabola) extrapolation. The user can choose to provide a time value t_i associated with the x_i, which are then used to compute the coefficients alpha_i (using Lagrange interpolation). Otherwise an equidistant distribution is assumed.

        //Least-squares fit through many points
        dg::Grid1d grid( 0., 2.*M_PI, 1, 200);
        // We have a travelling wave at various times
        dg::HVec y0 = dg::evaluate( [=](double x) { return cos( x -0.0);}, grid);
        dg::HVec y1 = dg::evaluate( [=](double x) { return cos( x -0.1);}, grid);
        dg::HVec y2 = dg::evaluate( [=](double x) { return cos( x -0.2);}, grid);
        dg::HVec y  = dg::evaluate( [=](double x) { return cos( x -0.3);}, grid);
        dg::HVec guess ( y);
        double norm = dg::blas1::dot( y,y);
 
        // Given (t0, y0), (t1, y1), (t2, y2) give a good guess for y at t3
 
        dg::Extrapolation<dg::HVec> extra(3,y);
 
        extra.update( 0.0, y0);
        extra.extrapolate( 0.3, guess);
        // From 1 point we expect simply the same point
        dg::blas1::axpby( 1., y0, -1., guess);
        double err = sqrt( dg::blas1::dot( guess, guess) /norm);
        INFO( "Difference 1 is "<<err);
        CHECK( err < 1e-14);
 
        // From 3 points we have a parabolic extrapolation
        extra.update( 0.1, y1);
        extra.update( 0.1, y1); // This one silently overwrites previous t1=0.1
        extra.update( 0.2, y2);
        extra.extrapolate( 0.3, guess);
        dg::blas1::axpby( 1., y, -1., guess);
        err = sqrt( dg::blas1::dot( guess, guess) /norm);
        INFO( "Difference 2 is "<<err);
        // The accuracy of polynomial extrapolation is not great ...
        CHECK( err < 1e-3);

Note: Since extrapolation with higher order polynomials is so prone to oscillations anything higher than linear rarely leads to anything useful. So best stick to constant or linear extrapolation; The derivative of the interpolating polynomial at a new point reduces to familiar finite difference formulas

Template Parameters

ContainerType

Any class for which a specialization of TensorTraits exists and which fulfills the requirements of the there defined data and execution policies derived from AnyVectorTag and AnyPolicyTag. Among others

dg::HVec (serial), dg::DVec (cuda / omp), dg::MHVec (mpi + serial) or dg::MDVec (mpi + cuda / omp)
std::vector<dg::DVec> (vector of shared device vectors), std::array<double, 4> (array of 4 doubles) or std::map < std::string, dg::DVec> ( a map of named vectors)
double (scalar) and other primitive types ...

If there are several ContainerTypes in the argument list, then TensorTraits must exist for all of them

See also: See The dg dispatch system for a detailed explanation of our type dispatch system; https://en.wikipedia.org/wiki/Extrapolation

Member Typedef Documentation

◆ container_type

template<class ContainerType >

using dg::Extrapolation< ContainerType >::container_type = ContainerType

◆ value_type

template<class ContainerType >

using dg::Extrapolation< ContainerType >::value_type = get_value_type<ContainerType>

Constructor & Destructor Documentation

◆ Extrapolation() [1/2]

template<class ContainerType >

dg::Extrapolation< ContainerType >::Extrapolation ( )

inline

Leave values uninitialized.

◆ Extrapolation() [2/2]

template<class ContainerType >

dg::Extrapolation< ContainerType >::Extrapolation	(	unsigned	max,
		const ContainerType &	copyable )

inline

Set maximum extrapolation order and allocate memory.

Parameters

max	maximum of vectors to use for extrapolation. Choose between 0 (no extrapolation) 1 (constant), 2 (linear) or 3 (parabola) extrapolation. Higher values currently default back to a linear extrapolation
copyable	the memory is allocated based on this vector

Member Function Documentation

◆ derive() [1/2]

template<class ContainerType >

template<class ContainerType0 >

void dg::Extrapolation< ContainerType >::derive ( ContainerType0 & dot_x ) const

inline

Evaluate first derivative of interpolating polynomial (equidistant version)

Parameters

dot_x (write only) contains derived value on output

Note: Assumes that time equals t0 such that a backward difference formula will be evaluated

Template Parameters

ContainerType0 must be usable with ContainerType in The dg dispatch system

◆ derive() [2/2]

template<class ContainerType >

template<class ContainerType0 >

void dg::Extrapolation< ContainerType >::derive	(	value_type	t,
		ContainerType0 &	dot_x ) const

inline

Evaluate first derivative of interpolating polynomial

Equivalent to constructing the interpolating polynomial, deriving it once and then evaluating it at the required point

Parameters

t	time at which derivative of interpolating polynomial is evaluated
dot_x	(write only) contains derived value on output

Note: If t is chosen as the latest time of update t0, then the result coincides with the backward difference formula of order max

Attention: If max==1, the result is 0 (derivative of a constant); If the update function has not been called enough times to fill all values the result depends: (i) never called => dot_x is zero (ii) called at least once => the interpolating polynomial is constructed with all available values

Template Parameters

ContainerType0 must be usable with ContainerType in The dg dispatch system

◆ exists()

template<class ContainerType >

bool dg::Extrapolation< ContainerType >::exists ( value_type t ) const

inline

Check if time exists in current points.

Note: The check is numerical and returns true if \(|t - t_i| < 10^{-14}\) for all i; If extrapolate() is called on an exisiting time, the corresponding vector is returned while update() will overwrite the existing vector

Parameters

t	The time to check for

Returns: true if time exists in current points, false else

◆ extrapolate() [1/2]

template<class ContainerType >

template<class ContainerType0 >

void dg::Extrapolation< ContainerType >::extrapolate ( ContainerType0 & new_x ) const

inline

Extrapolate value (equidistant version)

Parameters

new_x (write only) contains extrapolated value on output

Note: Assumes that extrapolation time equals last inserted time+1

Template Parameters

ContainerType0 must be usable with ContainerType in The dg dispatch system

◆ extrapolate() [2/2]

template<class ContainerType >

template<class ContainerType0 >

void dg::Extrapolation< ContainerType >::extrapolate	(	value_type	t,
		ContainerType0 &	new_x ) const

inline

Extrapolate value to given time.

Construt and evaluate the interpolating polynomial at a given point

Parameters

t	time to which to extrapolate (or at which interpolating polynomial is evaluated)
new_x	(write only) contains extrapolated value on output

Template Parameters

ContainerType0 must be usable with ContainerType in The dg dispatch system

Attention: If the update function has not been called enough times to fill all values the result depends: (i) never called => new_x is zero (ii) called at least once => the interpolating polynomial is constructed with all available values

◆ extrapolate_least_squares()

template<class ContainerType >

template<class MatrixType0 , class ContainerType0 , class ContainerType1 >

void dg::Extrapolation< ContainerType >::extrapolate_least_squares	(	MatrixType0 &&	A,
		const ContainerType0 &	b,
		ContainerType0 &	new_x,
		const ContainerType1 &	weights )

inline

EXPERIMENTAL Perform a least squares extrapolation.

This algorithm computes \( b_i = A x_i\), then solves the least squares problem \( \min_a || B a - b||_W \) with \( b_i\) the columns of \( B\) to compute \( x_{new} = \sum_i a_i x_i\).

Note: This is different from LeastSquaresExtrapolation if the matrix \( A\) is time-dependent.; So far this was only tested for toefl simulations with default parameters where we get mixed results. Depending on the max parameter one can see a slow down or a speed-up. Around 10 we observed a positive effect of about 10% acceleration

Parameters

A	the matrix
b	the right hand side
new_x	Contains initial guess based on a least squares extrapolation on output
weights	The weights/volume in which to minimize

◆ get_max()

template<class ContainerType >

unsigned dg::Extrapolation< ContainerType >::get_max ( ) const

inline

Current extrapolation count.

Note: This may not coincide with the max set in the constructor if values have not been updated yet

Returns: The current size of the extrapolation

◆ head()

template<class ContainerType >

const ContainerType & dg::Extrapolation< ContainerType >::head ( ) const

inline

return the current head (the one most recently inserted)

Returns: current head (undefined if max==0)

◆ set_max()

template<class ContainerType >

void dg::Extrapolation< ContainerType >::set_max	(	unsigned	max,
		const ContainerType &	copyable )

inline

Set maximum extrapolation order and allocate memory.

Parameters

max	maximum of vectors to use for extrapolation. Choose between 0 (no extrapolation) 1 (constant), 2 (linear) or 3 (parabola) extrapolation. Higher values currently default back to a linear extrapolation
copyable	the memory is allocated based on this vector

◆ tail() [1/2]

template<class ContainerType >

ContainerType & dg::Extrapolation< ContainerType >::tail ( )

inline

DEPRECATED write access to tail value ( the one that will be deleted in the next update, undefined if max==0)

◆ tail() [2/2]

template<class ContainerType >

const ContainerType & dg::Extrapolation< ContainerType >::tail ( ) const

inline

DEPRECATED read access to tail value ( the one that will be deleted in the next update, undefined if max==0)

◆ update() [1/2]

template<class ContainerType >

template<class ContainerType0 >

void dg::Extrapolation< ContainerType >::update ( const ContainerType0 & new_entry )

inline

insert a new entry

Parameters

new_entry the new entry

Note: Assumes new time equals last inserted time+1

Template Parameters

ContainerType0 must be usable with ContainerType in The dg dispatch system

◆ update() [2/2]

template<class ContainerType >

template<class ContainerType0 >

void dg::Extrapolation< ContainerType >::update	(	value_type	t_new,
		const ContainerType0 &	new_entry )

inline

insert a new entry, deleting the oldest entry or update existing entry

Parameters

t_new	the time for the new entry
new_entry	the new entry ( replaces value of existing entry if `t_new` already exists)

Template Parameters

ContainerType0 must be usable with ContainerType in The dg dispatch system

The documentation for this struct was generated from the following file:

extrapolation.h

Public Types

Public Member Functions

Detailed Description

Member Typedef Documentation

◆ container_type

◆ value_type

Constructor & Destructor Documentation

◆ Extrapolation() [1/2]

◆ Extrapolation() [2/2]

Member Function Documentation

◆ derive() [1/2]

◆ derive() [2/2]

◆ exists()

◆ extrapolate() [1/2]

◆ extrapolate() [2/2]

◆ extrapolate_least_squares()

◆ get_max()

◆ head()

◆ set_max()

◆ tail() [1/2]

◆ tail() [2/2]

◆ update() [1/2]

◆ update() [2/2]