Eigen::Spline< Scalar_, Dim_, Degree_ > Class Template Reference

A class representing multi-dimensional spline curves. More...

#include <Spline.h>

Public Types

enum  { Dimension = Dim_ }
 
enum  { Degree = Degree_ }
 
typedef Scalar_ Scalar
 
typedef SplineTraits< Spline >::PointType PointType
 The point type the spline is representing. More...
 
typedef SplineTraits< Spline >::KnotVectorType KnotVectorType
 The data type used to store knot vectors. More...
 
typedef SplineTraits< Spline >::ParameterVectorType ParameterVectorType
 The data type used to store parameter vectors. More...
 
typedef SplineTraits< Spline >::BasisVectorType BasisVectorType
 The data type used to store non-zero basis functions. More...
 
typedef SplineTraits< Spline >::BasisDerivativeType BasisDerivativeType
 The data type used to store the values of the basis function derivatives. More...
 
typedef SplineTraits< Spline >::ControlPointVectorType ControlPointVectorType
 The data type representing the spline's control points. More...
 

Public Member Functions

 Spline ()
 Creates a (constant) zero spline. For Splines with dynamic degree, the resulting degree will be 0. More...
 
template<typename OtherVectorType , typename OtherArrayType >
 Spline (const OtherVectorType &knots, const OtherArrayType &ctrls)
 Creates a spline from a knot vector and control points. More...
 
template<int OtherDegree>
 Spline (const Spline< Scalar, Dimension, OtherDegree > &spline)
 Copy constructor for splines. More...
 
const KnotVectorTypeknots () const
 Returns the knots of the underlying spline. More...
 
const ControlPointVectorTypectrls () const
 Returns the ctrls of the underlying spline. More...
 
PointType operator() (Scalar u) const
 Returns the spline value at a given site \(u\). More...
 
SplineTraits< Spline >::DerivativeType derivatives (Scalar u, DenseIndex order) const
 Evaluation of spline derivatives of up-to given order. More...
 
template<int DerivativeOrder>
SplineTraits< Spline, DerivativeOrder >::DerivativeType derivatives (Scalar u, DenseIndex order=DerivativeOrder) const
 Evaluation of spline derivatives of up-to given order. More...
 
SplineTraits< Spline >::BasisVectorType basisFunctions (Scalar u) const
 Computes the non-zero basis functions at the given site. More...
 
SplineTraits< Spline >::BasisDerivativeType basisFunctionDerivatives (Scalar u, DenseIndex order) const
 Computes the non-zero spline basis function derivatives up to given order. More...
 
template<int DerivativeOrder>
SplineTraits< Spline, DerivativeOrder >::BasisDerivativeType basisFunctionDerivatives (Scalar u, DenseIndex order=DerivativeOrder) const
 Computes the non-zero spline basis function derivatives up to given order. More...
 
DenseIndex degree () const
 Returns the spline degree. More...
 
DenseIndex span (Scalar u) const
 Returns the span within the knot vector in which u is falling. More...
 

Static Public Member Functions

static DenseIndex Span (typename SplineTraits< Spline >::Scalar u, DenseIndex degree, const typename SplineTraits< Spline >::KnotVectorType &knots)
 Computes the span within the provided knot vector in which u is falling. More...
 
static BasisVectorType BasisFunctions (Scalar u, DenseIndex degree, const KnotVectorType &knots)
 Returns the spline's non-zero basis functions. More...
 
static BasisDerivativeType BasisFunctionDerivatives (const Scalar u, const DenseIndex order, const DenseIndex degree, const KnotVectorType &knots)
 Computes the non-zero spline basis function derivatives up to given order. More...
 

Static Private Member Functions

template<typename DerivativeType >
static void BasisFunctionDerivativesImpl (const typename Spline< Scalar_, Dim_, Degree_ >::Scalar u, const DenseIndex order, const DenseIndex p, const typename Spline< Scalar_, Dim_, Degree_ >::KnotVectorType &U, DerivativeType &N_)
 

Private Attributes

KnotVectorType m_knots
 
ControlPointVectorType m_ctrls
 

Detailed Description

template<typename Scalar_, int Dim_, int Degree_>
class Eigen::Spline< Scalar_, Dim_, Degree_ >

A class representing multi-dimensional spline curves.

The class represents B-splines with non-uniform knot vectors. Each control point of the B-spline is associated with a basis function

\begin{align*} C(u) & = \sum_{i=0}^{n}N_{i,p}(u)P_i \end{align*}

Template Parameters
Scalar_The underlying data type (typically float or double)
Dim_The curve dimension (e.g. 2 or 3)
Degree_Per default set to Dynamic; could be set to the actual desired degree for optimization purposes (would result in stack allocation of several temporary variables).

Member Typedef Documentation

◆ BasisDerivativeType

template<typename Scalar_ , int Dim_, int Degree_>
typedef SplineTraits<Spline>::BasisDerivativeType Eigen::Spline< Scalar_, Dim_, Degree_ >::BasisDerivativeType

The data type used to store the values of the basis function derivatives.

◆ BasisVectorType

template<typename Scalar_ , int Dim_, int Degree_>
typedef SplineTraits<Spline>::BasisVectorType Eigen::Spline< Scalar_, Dim_, Degree_ >::BasisVectorType

The data type used to store non-zero basis functions.

◆ ControlPointVectorType

template<typename Scalar_ , int Dim_, int Degree_>
typedef SplineTraits<Spline>::ControlPointVectorType Eigen::Spline< Scalar_, Dim_, Degree_ >::ControlPointVectorType

The data type representing the spline's control points.

◆ KnotVectorType

template<typename Scalar_ , int Dim_, int Degree_>
typedef SplineTraits<Spline>::KnotVectorType Eigen::Spline< Scalar_, Dim_, Degree_ >::KnotVectorType

The data type used to store knot vectors.

◆ ParameterVectorType

template<typename Scalar_ , int Dim_, int Degree_>
typedef SplineTraits<Spline>::ParameterVectorType Eigen::Spline< Scalar_, Dim_, Degree_ >::ParameterVectorType

The data type used to store parameter vectors.

◆ PointType

template<typename Scalar_ , int Dim_, int Degree_>
typedef SplineTraits<Spline>::PointType Eigen::Spline< Scalar_, Dim_, Degree_ >::PointType

The point type the spline is representing.

◆ Scalar

template<typename Scalar_ , int Dim_, int Degree_>
typedef Scalar_ Eigen::Spline< Scalar_, Dim_, Degree_ >::Scalar

The spline curve's scalar type.

Member Enumeration Documentation

◆ anonymous enum

template<typename Scalar_ , int Dim_, int Degree_>
anonymous enum
Enumerator
Dimension 

The spline curve's dimension.

40 { Dimension = Dim_ };
@ Dimension
Definition: Spline.h:40

◆ anonymous enum

template<typename Scalar_ , int Dim_, int Degree_>
anonymous enum
Enumerator
Degree 

The spline curve's degree.

41 { Degree = Degree_ };
@ Degree
Definition: Spline.h:41

Constructor & Destructor Documentation

◆ Spline() [1/3]

template<typename Scalar_ , int Dim_, int Degree_>
Eigen::Spline< Scalar_, Dim_, Degree_ >::Spline ( )
inline

Creates a (constant) zero spline. For Splines with dynamic degree, the resulting degree will be 0.

66  : m_knots(1, (Degree == Dynamic ? 2 : 2 * Degree + 2)),
68  // in theory this code can go to the initializer list but it will get pretty
69  // much unreadable ...
70  enum { MinDegree = (Degree == Dynamic ? 0 : Degree) };
71  m_knots.template segment<MinDegree + 1>(0) = Array<Scalar, 1, MinDegree + 1>::Zero();
72  m_knots.template segment<MinDegree + 1>(MinDegree + 1) = Array<Scalar, 1, MinDegree + 1>::Ones();
73  }
KnotVectorType m_knots
Definition: Spline.h:218
ControlPointVectorType m_ctrls
Definition: Spline.h:219
const int Dynamic
Definition: Constants.h:25
double Zero
Definition: pseudosolid_node_update_elements.cc:35

References Eigen::Spline< Scalar_, Dim_, Degree_ >::Degree, Eigen::Dynamic, Eigen::Spline< Scalar_, Dim_, Degree_ >::m_knots, and oomph::PseudoSolidHelper::Zero.

◆ Spline() [2/3]

template<typename Scalar_ , int Dim_, int Degree_>
template<typename OtherVectorType , typename OtherArrayType >
Eigen::Spline< Scalar_, Dim_, Degree_ >::Spline ( const OtherVectorType &  knots,
const OtherArrayType &  ctrls 
)
inline

Creates a spline from a knot vector and control points.

Parameters
knotsThe spline's knot vector.
ctrlsThe spline's control point vector.
81 : m_knots(knots), m_ctrls(ctrls) {}
const ControlPointVectorType & ctrls() const
Returns the ctrls of the underlying spline.
Definition: Spline.h:98
const KnotVectorType & knots() const
Returns the knots of the underlying spline.
Definition: Spline.h:93

◆ Spline() [3/3]

template<typename Scalar_ , int Dim_, int Degree_>
template<int OtherDegree>
Eigen::Spline< Scalar_, Dim_, Degree_ >::Spline ( const Spline< Scalar, Dimension, OtherDegree > &  spline)
inline

Copy constructor for splines.

Parameters
splineThe input spline.
88 : m_knots(spline.knots()), m_ctrls(spline.ctrls()) {}

Member Function Documentation

◆ BasisFunctionDerivatives()

template<typename Scalar_ , int Dim_, int Degree_>
SplineTraits< Spline< Scalar_, Dim_, Degree_ > >::BasisDerivativeType Eigen::Spline< Scalar_, Dim_, Degree_ >::BasisFunctionDerivatives ( const Scalar  u,
const DenseIndex  order,
const DenseIndex  degree,
const KnotVectorType knots 
)
static

Computes the non-zero spline basis function derivatives up to given order.

The function computes

\begin{align*} \frac{d^i}{du^i} N_{i,p}(u), \hdots, \frac{d^i}{du^i} N_{i+p+1,p}(u) \end{align*}

with i ranging from 0 up to the specified order.

Parameters
uParameter \(u \in [0;1]\) at which the non-zero basis function derivatives are computed.
orderThe order up to which the basis function derivatives are computes.
degreeThe degree of the underlying spline
knotsThe underlying spline's knot vector.
472  {
473  typename SplineTraits<Spline>::BasisDerivativeType der;
474  BasisFunctionDerivativesImpl(u, order, degree, knots, der);
475  return der;
476 }
DenseIndex degree() const
Returns the spline degree.
Definition: Spline.h:271
static void BasisFunctionDerivativesImpl(const typename Spline< Scalar_, Dim_, Degree_ >::Scalar u, const DenseIndex order, const DenseIndex p, const typename Spline< Scalar_, Dim_, Degree_ >::KnotVectorType &U, DerivativeType &N_)
Definition: Spline.h:353

References constants::degree.

◆ basisFunctionDerivatives() [1/2]

template<typename Scalar_ , int Dim_, int Degree_>
SplineTraits< Spline< Scalar_, Dim_, Degree_ >, DerivativeOrder >::BasisDerivativeType Eigen::Spline< Scalar_, Dim_, Degree_ >::basisFunctionDerivatives ( Scalar  u,
DenseIndex  order 
) const

Computes the non-zero spline basis function derivatives up to given order.

The function computes

\begin{align*} \frac{d^i}{du^i} N_{i,p}(u), \hdots, \frac{d^i}{du^i} N_{i+p+1,p}(u) \end{align*}

with i ranging from 0 up to the specified order.

Parameters
uParameter \(u \in [0;1]\) at which the non-zero basis function derivatives are computed.
orderThe order up to which the basis function derivatives are computes.
453  {
454  typename SplineTraits<Spline<Scalar_, Dim_, Degree_> >::BasisDerivativeType der;
455  BasisFunctionDerivativesImpl(u, order, degree(), knots(), der);
456  return der;
457 }
SplineTraits< Spline >::BasisDerivativeType BasisDerivativeType
The data type used to store the values of the basis function derivatives.
Definition: Spline.h:56

References constants::degree.

◆ basisFunctionDerivatives() [2/2]

template<typename Scalar_ , int Dim_, int Degree_>
template<int DerivativeOrder>
SplineTraits<Spline, DerivativeOrder>::BasisDerivativeType Eigen::Spline< Scalar_, Dim_, Degree_ >::basisFunctionDerivatives ( Scalar  u,
DenseIndex  order = DerivativeOrder 
) const

Computes the non-zero spline basis function derivatives up to given order.

The function computes

\begin{align*} \frac{d^i}{du^i} N_{i,p}(u), \hdots, \frac{d^i}{du^i} N_{i+p+1,p}(u) \end{align*}

with i ranging from 0 up to the specified order.

Parameters
uParameter \(u \in [0;1]\) at which the non-zero basis function derivatives are computed.
orderThe order up to which the basis function derivatives are computes. Using the template version of this function is more efficieent since temporary objects are allocated on the stack whenever this is possible.

◆ BasisFunctionDerivativesImpl()

template<typename Scalar_ , int Dim_, int Degree_>
template<typename DerivativeType >
void Eigen::Spline< Scalar_, Dim_, Degree_ >::BasisFunctionDerivativesImpl ( const typename Spline< Scalar_, Dim_, Degree_ >::Scalar  u,
const DenseIndex  order,
const DenseIndex  p,
const typename Spline< Scalar_, Dim_, Degree_ >::KnotVectorType U,
DerivativeType &  N_ 
)
staticprivate
355  {
356  typedef Spline<Scalar_, Dim_, Degree_> SplineType;
357  enum { Order = SplineTraits<SplineType>::OrderAtCompileTime };
358 
359  const DenseIndex span = SplineType::Span(u, p, U);
360 
361  const DenseIndex n = (std::min)(p, order);
362 
363  N_.resize(n + 1, p + 1);
364 
367 
368  Matrix<Scalar, Order, Order> ndu(p + 1, p + 1);
369 
370  Scalar saved, temp; // FIXME These were double instead of Scalar. Was there a reason for that?
371 
372  ndu(0, 0) = 1.0;
373 
374  DenseIndex j;
375  for (j = 1; j <= p; ++j) {
376  left[j] = u - U[span + 1 - j];
377  right[j] = U[span + j] - u;
378  saved = 0.0;
379 
380  for (DenseIndex r = 0; r < j; ++r) {
381  /* Lower triangle */
382  ndu(j, r) = right[r + 1] + left[j - r];
383  temp = ndu(r, j - 1) / ndu(j, r);
384  /* Upper triangle */
385  ndu(r, j) = static_cast<Scalar>(saved + right[r + 1] * temp);
386  saved = left[j - r] * temp;
387  }
388 
389  ndu(j, j) = static_cast<Scalar>(saved);
390  }
391 
392  for (j = p; j >= 0; --j) N_(0, j) = ndu(j, p);
393 
394  // Compute the derivatives
395  DerivativeType a(n + 1, p + 1);
396  DenseIndex r = 0;
397  for (; r <= p; ++r) {
398  DenseIndex s1, s2;
399  s1 = 0;
400  s2 = 1; // alternate rows in array a
401  a(0, 0) = 1.0;
402 
403  // Compute the k-th derivative
404  for (DenseIndex k = 1; k <= static_cast<DenseIndex>(n); ++k) {
405  Scalar d = 0.0;
406  DenseIndex rk, pk, j1, j2;
407  rk = r - k;
408  pk = p - k;
409 
410  if (r >= k) {
411  a(s2, 0) = a(s1, 0) / ndu(pk + 1, rk);
412  d = a(s2, 0) * ndu(rk, pk);
413  }
414 
415  if (rk >= -1)
416  j1 = 1;
417  else
418  j1 = -rk;
419 
420  if (r - 1 <= pk)
421  j2 = k - 1;
422  else
423  j2 = p - r;
424 
425  for (j = j1; j <= j2; ++j) {
426  a(s2, j) = (a(s1, j) - a(s1, j - 1)) / ndu(pk + 1, rk + j);
427  d += a(s2, j) * ndu(rk + j, pk);
428  }
429 
430  if (r <= pk) {
431  a(s2, k) = -a(s1, k - 1) / ndu(pk + 1, r);
432  d += a(s2, k) * ndu(r, pk);
433  }
434 
435  N_(k, r) = static_cast<Scalar>(d);
436  j = s1;
437  s1 = s2;
438  s2 = j; // Switch rows
439  }
440  }
441 
442  /* Multiply through by the correct factors */
443  /* (Eq. [2.9]) */
444  r = p;
445  for (DenseIndex k = 1; k <= static_cast<DenseIndex>(n); ++k) {
446  for (j = p; j >= 0; --j) N_(k, j) *= r;
447  r *= p - k;
448  }
449 }
const unsigned n
Definition: CG3DPackingUnitTest.cpp:11
float * p
Definition: Tutorial_Map_using.cpp:9
SCALAR Scalar
Definition: bench_gemm.cpp:45
DenseIndex span(Scalar u) const
Returns the span within the knot vector in which u is falling.
Definition: Spline.h:279
SplineTraits< Spline >::BasisVectorType BasisVectorType
The data type used to store non-zero basis functions.
Definition: Spline.h:53
#define min(a, b)
Definition: datatypes.h:22
const Scalar * a
Definition: level2_cplx_impl.h:32
char char char int int * k
Definition: level2_impl.h:374
EIGEN_DEFAULT_DENSE_INDEX_TYPE DenseIndex
Definition: Meta.h:75
double U
Swimming speed.
Definition: two_d_variable_diff_adapt.cc:53
r
Definition: UniformPSDSelfTest.py:20
std::ptrdiff_t j
Definition: tut_arithmetic_redux_minmax.cpp:2

References a, j, k, min, n, p, UniformPSDSelfTest::r, RachelsAdvectionDiffusion::U, and oomph::PseudoSolidHelper::Zero.

◆ basisFunctions()

template<typename Scalar_ , int Dim_, int Degree_>
SplineTraits< Spline< Scalar_, Dim_, Degree_ > >::BasisVectorType Eigen::Spline< Scalar_, Dim_, Degree_ >::basisFunctions ( Scalar  u) const

Computes the non-zero basis functions at the given site.

Splines have local support and a point from their image is defined by exactly \(p+1\) control points \(P_i\) where \(p\) is the spline degree.

This function computes the \(p+1\) non-zero basis function values for a given parameter value \(u\). It returns

\begin{align*} N_{i,p}(u), \hdots, N_{i+p+1,p}(u) \end{align*}

Parameters
uParameter \(u \in [0;1]\) at which the non-zero basis functions are computed.
345  {
346  return Spline::BasisFunctions(u, degree(), knots());
347 }
static BasisVectorType BasisFunctions(Scalar u, DenseIndex degree, const KnotVectorType &knots)
Returns the spline's non-zero basis functions.
Definition: Spline.h:239

References Eigen::Spline< Scalar_, Dim_, Degree_ >::BasisFunctions(), and constants::degree.

◆ BasisFunctions()

template<typename Scalar_ , int Dim_, int Degree_>
Spline< Scalar_, Dim_, Degree_ >::BasisVectorType Eigen::Spline< Scalar_, Dim_, Degree_ >::BasisFunctions ( Scalar  u,
DenseIndex  degree,
const KnotVectorType knots 
)
static

Returns the spline's non-zero basis functions.

The function computes and returns

\begin{align*} N_{i,p}(u), \hdots, N_{i+p+1,p}(u) \end{align*}

Parameters
uThe site at which the basis functions are computed.
degreeThe degree of the underlying spline.
knotsThe underlying spline's knot vector.
241  {
242  const DenseIndex p = degree;
243  const DenseIndex i = Spline::Span(u, degree, knots);
244 
245  const KnotVectorType& U = knots;
246 
247  BasisVectorType left(p + 1);
248  left(0) = Scalar(0);
249  BasisVectorType right(p + 1);
250  right(0) = Scalar(0);
251 
252  VectorBlock<BasisVectorType, Degree>(left, 1, p) =
253  u - VectorBlock<const KnotVectorType, Degree>(U, i + 1 - p, p).reverse();
254  VectorBlock<BasisVectorType, Degree>(right, 1, p) = VectorBlock<const KnotVectorType, Degree>(U, i + 1, p) - u;
255 
256  BasisVectorType N(1, p + 1);
257  N(0) = Scalar(1);
258  for (DenseIndex j = 1; j <= p; ++j) {
259  Scalar saved = Scalar(0);
260  for (DenseIndex r = 0; r < j; r++) {
261  const Scalar tmp = N(r) / (right(r + 1) + left(j - r));
262  N[r] = saved + right(r + 1) * tmp;
263  saved = left(j - r) * tmp;
264  }
265  N(j) = saved;
266  }
267  return N;
268 }
int i
Definition: BiCGSTAB_step_by_step.cpp:9
static DenseIndex Span(typename SplineTraits< Spline >::Scalar u, DenseIndex degree, const typename SplineTraits< Spline >::KnotVectorType &knots)
Computes the span within the provided knot vector in which u is falling.
Definition: Spline.h:229
SplineTraits< Spline >::KnotVectorType KnotVectorType
The data type used to store knot vectors.
Definition: Spline.h:47
Scalar_ Scalar
Definition: Spline.h:39
@ N
Definition: constructor.cpp:22
Eigen::Matrix< Scalar, Dynamic, Dynamic, ColMajor > tmp
Definition: level3_impl.h:365

References constants::degree, i, j, N, p, UniformPSDSelfTest::r, Eigen::Spline< Scalar_, Dim_, Degree_ >::Span(), tmp, and RachelsAdvectionDiffusion::U.

Referenced by Eigen::Spline< Scalar_, Dim_, Degree_ >::basisFunctions().

◆ ctrls()

template<typename Scalar_ , int Dim_, int Degree_>
const ControlPointVectorType& Eigen::Spline< Scalar_, Dim_, Degree_ >::ctrls ( ) const
inline

Returns the ctrls of the underlying spline.

98 { return m_ctrls; }

References Eigen::Spline< Scalar_, Dim_, Degree_ >::m_ctrls.

◆ degree()

template<typename Scalar_ , int Dim_, int Degree_>
DenseIndex Eigen::Spline< Scalar_, Dim_, Degree_ >::degree

Returns the spline degree.

271  {
272  if (Degree_ == Dynamic)
273  return m_knots.size() - m_ctrls.cols() - 1;
274  else
275  return Degree_;
276 }

References Eigen::Dynamic.

◆ derivatives() [1/2]

template<typename Scalar_ , int Dim_, int Degree_>
SplineTraits< Spline< Scalar_, Dim_, Degree_ >, DerivativeOrder >::DerivativeType Eigen::Spline< Scalar_, Dim_, Degree_ >::derivatives ( Scalar  u,
DenseIndex  order 
) const

Evaluation of spline derivatives of up-to given order.

The function returns

\begin{align*} \frac{d^i}{du^i}C(u) & = \sum_{i=0}^{n} \frac{d^i}{du^i} N_{i,p}(u)P_i \end{align*}

for i ranging between 0 and order.

Parameters
uParameter \(u \in [0;1]\) at which the spline derivative is evaluated.
orderThe order up to which the derivatives are computed.
328  {
329  typename SplineTraits<Spline>::DerivativeType res;
330  derivativesImpl(*this, u, order, res);
331  return res;
332 }
cout<< "Here is the matrix m:"<< endl<< m<< endl;Matrix< ptrdiff_t, 3, 1 > res
Definition: PartialRedux_count.cpp:3
void derivativesImpl(const SplineType &spline, typename SplineType::Scalar u, DenseIndex order, DerivativeType &der)
Definition: Spline.h:299

References Eigen::derivativesImpl(), and res.

Referenced by check_global_interpolation_with_derivatives2d().

◆ derivatives() [2/2]

template<typename Scalar_ , int Dim_, int Degree_>
template<int DerivativeOrder>
SplineTraits<Spline, DerivativeOrder>::DerivativeType Eigen::Spline< Scalar_, Dim_, Degree_ >::derivatives ( Scalar  u,
DenseIndex  order = DerivativeOrder 
) const

Evaluation of spline derivatives of up-to given order.

The function returns

\begin{align*} \frac{d^i}{du^i}C(u) & = \sum_{i=0}^{n} \frac{d^i}{du^i} N_{i,p}(u)P_i \end{align*}

for i ranging between 0 and order.

Parameters
uParameter \(u \in [0;1]\) at which the spline derivative is evaluated.
orderThe order up to which the derivatives are computed. Using the template version of this function is more efficieent since temporary objects are allocated on the stack whenever this is possible.

◆ knots()

template<typename Scalar_ , int Dim_, int Degree_>
const KnotVectorType& Eigen::Spline< Scalar_, Dim_, Degree_ >::knots ( ) const
inline

Returns the knots of the underlying spline.

93 { return m_knots; }

References Eigen::Spline< Scalar_, Dim_, Degree_ >::m_knots.

Referenced by eval_spline3d_onbrks().

◆ operator()()

template<typename Scalar_ , int Dim_, int Degree_>
Spline< Scalar_, Dim_, Degree_ >::PointType Eigen::Spline< Scalar_, Dim_, Degree_ >::operator() ( Scalar  u) const

Returns the spline value at a given site \(u\).

The function returns

\begin{align*} C(u) & = \sum_{i=0}^{n}N_{i,p}P_i \end{align*}

Parameters
uParameter \(u \in [0;1]\) at which the spline is evaluated.
Returns
The spline value at the given location \(u\).
284  {
285  enum { Order = SplineTraits<Spline>::OrderAtCompileTime };
286 
287  const DenseIndex span = this->span(u);
288  const DenseIndex p = degree();
289  const BasisVectorType basis_funcs = basisFunctions(u);
290 
291  const Replicate<BasisVectorType, Dimension, 1> ctrl_weights(basis_funcs);
292  const Block<const ControlPointVectorType, Dimension, Order> ctrl_pts(ctrls(), 0, span - p, Dimension, p + 1);
293  return (ctrl_weights * ctrl_pts).rowwise().sum();
294 }
SplineTraits< Spline >::BasisVectorType basisFunctions(Scalar u) const
Computes the non-zero basis functions at the given site.
Definition: Spline.h:344

References constants::degree, and p.

◆ span()

template<typename Scalar_ , int Dim_, int Degree_>
DenseIndex Eigen::Spline< Scalar_, Dim_, Degree_ >::span ( Scalar  u) const

Returns the span within the knot vector in which u is falling.

Parameters
uThe site for which the span is determined.
279  {
280  return Spline::Span(u, degree(), knots());
281 }

References constants::degree, and Eigen::Spline< Scalar_, Dim_, Degree_ >::Span().

◆ Span()

template<typename Scalar_ , int Dim_, int Degree_>
DenseIndex Eigen::Spline< Scalar_, Dim_, Degree_ >::Span ( typename SplineTraits< Spline< Scalar_, Dim_, Degree_ > >::Scalar  u,
DenseIndex  degree,
const typename SplineTraits< Spline< Scalar_, Dim_, Degree_ > >::KnotVectorType knots 
)
static

Computes the span within the provided knot vector in which u is falling.

231  {
232  // Piegl & Tiller, "The NURBS Book", A2.1 (p. 68)
233  if (u <= knots(0)) return degree;
234  const Scalar* pos = std::upper_bound(knots.data() + degree - 1, knots.data() + knots.size() - degree - 1, u);
235  return static_cast<DenseIndex>(std::distance(knots.data(), pos) - 1);
236 }

References constants::degree.

Referenced by Eigen::Spline< Scalar_, Dim_, Degree_ >::BasisFunctions(), and Eigen::Spline< Scalar_, Dim_, Degree_ >::span().

Member Data Documentation

◆ m_ctrls

template<typename Scalar_ , int Dim_, int Degree_>
ControlPointVectorType Eigen::Spline< Scalar_, Dim_, Degree_ >::m_ctrls
private

Control points.

Referenced by Eigen::Spline< Scalar_, Dim_, Degree_ >::ctrls().

◆ m_knots

template<typename Scalar_ , int Dim_, int Degree_>
KnotVectorType Eigen::Spline< Scalar_, Dim_, Degree_ >::m_knots
private

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