<< . .

. 51
( : 95)



. . >>

uniformly converge to f and to its derivatives, as h tends to zero. The third
order derivative converges as well, provided that β is uniformly bounded.

Example 8.8 Figure 8.10 shows the cubic spline approximating the function in
the Runge™s example, and its ¬rst, second and third order derivatives, on a grid
of 11 equally spaced nodes, while in Table 8.5 the error s3 ’ f ∞ is reported as
a function of h together with the computed order of convergence p. The results
clearly demonstrate that p tends to 4 (the theoretical order) as h tends to zero.

8.6 Approximation by Splines 353

h 1 0.5 0.25 0.125 0.0625
s3 ’ f 0.022 0.0032 2.7741e-4 1.5983e-5 9.6343e-7

p “ 2.7881 3.5197 4.1175 4.0522
TABLE 8.5. Experimental interpolation error for Runge™s function using cubic
splines
1 0.8

0.9 (a)
0.6
(b)
0.8
0.4
0.7

0.2
0.6

0.5 0

0.4
’0.2

0.3
’0.4
0.2

’0.6
0.1

0 ’0.8
’5 ’4 ’3 ’2 ’1 0 1 2 3 4 5 ’5 ’4 ’3 ’2 ’1 0 1 2 3 4 5
5
1

4 (d)
(c)
0.5
3

2
0
1

0
’0.5

’1

’1
’2

’3
’1.5
’4

’5
’2
’5 ’4 ’3 ’2 ’1 0 1 2 3 4 5
’5 ’4 ’3 ’2 ’1 0 1 2 3 4 5



FIGURE 8.10. Interpolating spline (a) and its ¬rst (b), second (c) and third (d)
order derivatives (in solid line) for the function of Runge™s example (in dashed
line)

8.6.2 B-splines
Let us go back to splines of a generic degree k, and consider the B-spline
(or bell-spline) basis, referring to divided di¬erences introduced in Section
8.2.1.

De¬nition 8.2 The normalized B-spline Bi,k+1 of degree k relative to the
distinct nodes xi , . . . , xi+k+1 is de¬ned as
Bi,k+1 (x) = (xi+k+1 ’ xi )g[xi , . . . , xi+k+1 ], (8.51)
where
(t ’ x)k if x ¤ t,
g(t) = (t ’ x)k = (8.52)
+
0 otherwise.
354 8. Polynomial Interpolation



Substituting (8.18) into (8.51) yields the following explicit representation
k+1
(xj+i ’ x)k
+
Bi,k+1 (x) = (xi+k+1 ’ xi ) . (8.53)
k+1
j=0
(xi+j ’ xi+l )
l=0
l=j

From (8.53) it turns out that the active nodes of Bi,k+1 (x) are xi , . . . , xi+k+1
and that Bi,k+1 (x) is non null only within the interval [xi , xi+k+1 ].
Actually, it can be proved that it is the unique non null spline of min-
imum support relative to nodes xi , . . . , xi+k+1 [Sch67]. It can also be
(l) (l)
shown that Bi,k+1 (x) ≥ 0 [dB83] and |Bi,k+1 (xi )| = |Bi,k+1 (xi+k+1 )| for
l = 0, . . . , k ’1 [Sch81]. B-splines admit the following recursive formulation
([dB72], [Cox72])

1 if x ∈ [xi , xi+1 ],
Bi,1 (x) =
0 otherwise, (8.54)
x ’ xi xi+k+1 ’ x
Bi+1,k (x), k ≥ 1,
Bi,k+1 (x) = Bi,k (x) +
xi+k ’ xi xi+k+1 ’ xi+1
which is usually preferred to (8.53) when evaluating a B-spline at a given
point.

Remark 8.3 It is possible to de¬ne B-splines even in the case of partially
coincident nodes, by suitably extending the de¬nition of divided di¬erences.
This leads to a new recursive form of Newton divided di¬erences given by
(see for further details [Die93])
±
 f [x1 , . . . , xn ] ’ f [x0 , . . . , xn’1 ] if x < x < . . . < x

 0 1 n
xn ’ x0
f [x0 , . . . , xn ] =
 f (n+1) (x0 )

 if x0 = x1 = . . . = xn .
(n + 1)!
Assuming that m (with 1 < m < k + 2) of the k + 2 nodes xi , . . . , xi+k+1
are coincident and equal to », then (8.46) will contain a linear combination
k+1’j
of the functions (» ’ x)+ , for j = 1, . . . , m. As a consequence, the
B-spline can have continuous derivatives at » only up to order k ’ m and,
therefore, it is discontinuous if m = k + 1. It can be checked [Die93] that,
if xi’1 < xi = . . . = xi+k < xi+k+1 , then
±
 k
xi+k+1 ’ x

if x ∈ [xi , xi+k+1 ],
xi+k+1 ’ xi
Bi,k+1 (x) =

0 otherwise,
8.6 Approximation by Splines 355

while for xi < xi+1 = . . . = xi+k+1 < xi+k+2
±
 k
x ’ xi

if x ∈ [xi , xi+k+1 ],
xi+k+1 ’ xi
Bi,k+1 (x) =

0 otherwise.

Combining these formulae with the recursive relation (8.54) allows for con-
structing B-splines with coincident nodes.

Example 8.9 Let us examine the special case of cubic B-splines on equally
spaced nodes xi+1 = xi + h for i = 0, ..., n ’ 1. Equation (8.53) becomes

6h3 Bi,4 (x) =
±
 (x ’ xi ) , if x ∈ [xi , xi+1 ],
3



3
 h + 3h2 (x ’ x ) + 3h(x ’ x )2 ’ 3(x ’ x )3 , if x ∈ [x , x ],


 i+1 i+1 i+1 i+1 i+2


h3 + 3h2 (xi+3 ’ x) + 3h(xi+3 ’ x)2 ’ 3(xi+3 ’ x)3 , if x ∈ [xi+2 , xi+3 ],




 (xi+4 ’ x)3 ,
 if x ∈ [xi+3 , xi+4 ],





0 otherwise.

In Figure 8.11 the graph of Bi,4 is shown in the case of distinct nodes and of

partially coincident nodes.


1

0.9

0.8

0.7

0.6

0.5

0.4

0.3

0.2

0.1

0
’2 ’1 0 1 2


FIGURE 8.11. B-spline with distinct nodes (in solid line) and with three coin-
cident nodes at the origin (in dashed line). Notice the discontinuity of the ¬rst
derivative


Given n + 1 distinct nodes xj , j = 0, . . . , n, n ’ k linearly independent
B-splines of degree k can be constructed, though 2k degrees of freedom are
356 8. Polynomial Interpolation

still available to generate a basis for Sk . One way of proceeding consists of
introducing 2k ¬ctitious nodes

x’k ¤ x’k+1 ¤ . . . ¤ x’1 ¤ x0 = a,
(8.55)
b = xn ¤ xn+1 ¤ . . . ¤ xn+k

which the B-splines Bi,k+1 , with i = ’k, . . . , ’1 and i = n ’ k, . . . , n ’ 1,
are associated with. By doing so, any spline sk ∈ Sk can be uniquely written
as
n’1
sk (x) = ci Bi,k+1 (x). (8.56)
i=’k

The real numbers ci are the B-spline coe¬cients of sk . Nodes (8.55) are
usually chosen as coincident or periodic.
1. Coincident: this choice is suitable for enforcing the values attained
by a spline at the end points of its de¬nition interval. In such a case,
indeed, thanks to Remark 8.3 about B-splines with coincident nodes,
we get

sk (a) = c’k , sk (b) = cn’1 . (8.57)

2. Periodic, that is

x’i = xn’i ’ b + a, xi+n = xi + b ’ a, i = 1, . . . , k.

This choice is useful if the periodicity conditions (8.44) have to be
imposed.

Remark 8.4 (Inserting nodes) Using B-splines instead of cardinal spli-
nes is advantageous when handling, with a reduced computational e¬ort, a
given con¬guration of nodes for which a spline sk is known. In particular,
assume that the coe¬cients ci of sk (in form (8.56)) are available over the
nodes x’k , x’k+1 , . . . , xn+k , and that we wish to add to these a new node
x.
The spline sk ∈ Sk , de¬ned over the new set of nodes, admits the follow-
˜
ing representation with respect to a new B-spline basis Bi,k+1

n’1
sk (x) = di Bi,k+1 (x).
i=’k

The new coe¬cients di can be computed starting from the known coe¬-
cients ci using the following algorithm [Boe80]:
8.7 Splines in Parametric Form 357

let x ∈ [xj , xj+1 ); then, construct a new set of nodes {yi } such that

yi = xi for i = ’k, . . . , j, yj+1 = x,
yi = xi’1 for i = j + 2, . . . , n + k + 1;

de¬ne
±
for i = ’k, . . . , j ’ k,
1


y
j+1 ’ yi
for i = j ’ k + 1, . . . , j,
ωi =
 yi+k+1 ’ yi



0 for i = j + 1, . . . , n;

compute

di = ωi ci + (1 ’ ωi )ci i = ’k, ..., n ’ 1.

This algorithm has good stability properties and can be generalized to the
case where more than one node is inserted at the same time (see [Die93]).




8.7 Splines in Parametric Form
Using interpolating splines presents the following two drawbacks:

1. the resulting approximation is of good quality only if the function
f does not exhibit large derivatives (in particular, we require that
|f (x)| < 1 for every x). Otherwise, oscillating behaviors may arise
in the spline, as demonstrated by the example considered in Figure
8.12 which shows, in solid line, the cubic interpolating spline over the
following set of data (from [SL89])

xi 8.125 8.4 9 9.845 9.6 9.959 10.166 10.2
fi 0.0774 0.099 0.28 0.6 0.708 1.3 1.8 2.177

2. sk depends on the choice of the coordinate system. In fact, performing
a clockwise rotation of 36 degrees of the coordinate system in the
above example, would lead to the spline without spurious oscillations
reported in the boxed frame in Figure 8.12.
All the interpolation procedures considered so far depend on the cho-
sen Cartesian reference system, which is a negative feature if the
spline is used for a graphical representation of a given ¬gure (for in-
stance, an ellipse). Indeed, we would like such a representation to
be independent of the reference system, that is, to have a geometric
invariance property.
358 8. Polynomial Interpolation

2.5
’4.2

’4.4

2
’4.6

’4.8

1.5
’5

’5.2
6 7 8 9 10
1




0.5




0
8 8.5 9 9.5 10 10.5


FIGURE 8.12. Geometric noninvariance for an interpolating cubic spline s3 : the
set of data for s3 in the boxed frame is the same as in the main ¬gure, rotated
by 36 degrees. The rotation diminishes the slope of the interpolated curve and
eliminates any oscillation from s3 . Notice that resorting to a parametric spline
(dashed line) removes the oscillations in s3 without any rotation of the reference
system


A solution is provided by parametric splines, in which any component of the
curve, written in parametric form, is approximated by a spline function.
Consider a plane curve in parametric form P(t) = (x(t), y(t)), with t ∈
[0, T ], then take the set of the points in the plane of coordinates Pi =
(xi , yi ), for i = 0, . . . , n, and introduce a partition onto [0, T ]: 0 = t0 <
t 1 < . . . < tn = T .
Using the two sets of values {ti , xi } and {ti , yi } as interpolation data,
we obtain the two splines sk,x and sk,y , with respect to the independent
variable t, that interpolate x(t) and y(t), respectively. The parametric curve
Sk (t) = (sk,x (t), sk,y (t)) is called the parametric spline. Obviously, di¬erent
parameterizations of the interval [0, T ] yield di¬erent splines (see Figure
8.13).
A reasonable choice of the parameterization makes use of the length of
each segment Pi’1 Pi ,

(xi ’ xi’1 )2 + (yi ’ yi’1 )2 , i = 1, . . . , n.
li =
i
Setting t0 = 0 and ti = k=1 lk for i = 1, . . . , n, every ti represents the
cumulative length of the piecewise line that joins the points from P0 to
Pi . This function is called the cumulative length spline and approximates
satisfactorily even those curves with large curvature. Moreover, it can also
be proved (see [SL89]) that it is geometrically invariant.
Program 68 implements the construction of cumulative parametric cu-
bic splines in two dimensions (it can be easily generalized to the three-
8.7 Splines in Parametric Form 359



<< . .

. 51
( : 95)



. . >>