Medusa  1.1
Coordinate Free Mehless Method implementation
mm::RBFInterpolant< RBFType, vec_t > Class Template Reference

#include <RBFInterpolant_fwd.hpp>

Detailed Description

template<typename RBFType, typename vec_t>
class mm::RBFInterpolant< RBFType, vec_t >

Class representing a RBF Interpolant over a set of nodes of the form.

\[ u(p) = \sum_{i=1}^n \alpha_i \phi\left(\left\|\frac{p - p_i}{s}\right\|\right) + \sum_{i=1}^s \beta_i q_i\left(\frac{p - p_c}{s} \right), \]

where \(\phi\) is the RBF and \(q_i\) are monomials.

Note
This class does not compute the interpolant, that is done by RBFFD::getApproximant.
See also
RBFFD
Template Parameters
basis_tRadial basis function used. Must satisfy the Radial Basis Function concept.
vec_tVector type used.

Usage example:

double h = 0.1;
std::vector<Vec2d> pts = {{0, 0}, {0, h}, {0, -h}, {h, 0}, {-h, 0}};
Eigen::VectorXd values(pts.size()); values << 0.1, 0.1, 0.1, 0.2, 0.0; // linear function in x
RBFFD<Gaussian<double>, Vec2d, ScaleToFarthest> approx(1.0, 2);
auto appr = approx.getApproximant({0, 0}, pts, values);
double val = appr({h/2, 0});
double lap = appr({h/2, 0}, Lap<2>());

Definition at line 34 of file RBFInterpolant_fwd.hpp.

Public Member Functions

 RBFInterpolant (const rbf_t &rbf, const Monomials< vec_t > &mon, const vector_t &point, const std::vector< vector_t > &support, scalar_t scale, const Eigen::Matrix< scalar_t, Eigen::Dynamic, 1 > &coefficients)
 Construct a RBF interpolant with known coefficients. More...
 
scalar_t operator() (const vector_t &point) const
 Evaluate the interpolant at given point. More...
 
template<typename operator_t >
scalar_t operator() (const vector_t &point, const operator_t&op) const
 Evaluate an operator applied to the interpolant at given point. More...
 
const vector_tpoint () const
 Get the center point. More...
 
scalar_t scale () const
 Get the scale. More...
 
const Eigen::Matrix< scalar_t, Eigen::Dynamic, 1 > & coefficients () const
 Get the coefficient vector. More...
 

Public Types

typedef vec_t::scalar_t scalar_t
 Scalar type. More...
 
typedef vec_t vector_t
 Vector type. More...
 
typedef RBFType rbf_t
 Radial basis function type. More...
 

Private Attributes

rbf_t rbf_
 RBF used for interpolation. More...
 
Monomials< vec_t > mon_
 Augmenting monomials. More...
 
vector_t point_
 Center point. More...
 
std::vector< vector_tsupport_
 Local scaled stencil points. More...
 
scalar_t scale_
 Scale. More...
 
Eigen::Matrix< scalar_t, Eigen::Dynamic, 1 > coefficients_
 Coefficients (for scaled fn.) More...
 

Constructor & Destructor Documentation

◆ RBFInterpolant()

template<typename RBFType , typename vec_t >
mm::RBFInterpolant< RBFType, vec_t >::RBFInterpolant ( const rbf_t rbf,
const Monomials< vec_t > &  mon,
const vector_t point,
const std::vector< vector_t > &  support,
scalar_t  scale,
const Eigen::Matrix< scalar_t, Eigen::Dynamic, 1 > &  coefficients 
)

Construct a RBF interpolant with known coefficients.

The interpolant is of the form

\[ u(p) = \sum_{i=1}^n \alpha_i \phi\left(\left\|\frac{p - p_i}{s}\right\|\right) + \sum_{i=1}^s \beta_i q_i\left(\frac{p - p_c}{s} \right), \]

where \(\phi\) is the RBF and \(q_i\) are monomials.

Parameters
rbfRadial basis function \(\phi\), must satisfy the Radial Basis Function concept.
monMonomial basis \(q_i\).
pointCenter point \(p_c\) of the RBF interpolation (used for monomial shift).
supportNonscaled stencil of the RBFFD approximation.
scaleScale \(s\) used in RBFFD computation, see Scale function concept.
coefficientsThe \(n+s\) coefficients \(\alpha_i, \beta_i\) of the approximation.

Definition at line 16 of file RBFInterpolant.hpp.

Member Function Documentation

◆ coefficients()

template<typename RBFType , typename vec_t >
const Eigen::Matrix<scalar_t, Eigen::Dynamic, 1>& mm::RBFInterpolant< RBFType, vec_t >::coefficients ( ) const
inline

Get the coefficient vector.

Definition at line 80 of file RBFInterpolant_fwd.hpp.

◆ operator()() [1/2]

template<typename RBFType , typename vec_t >
vec_t::scalar_t mm::RBFInterpolant< RBFType, vec_t >::operator() ( const vector_t point) const

Evaluate the interpolant at given point.

Definition at line 30 of file RBFInterpolant.hpp.

◆ operator()() [2/2]

template<typename RBFType , typename vec_t >
template<typename operator_t >
scalar_t mm::RBFInterpolant< RBFType, vec_t >::operator() ( const vector_t point,
const operator_t&  op 
) const

Evaluate an operator applied to the interpolant at given point.

◆ point()

template<typename RBFType , typename vec_t >
const vector_t& mm::RBFInterpolant< RBFType, vec_t >::point ( ) const
inline

Get the center point.

Definition at line 76 of file RBFInterpolant_fwd.hpp.

◆ scale()

template<typename RBFType , typename vec_t >
scalar_t mm::RBFInterpolant< RBFType, vec_t >::scale ( ) const
inline

Get the scale.

Definition at line 78 of file RBFInterpolant_fwd.hpp.

Member Data Documentation

◆ coefficients_

template<typename RBFType , typename vec_t >
Eigen::Matrix<scalar_t, Eigen::Dynamic, 1> mm::RBFInterpolant< RBFType, vec_t >::coefficients_
private

Coefficients (for scaled fn.)

Definition at line 46 of file RBFInterpolant_fwd.hpp.

◆ mon_

template<typename RBFType , typename vec_t >
Monomials<vec_t> mm::RBFInterpolant< RBFType, vec_t >::mon_
private

Augmenting monomials.

Definition at line 42 of file RBFInterpolant_fwd.hpp.

◆ point_

template<typename RBFType , typename vec_t >
vector_t mm::RBFInterpolant< RBFType, vec_t >::point_
private

Center point.

Definition at line 43 of file RBFInterpolant_fwd.hpp.

◆ rbf_

template<typename RBFType , typename vec_t >
rbf_t mm::RBFInterpolant< RBFType, vec_t >::rbf_
private

RBF used for interpolation.

Definition at line 41 of file RBFInterpolant_fwd.hpp.

◆ scale_

template<typename RBFType , typename vec_t >
scalar_t mm::RBFInterpolant< RBFType, vec_t >::scale_
private

Scale.

Definition at line 45 of file RBFInterpolant_fwd.hpp.

◆ support_

template<typename RBFType , typename vec_t >
std::vector<vector_t> mm::RBFInterpolant< RBFType, vec_t >::support_
private

Local scaled stencil points.

Definition at line 44 of file RBFInterpolant_fwd.hpp.


The documentation for this class was generated from the following files:
mm::sh::lap
static const shape_flags lap
Indicates to calculate laplace shapes.
Definition: shape_flags.hpp:24
mm::Vec2d
Vec< double, 2 > Vec2d
Convenience typedef for 2d vector of doubles.
Definition: Vec_fwd.hpp:34