Numerical Analysis II Assistant lecture
Physics Dept. Shrooq M. Azzo
Third stage 2021-2022
Chapter One
Interpolation and Polynomial Approximation
It means insertion or filling up intermediate terms of the series.
Suppose we are given the following values of y = f(x) for a set of values
of x:
𝑥: 𝑥0 𝑥1 𝑥2 … . 𝑥𝑛
𝑦: 𝑦0 𝑦1 𝑦2 … . 𝑦𝑛
Thus the process of finding the value of y corresponding to any value of
𝑥 = 𝑥𝑖 between 𝑥0 and 𝑥𝑛 is called interpolation. Hence interpolation is
the technique of estimating the value of a function for any intermediate
value of the independent variable, while the process of computing the
value of the function outside the given range is called extrapolation.
polynomial
Let f(x) be a continuous function defined on some interval [a, b], and be
prescribed at n + 1 distinct tabular points 𝑥0 , 𝑥1 , … , 𝑥𝑛 such that 𝑎 = 𝑥0 <
𝑥1 < 𝑥2 < ⋯ < 𝑥𝑛 = 𝑏. The distinct tabular points𝑥0 , 𝑥1 , … , 𝑥𝑛 may be
non-equispaced or equispaced, that is 𝑥𝑘+1 − 𝑥𝑘 = ℎ, 𝑘 = 0,1, … , 𝑛 − 1.
The problem of polynomial approximation is to find a polynomial 𝑃𝑛 (𝑥),
of degree ≤ n, which fits the given data exactly, that is,
𝑃𝑛 (𝑥𝑖 ) = 𝑓(𝑥𝑖 ), 𝑖 = 0,1,2, … 𝑛, where 𝑃𝑛 (𝑥) = 𝑎0 + 𝑎1 𝑥 + 𝑎0 𝑥 2 + ⋯ +
𝑎𝑛 𝑥 𝑛 , 𝑛 ≠ 0 then,
If n=1 the interpolation is called linear interpolation.
If n=2 the interpolation is called Quadratic interpolation.
If n=3 the interpolation is called cubic interpolation.
Note: the degree of the polynomial must be less one then the number of
points.
1
Numerical Analysis II Assistant lecture
Physics Dept. Shrooq M. Azzo
Third stage 2021-2022
Example: from the table find 𝑓(0), 𝑓(6) using square interpolation at the
point (−1,2,3), (2,3,4)
𝑋 -1 2 3 4
𝑌 = 𝐹(𝑥) 3 1 2 4
2
Numerical Analysis II Assistant lecture
Physics Dept. Shrooq M. Azzo
Third stage 2021-2022
Lagrange Interpolating Polynomials
Let the data
X 𝑥0 𝑥1 𝑥2 …. 𝑥𝑛
F(x) 𝑓(𝑥0 ) 𝑓(𝑥1 ) 𝑓(𝑥2 ) …. 𝑓(𝑥𝑛 )
be given at distinct unevenly spaced points or non-uniform points
𝑥0 , 𝑥1 , … , 𝑥𝑛 . This data may also be given at evenly spaced points , nth
Lagrange Interpolating Polynomial defined as
𝑃𝑛 (𝑥) = 𝑓(𝑥0 )𝐿0 (𝑥) + 𝑓(𝑥1 )𝐿1 (𝑥) + ⋯ + 𝑓(𝑥𝑛 )𝐿𝑛 (𝑥)
𝑛
= ∑ 𝑓(𝑥𝑖 )𝐿𝑖 (𝑥)
𝑖=0
Where 𝑓(𝑥𝑖 ) = 𝑓𝑖 and 𝐿𝑖 (𝑥), i = 0, 1, 2, …,n are polynomials of degree n
given from
(𝑥 − 𝑥0 )(𝑥 − 𝑥1 ) … (𝑥 − 𝑥𝑖−1 )(𝑥 − 𝑥𝑖+1 ) … (𝑥 − 𝑥𝑛 )
𝐿𝑖 (𝑥) =
(𝑥𝑖 − 𝑥0 )(𝑥𝑖 − 𝑥1 ) … (𝑥𝑖 − 𝑥𝑖−1 )(𝑥𝑖 − 𝑥𝑖+1 ) … (𝑥𝑖 − 𝑥𝑛 )
n
(𝑥 − 𝑥𝑗 )
=∏
𝑗=0
(𝑥𝑖 − 𝑥𝑗 )
𝑗≠𝑖
Example: from the table find f(2),f(5) and f(7) using Lagrange formula
𝑋 0 1 3
𝐹(𝑥) -1 2 7
3
Numerical Analysis II Assistant lecture
Physics Dept. Shrooq M. Azzo
Third stage 2021-2022
Lagrange Interpolating Polynomial Algorithm
INPUT for i=1,2,…n, 𝑥𝑖 , 𝑦𝑖 = 𝑓(𝑥𝑖 ), maximum number of iterations n, 𝑥 ∗
(𝑥−𝑥𝑗 )
Step1 compute 𝐿𝑖 (𝑥) = ∏n𝑗=0
(𝑥𝑖 −𝑥𝑗 )
𝑗≠𝑖
Step2 compute 𝑝(𝑥 ∗ ) = ∑𝑛𝑖=0 𝑓(𝑥𝑖 )𝐿𝑖 (𝑥)
Print 𝑓(𝑥 ∗ ) and stop.
Finite difference
Suppose that the function 𝑦 = 𝑓(𝑥) is tabulated for the equally spaced
values 𝑥 = 𝑥0 , 𝑥0 + ℎ, 𝑥0 + 2ℎ, … , 𝑥0 + 𝑛ℎ,giving
𝑦 = 𝑦0 , 𝑦1 , 𝑦2 , … , 𝑦𝑛 . To determine the values of f(x) or f ′(x) for some
intermediate values of x, the following three types of differences are useful:
1. Forward differences: The differences 𝑦1 − 𝑦0 , 𝑦2 − 𝑦1 , 𝑦3 −
𝑦2 , … , 𝑦𝑛 − 𝑦𝑛−1 , when denoted by ∆𝑦0 , ∆𝑦1 , ∆𝑦2 , … , ∆𝑦𝑛−1
respectively, are called the first forward differences where ∆ is the
forward difference operator.
Thus the first forward differences are
∆𝑦𝑖 = 𝑦𝑖+1 − 𝑦𝑖
Similarly, the second forward differences are defined by
∆2 𝑦𝑖 = ∆𝑦𝑖+1 − ∆𝑦𝑖
= (𝑦𝑖+2 − 𝑦𝑖+1 ) − (𝑦𝑖+1 − 𝑦𝑖 )
= 𝑦𝑖+2 − 2𝑦𝑖+1 + 𝑦𝑖
Similarly,
∆3 𝑦𝑖 = 𝑦𝑖+3 − 3𝑦𝑖+2 + 3𝑦𝑖+1 − 𝑦𝑖
In general,∆𝑘 𝑦𝑖 = ∆𝑘−1 𝑦𝑖+1 − ∆𝑘−1 𝑦𝑖 defines the 𝐾 𝑡ℎ forward
differences.
4
Numerical Analysis II Assistant lecture
Physics Dept. Shrooq M. Azzo
Third stage 2021-2022
Forward difference table
x y ∆𝑦 ∆2 𝑦 ∆3 𝑦 ∆4 𝑦 ∆5 𝑦
𝑥0 𝑦0
∆𝑦0
𝑥1 𝑦1 ∆2 𝑦0
∆𝑦1 ∆3 𝑦0
𝑥2 𝑦2 ∆2 𝑦1 ∆4 𝑦0
∆𝑦2 ∆3 𝑦1 ∆5 𝑦0
𝑥3 𝑦3 ∆2 𝑦2 ∆4 𝑦1
∆𝑦3 ∆3 𝑦2
𝑥4 𝑦4 ∆2 𝑦3
∆𝑦4
𝑥5 𝑦5
Example: find the forward difference from the table
x 2 3 4 5 6
F(x) 5 10 17 26 37
5
Numerical Analysis II Assistant lecture
Physics Dept. Shrooq M. Azzo
Third stage 2021-2022
2. Backward differences: The differences 𝑦1 − 𝑦0 , 𝑦2 − 𝑦1 , 𝑦3 −
𝑦2 , … , 𝑦𝑛 − 𝑦𝑛−1 , when denoted by ∇𝑦0 , ∇𝑦1 , ∇𝑦2 , … , ∇𝑦𝑛
respectively, are called the first forward differences where ∇ is the
backward difference operator
Similarly, we define higher order backward differences as,
∇𝑦𝑖 = 𝑦𝑖 − 𝑦𝑖−1
∇2 𝑦𝑖 = ∇𝑦𝑖 − ∇𝑦𝑖−1
= (𝑦𝑖 − 𝑦𝑖−1 ) − (𝑦𝑖−1 − 𝑦𝑖−2 )
= 𝑦𝑖 − 2𝑦𝑖−1 + 𝑦𝑖−2
In general, ∇𝑘 𝑦𝑖 = ∇𝑘−1 𝑦𝑖 − ∇𝑘−1 𝑦𝑖−1 defines the 𝐾 𝑡ℎ backward
differences.
Backward difference table
x y ∇𝑦 ∇2 𝑦 ∇3 𝑦 ∇4 𝑦 ∇5 𝑦
𝑥0 𝑦0
𝛻𝑦1
𝑥1 𝑦1 ∇2 𝑦2
𝛻𝑦2 ∇3 𝑦3
𝑥2 𝑦2 ∇2 𝑦3 ∇4 𝑦4
𝛻𝑦3 ∇3 𝑦4 ∇5 𝑦5
𝑥3 𝑦3 ∇2 𝑦4 ∇4 𝑦5
𝛻𝑦4 ∇3 𝑦5
𝑥4 𝑦4 ∇2 𝑦5
𝛻𝑦5
𝑥5 𝑦5
Example: find the backward difference from the table
x 2 3 4 5 6
F(x) 5 10 17 26 37
6
Numerical Analysis II Assistant lecture
Physics Dept. Shrooq M. Azzo
Third stage 2021-2022
3. Central difference: The central difference operator 𝛿 is
defined by the relations
𝑦1 − 𝑦0 = 𝛿𝑦1⁄2 , 𝑦2 − 𝑦1 = 𝛿𝑦3⁄2 , … , 𝑦𝑛 − 𝑦𝑛−1 = 𝛿𝑦𝑛−1
2
Similarly, high order central differences are defined as
𝛿𝑦3⁄2 − 𝛿𝑦1⁄2 = 𝛿 2 𝑦1 , 𝛿𝑦5⁄2 − 𝛿𝑦3⁄2 = 𝛿 2 𝑦2
and so on.
These differences are shown as follows:
x y Δ𝑦 δ2 𝑦 δ3 𝑦 δ4 𝑦 δ5 𝑦
𝑥0 𝑦0
𝛿 𝑦1⁄2
𝑥1 𝑦1 δ2 𝑦1
𝛿 𝑦3⁄2 δ3 𝑦3⁄2
𝑥2 𝑦2 δ2 𝑦2 δ4 𝑦2
𝛿 𝑦5⁄2 δ3 𝑦5⁄2 δ5 𝑦5⁄2
𝑥3 𝑦3 δ2 𝑦3 δ4 𝑦3
𝛿 𝑦7⁄2 δ3 𝑦7⁄2
𝑥4 𝑦4 δ2 𝑦4
𝛿 𝑦9⁄2
𝑥5 𝑦5
OTHER DIFFERENCE OPERATORS
1. Shift operator E: Shift operator E is the operation of increasing
the argument x by h so that
𝐸𝑓(𝑥) = 𝑓(𝑥 + ℎ)
𝐸 2 𝑓(𝑥) = 𝑓(𝑥 + 2ℎ) and so on. The inverse operator, 𝐸 −1 , is defined
by
7
Numerical Analysis II Assistant lecture
Physics Dept. Shrooq M. Azzo
Third stage 2021-2022
𝐸 −1 𝑓(𝑥) = 𝑓(𝑥 − ℎ).
Also𝐸 𝑛 𝑦𝑥 = 𝑦𝑥+𝑛ℎ .
2. Averaging operator 𝝁: The averaging operator is defined by
1
𝜇𝑦𝑥 = [𝑦𝑥+1ℎ + 𝑦𝑥−1ℎ ]
2 2 2
In difference calculus, E is the fundamental operator and ∇ , Δ, δ, μ can
be expressed in terms of E.
Relation between operators
1. 𝐸 = 1 + ∆
2. 𝐸∆= ∆𝐸
3. 𝐸𝛻 = ∇𝐸 = 𝐸 − 1 = ∆
4. 𝐸 −1 = 1 − ∇
5. ∆2 = 𝐸 2 − 2𝐸 + 1 = (𝐸 − 1)2
6. 𝐸 = 𝑒ℎ𝐷 ,D is differential operation
1 1
−
7. 𝛿 = 𝐸 − 𝐸 2 2
𝑈
8. 𝛿 = 2𝑠𝑖𝑛ℎ , U=hD
2
1 𝛿2
9. 𝐸 = 1 + 𝛿 2 ± 𝛿 √1 + = 𝑒𝑢
2 4
1 𝛿2
10. ∆= 𝛿 2 ± 𝛿 √1 +
2 4
= 𝑒𝑢 − 1
1 𝛿2
11.∇= − 𝛿2 ± 𝛿√1 + = 1 − 𝑒𝑢
2 4
1 1
12.𝛿 = ∆(1 + ∆)−2 = ∇(1 − ∇)−2
1
−
13.∇= 𝛿 𝐸 2
8
Numerical Analysis II Assistant lecture
Physics Dept. Shrooq M. Azzo
Third stage 2021-2022
Difference of a polynomial
If 𝑥0 , 𝑥1 , 𝑥2 , … , 𝑥𝑛 are given set of observations with common
difference h and let 𝑦0 , 𝑦1 , 𝑦2 , … , 𝑦𝑛 are their corresponding values,
where𝑦 = 𝑓(𝑥) be the given function, Suppose that we want to
interpolate 𝑓(𝑥) near the point 𝑥0 , Set 𝑥 = 𝑥𝑚 = 𝑥0 + 𝑚ℎ , 𝑓(𝑥𝑚 ) = 𝑓𝑚 .
Then
𝑓𝑚 = 𝑓(𝑥𝑚 ) = 𝑓(𝑥0 + 𝑚ℎ)
= 𝐸 𝑚 𝑓(𝑥0 ) = 𝐸 𝑚 𝑓0
∵ 𝐸 = 1 + ∆⇒ 𝐸 𝑚 = (1 + ∆)𝑚
∴ 𝑓𝑚 = 𝑓(𝑥𝑚 ) = (1 + ∆)𝑚 𝑓0
𝑚 𝑚 𝑚
𝑓𝑚 = [1 + ( ) ∆ + ( ) ∆2 + ⋯ + ( ) ∆𝑛 ] 𝑓0
1 2 𝑛
𝑚(𝑚−1) 𝑚(𝑚−1)(𝑚−2)…(𝑚−𝑛+1)
= 𝑓0 + 𝑚∆𝑓0 + ∆2 𝑓0 + ⋯ + ∆𝑛 𝑓0
2! 𝑛!
This is an alternate form of the Newton’s forward difference
interpolation formula.
And Let x be any point near 𝑥𝑛 , Set 𝑥 = 𝑥𝑚 = 𝑥𝑛 + 𝑚ℎ xn, Then,
∵ 𝑓𝑚 = 𝑓(𝑥𝑚 ) = 𝐸 𝑚 𝑓0 = (𝐸 −1 )−𝑚 𝑓0
𝑓𝑚 = (1 − ∇)−𝑚 𝑓0 = (1 + (−∇))−𝑚 𝑓0
−𝑚
= (∑−𝑚 𝑘
𝑘=0 ( 𝑘 ) (−∇) ) 𝑓0
𝑚(𝑚 + 1) 2
𝑓𝑚 = 𝑓0 + 𝑚∇𝑓0 + ∇ 𝑓0 + ⋯
2!
𝑚(𝑚 + 1)(𝑚 + 2) … (𝑚 + 𝑛 − 1) 𝑛
+ ∇ 𝑓0
𝑛!
This is an alternate form of the Newton’s backward difference
interpolation formula.
We shall study now the central difference formulae most suited for
interpolation near the middle of a tabulated set, the we have three formula
to find it
9
Numerical Analysis II Assistant lecture
Physics Dept. Shrooq M. Azzo
Third stage 2021-2022
1) Evert's formula
𝑚+1 2 𝑚+2 4 𝑚 + 𝑘 2𝑘
𝑓𝑚 ≅ 𝑚𝑓1 + ( ) 𝛿 𝑓1 + ( ) 𝛿 𝑓1 + ⋯ + ( ) 𝛿 𝑓1 + ⋯
3 5 2𝑘 + 1
𝑞+1 2 𝑞+2 2
+ 𝑞𝑓0 + ( ) 𝛿 𝑓0 + ( ) 𝛿 𝑓0 + ⋯
3 5
𝑞+𝑘
+( ) 𝛿 2𝑘 𝑓0
2𝑘 + 1
𝑥𝑚 −𝑥0
Where, 𝑚 ,𝑞 = 1 − 𝑚
ℎ
2) Bessel's formula
1 𝑚 1 1 𝑚
𝑓𝑚 = 𝜇𝑓1 + (𝑚 − ) 𝛿𝑓1 + ( ) 𝜇𝛿 2 𝑓1 + (𝑚 − ) ( ) + ⋯
2 2 2 2 2 3 2 2
𝑚+𝑘−1
+( ) 𝜇𝛿 2𝑘 𝑓1
2𝑘 2
1 1 𝑚 + 𝑘 − 1 2𝑘+1
+ (𝑚 − ) ( )𝛿 𝑓1
2𝑘 + 1 2 2𝑘 2
1
Where,𝜇𝛿𝑓𝑘 = [𝛿𝑓𝑘+1 + 𝛿𝑓𝑘−1 ]
2 2 2
3) Stirling's formula
𝑚2 2 𝑚(𝑚2 − 1) 3 𝑚2 (𝑚2 − 1) 4
𝑓𝑚 = 𝑓0 + 𝑚𝜇𝑓0 + 𝛿 𝑓0 + 𝜇𝛿 𝑓0 + 𝛿 𝑓0
2 3! 4!
+⋯
Example:- find 𝑓(4.5), 𝑓(9), 𝑓(6.4) from the table
x 4 6 8 10
F(x) 1 3 8 20
Sol.:- because 𝑥𝑚 = 4.5 we take 𝑥0 = 4
𝑥𝑚 −𝑥0 4.5−4 1
𝑚= = =
ℎ 2 4
The difference table is
10
Numerical Analysis II Assistant lecture
Physics Dept. Shrooq M. Azzo
Third stage 2021-2022
𝑥𝑖 𝑓𝑖 ∆𝑓𝑖 ∆2 𝑓𝑖 ∆3 𝑓𝑖
4 1
2
6 3 3
5 4
8 8 7
12
10 20
𝑚(𝑚 − 1) 𝑚(𝑚 − 1)(𝑚 − 2)
𝑓(4.5) ≅ 𝑓0 + 𝑚∆𝑓0 + ∆2 𝑓0 + ∆3 𝑓0
2! 3!
1 1 1 1 1
1 ( −1) ( −1)( −2)
= 1 + 4 ∗ 2 + 4 42 ∗3+ 4 4 4
∗ 4 = 1.4375
6
To find 𝑓(9), because 𝑥𝑚 = 9 we take 𝑥0 = 10
𝑥𝑚 −𝑥0 9−10 1
𝑚= = =−
ℎ 2 2
𝑚(𝑚 + 1) 𝑚(𝑚 + 1)(𝑚 + 2)
𝑓(9) ≅ 𝑓0 + 𝑚∇𝑓0 + ∇2 𝑓0 + ∇3 𝑓0
2! 3!
1 1 1 1 1
1 − (− +1) − (− +1)(− +2)
2 2 2
= 20 + (− ) ∗ 12 + 2 2
∗7+ ∗ 4 = 12.875
2 2 6
To find 𝑓(6.4), because 𝑥𝑚 = 6.4we take 𝑥0 = 6, 𝑥1 = 8
𝑥𝑚 −𝑥0 6.4−6
𝑚= = = 0.2
ℎ 2
𝑞 = 1 − 𝑚 = 0.8
𝑓(𝑥𝑚 ) = 𝑓(6.4)
𝑚+1 2 𝑚+2 4
≅ 𝑚𝑓1 + ( ) 𝛿 𝑓1 + ( ) 𝛿 𝑓1 + ⋯ + 𝑞𝑓0
3 5
𝑞+1 2 𝑚+2 4
+( ) 𝛿 𝑓0 + ( ) 𝛿 𝑓0
3 5
𝑚+1 2 𝑞+1 2
= 𝑚𝑓1 + ( ) 𝛿 𝑓1 + 𝑞𝑓0 + ( ) 𝛿 𝑓0
3 3
0.2 + 1 0.8 + 1
= 0.2 ∗ 8 + ( ) ∗ 7 + 0.8 ∗ 3 + ( )∗3
3 3
1.2(1.2−1)(1.2−2)(1.2−3)! 1.8(1.8−1)(1.8−2)(1.8−3)!
= 1.6 + ∗ 7 + 2.4 + ∗3
3!(1.2−3)! 3!(1.8−3)!
11
Numerical Analysis II Assistant lecture
Physics Dept. Shrooq M. Azzo
Third stage 2021-2022
= 3.632
12