Transcript Document
Chapter 15
Curve Fitting : Splines
Gab Byung Chae
Oscillation in a higher order
interpolation
The alternative : spline
Notation used to derive splines
n-1 intervals in n data points
Splines in
Linear spline
Quadratic spline
Cubic spline
st
1 ,
nd
2 ,
rd
3
order
15.2 Linear Splines
linear splines passing thru
( x1 , f 1 ),...., ( x n , f n ),
x1 x 2 ... x n 1 x n
f 2 f1
s
(
x
)
f
( x x1 ), x1 x x 2
1
1
x 2 x1
f i 1 f i
s ( x ) si ( x ) f i
( x x i ), x i x x i 1
x i 1 x i
f n f n 1
s ( x ) f n 1
( x x n 1 ), x n 1 x x n
n 1
x n x n 1
Example : (3,2.5),(4.5,1),(7,2.5),(9,.5) -> (5,?)
15.3 Quadratic Splines
The objective in quadratic splines is to derive
a second-order polynomial for each interval
between data points
The polynomial :
si(x) =ai + bi(x-xi)+ ci(x-xi)2
For n data points, there are n-1 intervals and,
consequently, 3(n-1) unknown constants to
evaluate. 3(n-1) equations or conditions
are required to evaluate the unknowns.
s1 ( x ) a1 b1 ( x x1 ) c1 ( x x1 ) 2 , x1 x x 2
2
s ( x ) s i ( x ) a i bi ( x x i ) c i ( x x i ) , x i x x i 1
s n 1 ( x ) a n 1 b n 1 ( x x n 1 ) c n 1 ( x x n 1 ) 2 , x n 1 x x n
Let hi x i 1 x i , for all i=1,,..,n-1
1. Continuity condition : the function must pass
through all the points.
si ( xi ) f i ai f i
--- (15.6)
for all i=1,,..,n-1 : (n-1) equations
2. The function values of adjacent polynomials
must be equal at the knots.
s i ( x i 1 ) f i 1 f i bi ( x i 1 x i ) c i ( x i 1 x i ) f i 1
2
2
f i bi hi c i hi f i 1
--- (15.8)
for all i=1,,..,n-1 : (n-1) equations
Let hi x i 1 x i , for all i=1,,..,n-1
3. The first derivatives at the interior nodes must
be equal.
(n-2) eq.s
s i ' ( x i 1 ) s i 1 ' ( x i 1 ) bi 2 c i hi bi 1
--- (15.9)
4. Assume that the second derivative is zero at
the first point.
s1 " ( x1 ) 0 c1 0
1 equation
Example : (3,2.5),(4.5,1),(7,2.5),(9,.5) -> (5,?)
Example 15.2
Problem] Fit quadratic splines to the data in
Table 15.1. Use the results to estimate the
value at x=5.
Solution] Equation (15.8) is written for i=1
through 3 (with c1=0) to give
f1+ b1h1 = f2
f2+ b2h2 + c2h22 = f3
f3+ b3h3 + c3h32 = f4
Equation (15.9) creates 2 conditions with
c1=0 :
b1 = b2
b2+ 2c2h2 = b3
f1
f2
f3
f4
=
=
=
=
2.5
1.0
2.5
0.5
h1 = 4.5-3.0 = 1.5
h2 = 7.0-4.5 = 2.5
h3 = 9.0-7.0 = 2.0
2.5 + 1.5 b1
= 1.0
1.0
+ 2.5 b2
+ 6.25 c2
= 2.5
2.5
+ 2.0 b3
+ 4c3 = 0.5
b1 - b2
=0
b2 - b3
+ 5c2
=0
b1 = -1
b2 = -1
b3= 2.2
c2 = 0.64
c3 = -1.6
s1(x) = 2.5 –(x-3)
s2(x) = 1.0 –(x-4.5)+0.64(x-4.5)2
s3(x) = 2.5 + 2.2(x-7.0)-1.6(5-7.0)2
x=5 이므로
s2(5) = 1.0 –(5-4.5)+0.64(5-4.5)2 =0.66
15.4 Cubic Splines
s1 ( x ) a 1 b1 ( x x1 ) c1 ( x x1 ) 2 d 1 ( x x1 ) 3 , x1 x x 2
2
3
s ( x ) s i ( x ) a i bi ( x x i ) c i ( x x i ) d i ( x x i ) , x i x x i 1
s n 1 ( x ) a n 1 b n 1 ( x x n 1 ) c n 1 ( x x n 1 ) 2 d n 1 ( x x n 1 ) 3 , x n 1 x x n
si ( xi ) f i a i f i
a.
2
(n-1) eq.s
3
s i ( x i 1 ) f i 1 f i bi hi c i hi d i hi f i 1
2
s i ' ( x i 1 ) s i 1 ' ( x i 1 ) bi 2 c i hi 3 d i hi bi 1
hi x i 1 x i ,
c.
(n-1) eq.s
(n-2) eq.s
d. (n-2) eq.s
2 c n 1 6 d n 1 h n 1 0 e. 2 eq.s
s i " ( x i 1 ) s i 1 " ( x i 1 ) c i 3 d i hi c i 1
s1 " ( x1 ) s n 1 " ( x n ) 0 2 c1
b.
Solving Cubic Splines
f. c i 3 d i h i c i 1 ( Eq . d ) d i ( c i 1 c i ) / 3 h i
2
3
g . f i bi hi c i hi d i hi f i bi hi
bi
f i 1 f i
hi
hi
3
hi
2
3
( 2 c i c i 1 ) f i 1 ( Eqs . b & f )
( 2 c i c i 1 )
2
h . b i 2 c i h i 3 d i h i b i h i ( c i c i 1 ) b i 1 ( Eqs . c & f )
i . h i 1 c i 1 2 ( h i 1 h i ) c i h i c i 1 3
j. 2 c1 0 ,
1
h
1
0
0
0
0
f i 1 f i
2 c n 1 6 d n 1 h n 1 2 c n 0
0
0
0
3
hi
f i f i 1
( Eqs . g & h )
h i 1
(e, f )
0
2 ( h1 h 2 )
h2
0
0
0
0
0
0
hn 2
2 ( h n 2 h n 1 )
0
0
0
0
0 c1
3
c2
0
0 :
0 :
fn
h n 1 c n 1
3
xn
1 c n
f3 f2
x3 x 2
0
3
f 2 f1
x 2 x1
:
f n 1
x n 1
:
3
0
f n 1 f n 2
x n 1 x n 2
Example 15.3
Fit cubic splines to the data . Utilize the
results to estimate the value at x =5
Solution : employ Eq.(15.27)
1
h
1
0
0
f1
f2
f3
f4
0
0
2 ( h1 h 2 )
h2
h2
2 ( h1 h 2 )
0
0
=
=
=
=
2.5
1.0
2.5
0.5
0 c1
0
0 c 2 3 ( f [ x 3 , x 2 ] f [ x 2 , x1 ])
h 3 c 3 3 ( f [ x 4 , x 3 ] f [ x 3 , x 2 ])
1 c4
0
h1 = 4.5-3.0 = 1.5
h2 = 7.0-4.5 = 2.5
h3 = 9.0-7.0 = 2.0
So we have
1
1 .5
0
0
0
0
8
2 .5
2 .5
9
0
0
0 c1 0
0 c 2 4 .8
2 c 3 4 .8
1 c 4 0
Using MATLAB the results are
c1 = 0
c2 = 0.839543726
c3 = -0.766539924
c4 = 0
Using Eq. (15.21) and
b1 = -1.419771863
b2 = -0.160456274
b3 = 0.022053232
(15.18),
d1 = 0.186565272
d2 = -0.214144487
d3 = 0.127756654
s1 ( x ) 2 . 5 1 . 419771863 ( x 3 ) 0 . 186565272 ( x 3 )
3
s 2 ( x ) 1 . 0 0 . 160456274 ( x 4 . 5 ) 0 . 839543726 ( x 4 . 5 )
0 . 214144487 ( x 4 . 5 )
3
s 3 ( x ) 2 . 5 0 . 022053232 ( x 7 . 0 ) 0 . 766539924 ( x 7 . 0 )
0 . 127756654 ( x 7 . 0 )
2
2
3
s 2 ( 5 ) 1 . 0 0 . 160456274 ( 5 4 . 5 ) 0 . 839543726 ( 5 4 . 5 )
0 . 214144487 ( 5 4 . 5 ) 1 . 102889734
3
2
15.4.2 End Conditions
Clamped End Condition : Specifying the first
derivatives at the first and last nodes.
Not-a-Knot End Condition : Force continuity
of the third derivative at the second and the
next-to-last knots. – same cubic functions
will apply to each of the first and last two
adjacent segments,
15.5 Piecewise interpolation in
Matlab
15.5.1 MATLAB function : spline
Syntax :
yy = spline(x, y, xx)
Example 15.4 Splines in MATLAB
Use MATLAB to fit 9 equally spaced data points
sampled from this function in the interval [-1, 1].
a) a not-a-knot spline b) a clamped spline with end
slopes of f1’=1 and f’n-1 = -4
f ( x)
1
1 25 x
2
Solution :
a)
>> x = linspace(-1,1,9);
>> y = 1./(1+25*x.^2);
>> xx = linspace(-1,1);
>> yy = spline(x,y,xx);
>> yr = 1./(1+25*xx.^2);
>> plot(x,y,’o’,xx,yy,xx,yr,’—’)
Cubic spline of the Runge
function
f ( x)
1
1 25 x
2
b)
>> yc =[1 y -4];
>> yyc =spline(x,yx,xx);
>> plot(x,y,’o’,xx,yyc,xx,yr,’—’)
15.5.2 MATLAB function : interp1
Implement a number of different types of
piecewise one-dimensional interpolation.
Syntax :
yi = interp1(x,y,xi, ‘method’)
yi = a vector containing the results of the
interpolation as evaluated at the points in the
vector xi , and ‘method’ = the desired method.
method
Nearest – nearest neighbor interpolation.
Linear – linear interpolation
Spline – piecewise cubic spline
interpolation(identical to the spline function)
pchip and cubic – piecewise cubic Hermite
interpolation.
The default is linear.
MATLAB’s methods
>> yy=spline(x,y,xx); cubic spilne
>> yy=interp1(x,y,xx): (piecewise) linear spline
>> yy=interp1(x,y,xx,’spline’): cubic spline
>> yy=interp1(x,y,xx,’nearest’) : nearest
neighbor method
>> yy=interp1(x,y,xx,’pchip’): piecewise cubic
Hermite interpolation
>> zz=interp2(x,y,z,xx,yy): 2 dimensional
interpolation
Example 15.5 Trade-offs Using
interpl
Linear, nearest, cubic spline, piecewise cubic
Hermite interpolation
>>
>>
>>
>>
>>
t = [ 0 20 40 56 68 80 84 96 104 110];
v = [0 20 20 38 80 80 100 100 125 125];
tt = linspace(0,110);
vl=interp1(t,v,tt); %(piecewise) linear spline
plot(t,v,’o’,tt,vl)
>> vn=interp1(t,v,tt,’nearest’); % nearest
>> plot(t,v,’o’,tt,vn)
>> vs=interp1(t,v,tt,’spline’); % cubic spline
>> plot(t,v,’o’,tt,vs)
>> vh=interp1(t,v,tt,’pchip’); % piecewise cubic
Hermite interpolation
>> plot(t,v,’o’,tt,vh)