This page was generated from doc/euclidean/natural-non-uniform.ipynb. Interactive online version: .

Non-Uniform Natural Splines

The derivation is similar to the uniform case, but this time the parameter intervals can have arbitrary values.

import sympy as sp
sp.init_printing(order='grevlex')

from utility import NamedExpression

t = sp.symbols('t')

Just like in the uniform case, we are considering two adjacent spline segments, but this time we must allow arbitrary parameter values:

t3, t4, t5 = sp.symbols('t3:6')

b_monomial = sp.Matrix([t**3, t**2, t, 1]).T
b_monomial

$\displaystyle \left[\begin{matrix}t^{3} & t^{2} & t & 1\end{matrix}\right]$

coefficients3 = sp.symbols('a:dbm3')[::-1]
coefficients4 = sp.symbols('a:dbm4')[::-1]

b_monomial.dot(coefficients3)

$\displaystyle \boldsymbol{d}_{3} t^{3} + \boldsymbol{c}_{3} t^{2} + \boldsymbol{b}_{3} t + \boldsymbol{a}_{3}$

p3 = NamedExpression(
    'pbm3',
    b_monomial.dot(coefficients3).subs(t, (t - t3)/(t4 - t3)))
p4 = NamedExpression(
    'pbm4',
    b_monomial.dot(coefficients4).subs(t, (t - t4)/(t5 - t4)))
display(p3, p4)

$\displaystyle \boldsymbol{p}_{3} = \frac{\boldsymbol{d}_{3} \left(t - t_{3}\right)^{3}}{\left(- t_{3} + t_{4}\right)^{3}} + \frac{\boldsymbol{c}_{3} \left(t - t_{3}\right)^{2}}{\left(- t_{3} + t_{4}\right)^{2}} + \frac{\boldsymbol{b}_{3} \left(t - t_{3}\right)}{- t_{3} + t_{4}} + \boldsymbol{a}_{3}$

$\displaystyle \boldsymbol{p}_{4} = \frac{\boldsymbol{d}_{4} \left(t - t_{4}\right)^{3}}{\left(- t_{4} + t_{5}\right)^{3}} + \frac{\boldsymbol{c}_{4} \left(t - t_{4}\right)^{2}}{\left(- t_{4} + t_{5}\right)^{2}} + \frac{\boldsymbol{b}_{4} \left(t - t_{4}\right)}{- t_{4} + t_{5}} + \boldsymbol{a}_{4}$

pd3 = p3.diff(t)
pd4 = p4.diff(t)
display(pd3, pd4)

$\displaystyle \frac{d}{d t} \boldsymbol{p}_{3} = \frac{3 \boldsymbol{d}_{3} \left(t - t_{3}\right)^{2}}{\left(- t_{3} + t_{4}\right)^{3}} + \frac{\boldsymbol{c}_{3} \left(2 t - 2 t_{3}\right)}{\left(- t_{3} + t_{4}\right)^{2}} + \frac{\boldsymbol{b}_{3}}{- t_{3} + t_{4}}$

$\displaystyle \frac{d}{d t} \boldsymbol{p}_{4} = \frac{3 \boldsymbol{d}_{4} \left(t - t_{4}\right)^{2}}{\left(- t_{4} + t_{5}\right)^{3}} + \frac{\boldsymbol{c}_{4} \left(2 t - 2 t_{4}\right)}{\left(- t_{4} + t_{5}\right)^{2}} + \frac{\boldsymbol{b}_{4}}{- t_{4} + t_{5}}$

equations = [
    p3.evaluated_at(t, t3).with_name('xbm3'),
    p3.evaluated_at(t, t4).with_name('xbm4'),
    p4.evaluated_at(t, t4).with_name('xbm4'),
    p4.evaluated_at(t, t5).with_name('xbm5'),
    pd3.evaluated_at(t, t3).with_name('xbmdot3'),
    pd3.evaluated_at(t, t4).with_name('xbmdot4'),
    pd4.evaluated_at(t, t4).with_name('xbmdot4'),
    pd4.evaluated_at(t, t5).with_name('xbmdot5'),
]

We introduce a few new symbols to simplify the display, but we keep calculating with \(t_i\):

deltas = {
    t3: 0,
    t4: sp.Symbol('Delta3'),
    t5: sp.Symbol('Delta3') + sp.Symbol('Delta4'),
}

for e in equations:
    display(e.subs(deltas))

$\displaystyle \boldsymbol{x}_{3} = \boldsymbol{a}_{3}$

$\displaystyle \boldsymbol{x}_{4} = \boldsymbol{a}_{3} + \boldsymbol{b}_{3} + \boldsymbol{c}_{3} + \boldsymbol{d}_{3}$

$\displaystyle \boldsymbol{x}_{4} = \boldsymbol{a}_{4}$

$\displaystyle \boldsymbol{x}_{5} = \boldsymbol{a}_{4} + \boldsymbol{b}_{4} + \boldsymbol{c}_{4} + \boldsymbol{d}_{4}$

$\displaystyle \dot{\boldsymbol{x}}_{3} = \frac{\boldsymbol{b}_{3}}{\Delta_{3}}$

$\displaystyle \dot{\boldsymbol{x}}_{4} = \frac{\boldsymbol{b}_{3}}{\Delta_{3}} + \frac{2 \boldsymbol{c}_{3}}{\Delta_{3}} + \frac{3 \boldsymbol{d}_{3}}{\Delta_{3}}$

$\displaystyle \dot{\boldsymbol{x}}_{4} = \frac{\boldsymbol{b}_{4}}{\Delta_{4}}$

$\displaystyle \dot{\boldsymbol{x}}_{5} = \frac{\boldsymbol{b}_{4}}{\Delta_{4}} + \frac{2 \boldsymbol{c}_{4}}{\Delta_{4}} + \frac{3 \boldsymbol{d}_{4}}{\Delta_{4}}$

coefficients = sp.solve(equations, coefficients3 + coefficients4)

for c, e in coefficients.items():
    display(NamedExpression(c, e.subs(deltas)))

$\displaystyle \boldsymbol{d}_{3} = \Delta_{3} \dot{\boldsymbol{x}}_{3} + \Delta_{3} \dot{\boldsymbol{x}}_{4} + 2 \boldsymbol{x}_{3} - 2 \boldsymbol{x}_{4}$

$\displaystyle \boldsymbol{c}_{3} = - 2 \Delta_{3} \dot{\boldsymbol{x}}_{3} - \Delta_{3} \dot{\boldsymbol{x}}_{4} - 3 \boldsymbol{x}_{3} + 3 \boldsymbol{x}_{4}$

$\displaystyle \boldsymbol{b}_{3} = \Delta_{3} \dot{\boldsymbol{x}}_{3}$

$\displaystyle \boldsymbol{a}_{3} = \boldsymbol{x}_{3}$

$\displaystyle \boldsymbol{d}_{4} = - \Delta_{3} \dot{\boldsymbol{x}}_{4} - \Delta_{3} \dot{\boldsymbol{x}}_{5} + \dot{\boldsymbol{x}}_{4} \left(\Delta_{3} + \Delta_{4}\right) + \dot{\boldsymbol{x}}_{5} \left(\Delta_{3} + \Delta_{4}\right) + 2 \boldsymbol{x}_{4} - 2 \boldsymbol{x}_{5}$

$\displaystyle \boldsymbol{c}_{4} = 2 \Delta_{3} \dot{\boldsymbol{x}}_{4} + \Delta_{3} \dot{\boldsymbol{x}}_{5} - 2 \dot{\boldsymbol{x}}_{4} \left(\Delta_{3} + \Delta_{4}\right) - \dot{\boldsymbol{x}}_{5} \left(\Delta_{3} + \Delta_{4}\right) - 3 \boldsymbol{x}_{4} + 3 \boldsymbol{x}_{5}$

$\displaystyle \boldsymbol{b}_{4} = - \Delta_{3} \dot{\boldsymbol{x}}_{4} + \dot{\boldsymbol{x}}_{4} \left(\Delta_{3} + \Delta_{4}\right)$

$\displaystyle \boldsymbol{a}_{4} = \boldsymbol{x}_{4}$

pdd3 = pd3.diff(t)
pdd4 = pd4.diff(t)
display(pdd3, pdd4)

$\displaystyle \frac{d^{2}}{d t^{2}} \boldsymbol{p}_{3} = \frac{3 \boldsymbol{d}_{3} \left(2 t - 2 t_{3}\right)}{\left(- t_{3} + t_{4}\right)^{3}} + \frac{2 \boldsymbol{c}_{3}}{\left(- t_{3} + t_{4}\right)^{2}}$

$\displaystyle \frac{d^{2}}{d t^{2}} \boldsymbol{p}_{4} = \frac{3 \boldsymbol{d}_{4} \left(2 t - 2 t_{4}\right)}{\left(- t_{4} + t_{5}\right)^{3}} + \frac{2 \boldsymbol{c}_{4}}{\left(- t_{4} + t_{5}\right)^{2}}$

sp.Eq(pdd3.expr.subs(t, t4), pdd4.expr.subs(t, t4))

$\displaystyle \frac{3 \boldsymbol{d}_{3} \left(- 2 t_{3} + 2 t_{4}\right)}{\left(- t_{3} + t_{4}\right)^{3}} + \frac{2 \boldsymbol{c}_{3}}{\left(- t_{3} + t_{4}\right)^{2}} = \frac{2 \boldsymbol{c}_{4}}{\left(- t_{4} + t_{5}\right)^{2}}$

_.subs(coefficients).subs(deltas).simplify()

$\displaystyle \frac{2 \left(\Delta_{3} \dot{\boldsymbol{x}}_{3} + 2 \Delta_{3} \dot{\boldsymbol{x}}_{4} + 3 \boldsymbol{x}_{3} - 3 \boldsymbol{x}_{4}\right)}{\Delta_{3}^{2}} = \frac{2 \left(- 2 \Delta_{4} \dot{\boldsymbol{x}}_{4} - \Delta_{4} \dot{\boldsymbol{x}}_{5} - 3 \boldsymbol{x}_{4} + 3 \boldsymbol{x}_{5}\right)}{\Delta_{4}^{2}}$

Like in the uniform case, we can generalize by renaming index \(4\) to \(i\):

\begin{equation*} \frac{1}{\Delta_{i-1}} \dot{\boldsymbol{x}}_{i-1} + \left(\frac{2}{\Delta_{i-1}} + \frac{2}{\Delta_i}\right) \dot{\boldsymbol{x}}_i + \frac{1}{\Delta_i} \dot{\boldsymbol{x}}_{i+1} = \frac{3 (\boldsymbol{x}_i - \boldsymbol{x}_{i-1})}{{\Delta_{i-1}}^2} + \frac{3 (\boldsymbol{x}_{i+1} - \boldsymbol{x}_i)}{{\Delta_i}^2} \end{equation*}

We are not showing the full matrix here, because it would be quite a bit more complicated and less helpful than in the uniform case.

End Conditions

Like in the uniform case, we can come up with a few end conditions in order to define the missing matrix rows.

The Python class splines.Natural uses “natural” end conditions by default.

“Natural” end conditions are derived in a separate notebook, yielding these expressions:

\begin{align*} 2 \Delta_0 \dot{\boldsymbol{x}}_0 + \Delta_0 \dot{\boldsymbol{x}}_1 &= 3 (\boldsymbol{x}_1 - \boldsymbol{x}_0) \\ \Delta_{N-2} \dot{\boldsymbol{x}}_{N-2} + 2 \Delta_{N-2} \dot{\boldsymbol{x}}_{N-1} &= 3 (\boldsymbol{x}_{N-1} - \boldsymbol{x}_{N-2}) \end{align*}

Other end conditions can be derived as shown in the notebook about uniform “natural” splines.