Bembel
Spline

The Spline module provides basic routines related to spline function and local polynomials. More...

Classes

class  Bembel::Basis::PSpecificBasisHandler< P, Scalar >
 The functions above have a fixed compile time polynomial degree. The PSpecificBasis handler is used to convert this to a runtime p through template recursion. More...
 
struct  Bembel::static_assertKlessseqN< condition >
 Here we hardcode the binomial coefficients through template recursion and perform bounds checking, just in case... More...
 
class  Bembel::Basis::PSpecificShapeFunctionHandler< P >
 These routines implement a template recursion that allows to choose a compile time instantiation of a basis-evaluation routine with a runtime p. To replace the underlying basis, only these routines should be changed. More...
 

Functions

template<int P, typename Scalar >
void Bembel::Basis::phi_ (Eigen::Matrix< Scalar, Eigen::Dynamic, 1 > *c, Scalar w, double x)
 evaluates the 1D basis at x weighted with a quadrature weight w
 
template<int P, typename Scalar >
void Bembel::Basis::phi_dx_ (Eigen::Matrix< Scalar, Eigen::Dynamic, 1 > *c, Scalar w, double x)
 evaluates the derivative of phi
 
template<int P, typename Scalar >
void Bembel::Basis::phiphi_ (Eigen::Matrix< Scalar, Eigen::Dynamic, 1 > *c, Scalar w, Eigen::Vector2d a)
 evaluates the 2D tensor product basis at a point a in [0,1]^2
 
template<int P, typename Scalar >
void Bembel::Basis::phiphi_dx_ (Eigen::Matrix< Scalar, Eigen::Dynamic, 1 > *c, Scalar w, Eigen::Vector2d a)
 evaluates the x-derivative of the 2D tensor product basis at a point a in [0,1]^2
 
template<int P, typename Scalar >
void Bembel::Basis::phiphi_dy_ (Eigen::Matrix< Scalar, Eigen::Dynamic, 1 > *c, Scalar w, Eigen::Vector2d a)
 evaluates the y-derivative of the 2D tensor product basis at a point a in [0,1]^2
 
template<int P, typename Scalar >
void Bembel::Basis::Phi_times_Phi_ (Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic > *c, Scalar w, Eigen::Vector2d xi, Eigen::Vector2d eta)
 evaluates the interaction of two phiphis, one at xi and one at eta. Used for e.g. gram matrices.
 
template<int P, typename Scalar >
void Bembel::Basis::Div_Phi_times_Div_Phi_ (Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic > *c, Scalar weight, Eigen::Vector2d xi, Eigen::Vector2d eta)
 same as above, just using the divergence
 
template<int N>
constexpr double Bembel::Basis::BernsteinX (double evaluation_point) noexcept
 Template recursion to produce Bernstein polynomials. This is only limited by the binomial coefficient, see Pascal.hpp.
 
template<typename T >
Eigen::SparseMatrix< T > Bembel::Spl::MakeProjection (const std::vector< T > &x_knots, const std::vector< T > &y_knots, const std::vector< T > &x_unique_knots, const std::vector< T > &y_unique_knots, const int polynomial_degree_x, const int polynomial_degree_y) noexcept
 implements Bezier extraction via the solution of an interpolation problem. More...
 
template<typename T >
Eigen::Matrix< T, -1, -1 > Bembel::Spl::DeBoor (Eigen::Matrix< T, -1, -1 > const &control_points, const std::vector< double > &knot_vector, const std::vector< double > &evaluation_points) noexcept
 "By the book" implementations of the Cox-DeBoor formula. Inefficient, do not use at bottlenecks. More...
 
template<typename T >
Eigen::Matrix< T, -1, -1 > Bembel::Spl::DeBoorDer (Eigen::Matrix< T, -1, -1 > const &control_points, std::vector< double > const &knot, std::vector< double > const &evaluation_points) noexcept
 A "by the book" implementation of the derivatives and TP-algos based on the DeBoor Recursion.
 
std::vector< double > Bembel::Spl::MakeBezierKnotVector (int polynomial_degree) noexcept
 Here, routines for the creation and processing of knot vectors are defined.
 
template<typename T >
std::vector< T > Bembel::Spl::MakeInterpolationPoints (const std::vector< T > &uniq, const std::vector< T > &mask) noexcept
 Creates a T-vector of equidistant itnerpolation points for a given knot vector without any repetitions.
 
Eigen::Matrix< double, -1, -1 > Bembel::Spl::GetInterpolationMatrix (int polynomial_degree, const std::vector< double > &mask)
 returns the coefficients to represent a function in the Bernstein basis on [0,1].
 
template<typename T >
void Bembel::Spl::GetCoefficients (const int incr, const std::vector< double > &mask, const std::vector< T > &values, T *coefs)
 this solves a generic interpolation problem.
 
template<typename T >
std::vector< T > Bembel::Spl::GetCoefficients (const int increments, const std::vector< double > &mask, const std::vector< T > &values)
 this solves a generic interpolation problem.
 

Detailed Description

The Spline module provides basic routines related to spline function and local polynomials.

Function Documentation

◆ DeBoor()

template<typename T >
Eigen::Matrix<T, -1, -1> Bembel::Spl::DeBoor ( Eigen::Matrix< T, -1, -1 > const &  control_points,
const std::vector< double > &  knot_vector,
const std::vector< double > &  evaluation_points 
)
noexcept

"By the book" implementations of the Cox-DeBoor formula. Inefficient, do not use at bottlenecks.

ISO C++ forbids variable length array

Definition at line 22 of file DeBoor.hpp.

◆ MakeProjection()

template<typename T >
Eigen::SparseMatrix<T> Bembel::Spl::MakeProjection ( const std::vector< T > &  x_knots,
const std::vector< T > &  y_knots,
const std::vector< T > &  x_unique_knots,
const std::vector< T > &  y_unique_knots,
const int  polynomial_degree_x,
const int  polynomial_degree_y 
)
noexcept

implements Bezier extraction via the solution of an interpolation problem.

The inverse matrix, for the solition of the linear system to come.

For each basisfunction phi_j solves for the right coeef such that \Sum

Definition at line 23 of file Bezierextraction.hpp.