SSJ  3.3.1
Stochastic Simulation in Java
Public Member Functions | Static Public Member Functions | List of all members
BSpline Class Reference

Represents a B-spline with control points at \((X_i, Y_i)\). More...

Inheritance diagram for BSpline:
[legend]
Collaboration diagram for BSpline:
[legend]

Public Member Functions

 BSpline (final double[] x, final double[] y, final int degree)
 Constructs a new uniform B-spline of degree degree with control points at (x[i], y[i]). More...
 
 BSpline (final double[] x, final double[] y, final double[] knots)
 Constructs a new uniform B-spline with control points at (x[i], y[i]), and knot vector given by the array knots. More...
 
double [] getX ()
 Returns the \(X_i\) coordinates for this spline. More...
 
double [] getY ()
 Returns the \(Y_i\) coordinates for this spline. More...
 
double getMaxKnot ()
 Returns the knot maximal value. More...
 
double getMinKnot ()
 Returns the knot minimal value. More...
 
double [] getKnots ()
 Returns an array containing the knot vector \((t_0, t_{m-1})\). More...
 
BSpline derivativeBSpline ()
 Returns the derivative B-spline object of the current variable. More...
 
BSpline derivativeBSpline (int i)
 Returns the \(i\)th derivative B-spline object of the current variable; \(i\) must be less than the degree of the original B-spline. More...
 
double evaluate (final double u)
 Returns the value of the function evaluated at \(x\). More...
 
double integral (double a, double b)
 Computes (or estimates) the integral of the function over the interval \([a, b]\). More...
 
double derivative (double u)
 Computes (or estimates) the first derivative of the function at point x. More...
 
double derivative (double u, int n)
 Computes (or estimates) the \(n\)th derivative of the function at point x. More...
 
double evalX (double u)
 
double evalY (double u)
 

Static Public Member Functions

static BSpline createInterpBSpline (double[] x, double[] y, int degree)
 Returns a B-spline curve of degree degree interpolating the \((x_i, y_i)\) points [42] . More...
 
static BSpline createApproxBSpline (double[] x, double[] y, int degree, int hp1)
 Returns a B-spline curve of degree degree smoothing \((x_i, y_i)\), for \(i=0,…,n\) points. More...
 

Detailed Description

Represents a B-spline with control points at \((X_i, Y_i)\).

Let \(\mathbf{P_i}=(X_i, Y_i)\), for \(i=0,…,n-1\), be a control point and let \(t_j\), for \(j=0,…,m-1\) be a knot. A B-spline [42]  of degree \(p=m-n-1\) is a parametric curve defined as

\[ \mathbf{P(t)} = \sum_{i=0}^{n-1} N_{i, p}(t) \mathbf{P_i},\mbox{ for }t_p\le t\le t_{m-p-1}. \]

Here,

\begin{align*} N_{i, p}(t) & = \frac{t-t_i}{t_{i+p} - t_i}N_{i, p-1}(t) + \frac{t_{i+p+1} - t}{t_{i+p+1} - t_{i+1}}N_{i+1, p-1}(t) \\ N_{i, 0}(t) & = \left\{ \begin{array}{ll} 1 & \mbox{ for }t_i\le t\le t_{i+1}, \\ 0 \mbox{ elsewhere}. \end{array} \right. \end{align*}

This class provides methods to evaluate \(\mathbf{P(t)}=(X(t), Y(t))\) at any value of \(t\), for a B-spline of any degree \(p\ge1\). Note that the evaluate(double) method of this class can be slow, since it uses a root finder to determine the value of \(t^*\) for which \(X(t^*)=x\) before it computes \(Y(t^*)\).

Constructor & Destructor Documentation

◆ BSpline() [1/2]

BSpline ( final double []  x,
final double []  y,
final int  degree 
)

Constructs a new uniform B-spline of degree degree with control points at (x[i], y[i]).

The knots of the resulting B-spline are set uniformly from x[0] to x[n-1].

Parameters
xthe values of \(X\).
ythe values of \(Y\).
degreethe degree of the B-spline.

◆ BSpline() [2/2]

BSpline ( final double []  x,
final double []  y,
final double []  knots 
)

Constructs a new uniform B-spline with control points at (x[i], y[i]), and knot vector given by the array knots.

Parameters
xthe values of \(X\).
ythe values of \(Y\).
knotsthe knots of the B-spline.

Member Function Documentation

◆ createApproxBSpline()

static BSpline createApproxBSpline ( double []  x,
double []  y,
int  degree,
int  hp1 
)
static

Returns a B-spline curve of degree degree smoothing \((x_i, y_i)\), for \(i=0,…,n\) points.

The precision depends on the parameter \(hp1\): \(1 \le\mathtt{degree} \le hp1<n\), which represents the number of control points used by the new B-spline curve, minimizing the quadratic error

\[ L = \sum_{i=0}^n\left( \frac{Y_i - S_i(X_i)}{W_i}\right)^2. \]

This method uses the uniformly spaced method for interpolating points with a B-spline curve and a uniformed clamped knot vector, as described in http://www.cs.mtu.edu/~shene/COURSES/cs3621/NOTES/.

Parameters
xthe values of \(X\).
ythe values of \(Y\).
degreethe degree of the B-spline.
hp1the desired number of control points.
Returns
the B-spline curve.

◆ createInterpBSpline()

static BSpline createInterpBSpline ( double []  x,
double []  y,
int  degree 
)
static

Returns a B-spline curve of degree degree interpolating the \((x_i, y_i)\) points [42] .

This method uses the uniformly spaced method for interpolating points with a B-spline curve, and a uniformed clamped knot vector, as described in http://www.cs.mtu.edu/~shene/COURSES/cs3621/NOTES/.

Parameters
xthe values of \(X\).
ythe values of \(Y\).
degreethe degree of the B-spline.
Returns
the B-spline curve.

◆ derivative() [1/2]

double derivative ( double  x)

Computes (or estimates) the first derivative of the function at point x.

Parameters
xthe point to evaluate the derivative to.
Returns
the value of the derivative.

Implements MathFunctionWithFirstDerivative.

◆ derivative() [2/2]

double derivative ( double  x,
int  n 
)

Computes (or estimates) the \(n\)th derivative of the function at point x.

For \(n=0\), this returns the result of umontreal.ssj.functions.MathFunction.evaluate(double).

Parameters
xthe point to evaluate the derivate to.
nthe order of the derivative.
Returns
the resulting derivative.
Exceptions
IllegalArgumentExceptionif n is negative or 0.

Implements MathFunctionWithDerivative.

◆ derivativeBSpline() [1/2]

BSpline derivativeBSpline ( )

Returns the derivative B-spline object of the current variable.

Using this function and the returned object, instead of the derivative method, is strongly recommended if one wants to compute many derivative values.

Returns
the derivative B-spline of the current variable.

◆ derivativeBSpline() [2/2]

BSpline derivativeBSpline ( int  i)

Returns the \(i\)th derivative B-spline object of the current variable; \(i\) must be less than the degree of the original B-spline.

Using this function and the returned object, instead of the derivative method, is strongly recommended if one wants to compute many derivative values.

Parameters
ithe degree of the derivative.
Returns
the ith derivative.

◆ evaluate()

double evaluate ( final double  x)

Returns the value of the function evaluated at \(x\).

Parameters
xvalue at which the function is evaluated
Returns
function evaluated at x

Implements MathFunction.

◆ getKnots()

double [] getKnots ( )

Returns an array containing the knot vector \((t_0, t_{m-1})\).

Returns
the knot vector.

◆ getMaxKnot()

double getMaxKnot ( )

Returns the knot maximal value.

Returns
the \(Y_i\) coordinates.

◆ getMinKnot()

double getMinKnot ( )

Returns the knot minimal value.

Returns
the \(Y_i\) coordinates.

◆ getX()

double [] getX ( )

Returns the \(X_i\) coordinates for this spline.

Returns
the \(X_i\) coordinates.

◆ getY()

double [] getY ( )

Returns the \(Y_i\) coordinates for this spline.

Returns
the \(Y_i\) coordinates.

◆ integral()

double integral ( double  a,
double  b 
)

Computes (or estimates) the integral of the function over the interval \([a, b]\).

Parameters
athe starting point of the interval.
bthe ending point of the interval.
Returns
the value of the integral.

Implements MathFunctionWithIntegral.


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