KASKADE 7 development version
|
A convenience class simplifying the implementation of different SDCTimeGrid derived classes. More...
#include <sdc.hh>
A convenience class simplifying the implementation of different SDCTimeGrid derived classes.
It stores the collocation points as well as integration and differentiation matrices, and provides access to them by implementing the pure virtual methods of the base class.
Public Types | |
typedef Dune::DynamicVector< double > | RealVector |
The type used for real vectors. More... | |
typedef DynamicMatrix< double > | RealMatrix |
The type used for real (dense) matrices. More... | |
Public Member Functions | |
virtual RealVector const & | points () const |
Time points in the time step. More... | |
virtual RealMatrix const & | integrationMatrix () const |
Integration matrix \( S \). More... | |
virtual RealMatrix const & | differentiationMatrix () const |
Differentiation matrix \( D \). More... | |
double | point (int k) const |
Time points in the time step. More... | |
virtual void | refine (RealMatrix &p)=0 |
Perform refinement of the grid, filling the prolongation matrix. More... | |
virtual RealMatrix | interpolate (RealVector const &x) const =0 |
Compute interpolation coefficients. More... | |
Protected Attributes | |
RealVector | pts |
RealMatrix | integ |
RealMatrix | diff |
|
inherited |
|
inherited |
|
inlinevirtual |
Differentiation matrix \( D \).
This computes a matrix such that polynomials given by values at \( t_0, \dots,t_n\) can easily be differentiated.
The Lagrangian interpolation functions \( L_k \) are defined by \( L_k(t_i) = \delta_{ik} , \quad i=0,\dots,n\). The matrix \( D \in \mathbb{R}^{n+1\times n+1}\) contains the values
\[ D_{ik} = \dot L_k(t_i) \]
This way, if \( u \) is defined in terms of its function values \( v_i = u(t_i) \), its derivatives can be evaluated by a matrix-vector multiplication:
\[ \dot u(\tau_i) = (Dv)_i \]
Implements Kaskade::SDCTimeGrid.
|
inlinevirtual |
Integration matrix \( S \).
This computes a matrix such that functions given by values at \( t_0,\dots,t_n\) can be easily integrated.
Interpolation is based on the nodes \(t_0\dots,t_n\). Depending on the actual implementation, the initial point \( t_0 \) might be ignored (i.e. have a quadrature weight of 0).
The Lagrangian interpolation functions \( L_k \) are defined by \( L_k(t_i) = \delta_{ik},\quad i=1,\dots,n \). The matrix \( S \in \mathbb{R}^{n\times n+1}\) contains the values
\[ S_{ik} = \int_{\tau=t_i}^{t_{i+1}} L_k(\tau) \, d\tau \]
(the leading column being zero). This way, if \( u \) is defined in terms of its function values \( v_i = u(t_i),\quad i=1,\dots,n \), the integrals can be evaluated by a matrix-vector multiplication:
\[ \int_{\tau=t_i}^{t_{i+1}} u(\tau) \, d\tau = (Sv)_i \]
Implements Kaskade::SDCTimeGrid.
|
pure virtualinherited |
Compute interpolation coefficients.
Returns a matrix \( w \in \mathbb{R}^{m+1\times n+1} \), such that the interpolation polynomial \( p \) to the values \( y_i \) at grid points \( t_i \) can be evaluated as
\[ p(x_i) = \sum_{j=0}^n w_{ij} y_j, \quad i=0,\dots,m. \]
Implemented in Kaskade::LobattoTimeGrid, Kaskade::RadauTimeGrid, and Kaskade::GaussTimeGrid.
|
inlineinherited |
|
inlinevirtual |
Time points in the time step.
The time step \( [t_0, t_n] \) contains \( n+1 \) time points \( t_i \), including the end points. Those are provided here. The time points are stored in increasing order.
Implements Kaskade::SDCTimeGrid.
|
pure virtualinherited |
Perform refinement of the grid, filling the prolongation matrix.
If the function representation is not sufficiently accurate, a finer grid of time points can be tried. This method refines the grid to \( m+1 \) time points \( s_i \), \( m>n \), and fills a prolongation matrix \( P\in \mathbb{R}^{m+1\times n+1} \) such that with \( v_i = u(t_i) \) and \( w=Pv \) it holds that \( w_i = u(s_i) \).
Different derived classes may implement this in different ways, or provide their own, more flexible ways of refining the grid.
[out] | p | the prolongation matrix |
Implemented in Kaskade::LobattoTimeGrid, Kaskade::RadauTimeGrid, and Kaskade::GaussTimeGrid.
|
protected |
Definition at line 192 of file sdc.hh.
Referenced by differentiationMatrix().
|
protected |
Definition at line 191 of file sdc.hh.
Referenced by integrationMatrix().
|
protected |