Numerical Analysis

Vijay Kumar

Department of Mathematics and Statistics

DDU Gorakhpur University, Gorakhpur
[email protected]

Sept. 2018
Numerical Analysis Least Squares Interpolation Formulae Integration

Outline

Method of Least Squares

Fitting of Straight line
Fitting of Non-linear curves
Interpolation
Different Operators
Interpolation formulae for Equal Intervals
1 Newton’s forward
2 Newton’s backward
Interpolation formulae for Unequal Intervals
1 Lagrange’s
2 Newton’s divided

Numerical Integration
1 Trapezoidal rule.
2 Simpson’s 1/3 rule.
3 Simpson’s 3/8 rule.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 2 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration

...

Why Numerical Methods ?

The complex mathematical problems faced by
scientists and engineers rarely can be solved by
analytical approaches, so numerical methods
are often necessary.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 3 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration

Computational problems:

Interpolation and curve fitting

Numerical integration
Numerical differentiation
Roots of an equation
Numerical optimization
Matrix algebra
Solution of systems of linear algebraic equations
Solution of ordinary differential equations
Solution of partial differential equations
Computationally intensive Monte Carlo methods etc.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 4 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration

Numerical Analysis:

The statement
y = f (x) ; x0 ≤ x ≤ xn ,
means: corresponding to every value of x in the range x0 ≤ x ≤ xn , there exists
one or more values of y. Assuming that f (x) is single valued and continuous
and that it is known explicitly, then the values of f (x) corresponding to certain
given values of x, say x0 ≤ x ≤ xn can easily be computed and tabulated.

The central problem of numerical analysis is the converse one:

Given the set of tabular values (x0 , y0 ) , (x1 , y1 ) , (x2 , y2 ) , . . . , (xn , yn )

satisfying the relation y = f (x) where the nature of f (x) is not known, it is
required to find a simpler function φ (x), such that

f (x) and φ (x) agree at the set of tabulated points.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 5 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration

Numerical Analysis:

Numerical analysis always utilizes a discrete set of points

to represent functions.
Numerical method allows operations such as integration
and differentiation to be performed using discrete points.

No numerical method is completely trouble free in all

situations!
No numerical method is error free!
No numerical method is optimal for all types/forms of
an equation!

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 6 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Straight line Polynomial Nonlinear curves

Method of Least Squares:

Let the set of data points be (xi , yi ) ; i = 1, 2, . . . , n and let the curve given
by y = f (x) be fitted to this data.
At x = xi , the observed value of the ordinate is yi and corresponding value on
the fitting curve is f (xi ).
If ei is the residual of approximation at x = xi , then we have

ei = yi − f ( xi )

If we write
n
X n
X 2
S= e2i = (yi − f ( xi ))
i=1 i=1

then the method of least squares consists in minimizing S, i.e. the sum of
squares of residuals.
We shall study the least square fitting to the given data (xi , yi ) ; i = 1, 2, . . . , n
for different curves (linear and nonlinear etc.).

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 7 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Straight line Polynomial Nonlinear curves

Fitting of straight line y = a + b x:

Define the vertical difference

Let the set of data points be (xi , yi ); between P and the estimated line by
i = 1, 2, . . . , n. Suppose we wish to fit
the straight line y = a + b x to the given ei = (yi − (a + b xi ))
data set.
At x = xi , the estimated value of yi as y
obtained by the fitted straight line will be P(xi , yi)
ŷi = a + b xi ; i = 1, 2, . . . , n x
Residual a +b
Take any point P with coordinates, say y=

(xi , yi ). Drop a perpendicular from P M(xi , a + b xi)

to the x-axis. Intersecting our estimated
line at M with coordinates (xi , a + b xi ).
Then
O xi T x
OT = xi , P T = yi M T = a + b xi

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 8 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Straight line Polynomial Nonlinear curves

Scatter diagram:

P(xi , yi)
x
Residual a +b
y=

M(xi , a + b xi)

O xi T x

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 9 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Straight line Polynomial Nonlinear curves

Fitting of straight line contd...

These residuals, or deviations, from the estimated line will be positive or

negative as the actual point lies above or below the line.
If they are squared and summed, the resultant quantity must be nonnegative
and will vary directly with the spread of the points from the line.
Different pairs of values of a and b will give different lines and hence different
values for the sum of the squared residuals about the line. Thus we have
2
e2i = (yi − a − b xi )
n n
2
e2i =
P P
E= (yi − a − b xi )
i=1 i=1

The principle of least squares is that the a and b values should be chosen so
n
e2i as small as possible. A necessary condition is that the
P
as to make E =
i=1
partial derivatives of the sum with respect to a and b should both be zero. We
thus have

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 10 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Straight line Polynomial Nonlinear curves

Fitting of straight line contd...

n
∂E
P
∂a = −2 (yi − a − b xi ) = 0
i=1
n
∂E
P
∂b = −2 (yi − a − b xi ) xi = 0
i=1

Simplifying these two equations gives the normal equations for straight line :
n
P n
P
yi = n a + b xi
i=1 i=1
n n n
x2i
P P P
xi yi = a xi + b
i=1 i=1 i=1

Since the xi and yi are known quantities, the normal equations, can be solved.
n n n n
x2i and
P P P P
When these values yi , xi , xi yi computed from the data
i=1 i=1 i=1 i=1
are inserted, they give two simultaneous equations, which may be solved for a
and b.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 11 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Straight line Polynomial Nonlinear curves

Fitting of Nonlinear curves:

Fitting of second degree polynomial(Parabola):

Let the set of data points be (xi , yi ) ; i = 1, 2, . . . , n. Suppose we wish to fit
the second degree polynomial (parabola) y = a + b x + c x2 to the given data.
At x = xi , the estimated value of yi as obtained by the fitted parabola will be
ŷi = a + b xi + c x2i ; i = 1, 2, . . . , n



ei = yi − a + b xi + c x2i
n n
2
e2i = yi − a − b xi − c x2i
P P
E=
i=1 i=1

Differentiating w.r.t. a, b and c, we have

n
∂E


yi − a − b xi − c x2i = 0
P
∂a = −2
i=1
n
∂E


yi − a − b xi − c x2i xi = 0
P
∂b = −2
i=1
n
∂E


yi − a − b xi − c x2i x2i = 0
P
∂c = −2
i=1

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 12 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Straight line Polynomial Nonlinear curves

Fitting of Parabola contd...

The normal equations are given by

n n n
x2i
P P P
yi = n a + b xi + c
i=1 i=1 i=1
n n n n
x2i + c x3i
P P P P
xi yi = a xi + b
i=1 i=1 i=1 i=1
n n n n
x2i yi = a x2i + b x3i + c x4i
P P P P
i=1 i=1 i=1 i=1

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 13 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Straight line Polynomial Nonlinear curves

Semi-log transformation y = a eb x

Let the set of data points be (xi , yi ) ; i = 1, 2, . . . , n. Suppose we wish to fit

the curve
y = a eb x
Taking log, we have
log y = log a + (b log10 e) x
u = A + B x where u = log y ; A = log a and B = b log10 e
which a straight line. Therefore, the normal equations are given by
n
P n
P
ui = n A + B xi
i=1 i=1
n n n
x2i
P P P
xi ui = A xi + B
i=1 i=1 i=1

Solve for A and B. Hence

a = antilog (A) and b = B/ (b log10 e)

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 14 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Straight line Polynomial Nonlinear curves

Double-log transformation y = a xb
Let the set of data points be (xi , yi ) ; i = 1, 2, . . . , n. Suppose we wish to fit
the curve y = a xb . Taking log of both sides, we have

log y = log a + b log x

u = A + bv
where u = log y ; A = log a and v = log x. Therefore, the normal
equations are given by
n
P n
P
ui = n A + b vi
i=1 i=1
n n n
vi2
P P P
vi ui = A vi + b
i=1 i=1 i=1

n n n n
vi2 and
P P P P
Compute ui , vi , vi ui , substitute these values in normal
i=1 i=1 i=1 i=1
equations and solve for A and b. Hence

a = antilog (A)
V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 15 / 84
Numerical Analysis Least Squares Interpolation Formulae Integration Straight line Polynomial Nonlinear curves

Semi-log transformation y = a bx
Let the set of data points be (xi , yi ) ; i = 1, 2, . . . , n. Suppose we wish to fit
the curve
y = a bx
Taking log, we have

log y = log a + x log b or u=A+Bx

where u = log y ; A = log a and B = log b. Therefore, the normal

equations are given by
n
P n
P
ui = n A + B xi
i=1 i=1
n n n
x2i
P P P
xi ui = A xi + B
i=1 i=1 i=1

n n n n
x2i and
P P P P
Compute ui , xi , xi ui , substitute these values in normal
i=1 i=1 i=1 i=1
equations and solve for A and B. Hence

a = antilog (A)
V. Kumar, DDU Gorakhpur University and Methods
B.Sc.-I : Numerical b = antilog (B)
Sept. 2018 16 / 84
Numerical Analysis Least Squares Interpolation Formulae Integration Straight line Polynomial Nonlinear curves

b
Reciprocal transformation y = a + x

Let the set of data points be (xi , yi ) ; i = 1, 2, . . . , n. Suppose we wish to fit

the curve y = a + (b/x). At x = xi , the estimated value of yi as obtained by
the fitted curve will be ŷi = a + xbi ; i = 1, 2, . . . , n. The sum of residuals is
n n  2
e2i = yi − a − xbi
P P
given by E =
i=1 i=1
Differentiating w.r.t. a and b and equating to zero, we have
n  
∂E b
P
∂a = −2 yi − a − xi =0
i=1 
n  
∂E b 1
P
∂b = −2 yi − a − xi xi =0
i=1

The normal equations are given by

n
P n
P
yi = na + b (1/xi )
i=1 i=1
n n n

1/x2i
P P P
(yi /xi ) = a (1/xi ) + b
i=1 i=1 i=1

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 17 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Straight line Polynomial Nonlinear curves

b
Fitting of the curve y = a + x2

Let the set of data points be (xi , yi ) ; i = 1, 2, . . . , n. Suppose we wish to

fit the given curve. At x = xi the curve will give the estimated value of yi as
ŷi = f (xi ) = a + (b/x2i ) ; i = 1, 2, . . . , n which leads to residual sum of squares
n n  2
X X b
E= e2i = yi − a − 2
i=1 i=1
xi

Differentiating w.r.t. a and b and equating to zero, we have

n  
∂E
yi − a − xb2 = 0
P
∂a = −2 i
i=1 
n  
∂E b 1
P
∂b = −2 yi − a − x2 x2
=0
i i
i=1

The normal equations are given by

n n

1/x2i
P P
yi = na + b
i=1 i=1
n
n
n

yi /x2i 1/x2i + b 1/x4i
P P P
=a
i=1 i=1 i=1
V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 18 / 84
Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Interpolation : Introduction

Interpolation is the technique of estimating the value of a function for any

intermediate value of the independent variable.

The process of computing or finding the value of a function for any value of
the independent variable outside the given range is called extrapolation.

Interpolation is the method of computing the value of the function y = f (x)

for any given value of the independent variable x when a set of values of
y = f (x) for certain values of x are known or given.

If (xi , yi ), i = 0, 1, 2, . . . , n are the set of (n + 1) given data points of the

function y = f (x), then the process of finding the value of y corresponding to
any value of x 6= xi between x0 and xn , is called interpolation.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 19 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Introduction contd...:
Given the set of tabular values (x0 , y0 ) , (x1 , y1 ) , (x2 , y2 ) , . . . , (xn , yn )
satisfying the relation y = f (x) where the nature of f (x) is not known, it
is required to find a simpler function φ (x), such that

f (x) and φ (x) agree at the

set of tabulated points. y Mn
f(x): unknown function
A polynomial of degree n passes
φ n(x) : Approximate function Mi
through (n + 1) distinct points.
The justification for )
M1 y=f(x
approximation of an unknown M2
function y = f (x) by means M0
of a polynomial φ (x) is due φ n(x)
to Weierstrass’s, known as
”Fundamental theorem of y0 y1 y2 yi yn
numerical analysis” or
”Weierstrass’s theorem” a=x0 x1 x2 xi b=xn x

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 20 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Weierstrass’s theorem

Theorem
If f (x) is a continuous in x ∈ [a, b], then for a given ε > 0, there exists a
polynomial φn (x) of degree n, such that

|f (x) − φn (x)| < ε for all x ∈ [a, b]

This means that it is possible to find a polynomial φn (x) whose graph remains
within the region bounded by

y = f (x) − ε and y = f (x) + ε, for all x ∈ [a, b] ,

however, small ε may be.

If f (x) is a continuous in x ∈ [a, b], then we can approximate the function

f (x) by a polynomial φ (x) of sufficiently high degree.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 21 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Introduction contd...

If the function f (x) is known explicitly, then the value of y corresponding to

any value of x can easily be obtained.

On the other hand, if the function f (x) is not known, then it is very hard to
find the exact form of f (x) with the tabulated values (xi , yi ). In such cases,
the function f (x) can be replaced by a simpler, function, say, φ (x), which
has the same values as f (x) for x0 , . . . , xn .

The function φ (x) is called the interpolating or smoothing function and any
other value can be computed from φ (x).

If φ (x) is a polynomial, then φ (x) is called the interpolating polynomial and

the process of computing the intermediate values of y = f (x) is called the
polynomial interpolation.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 22 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Introduction contd...

There are several definitions available for the term interpolation.

Hiral defines interpolation as the estimation of a most likely
estimate in the given conditions. It is the technique of estimating
a past figure.
Theile’s definition of interpolation is ”Interpolation is the art of
reading between the lines of a table”.
Harper’s definition is ”Interpolation consists in reading a value
which lies between two extreme points”.
In the study of interpolation, we make the following assumptions:
there are no sudden jumps in the values of the dependent variable for the
period under consideration
the rate of change of figures from one period to another is uniform.
Here, we study the technique of interpolation based on the calculus of finite
differences. The following important interpolation formulae obtained or derived
based on forward, backward and central differences of a function are presented.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 23 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Interpolation Formulae:

Equal intervals
1 Newton’s binomial expansion formula
2 Newton’s forward interpolation formula
3 Newton’s backward interpolation formula
4 Gauss’s forward interpolation formula
5 Gauss’s backward interpolation formula
6 Bessel’s formula
7 Stirling’s formula
8 Laplace-Everett’s formula
Unequal intervals
1 Lagrange’s formula
2 Newton’s divided difference formula

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 24 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Finite Difference Operators:

Consider a function y = f (x) defined on (a, b). x and y are the independent
and dependent variables respectively.

If the points x0 , x1 , . . . , xn are taken at equidistance i.e., xi = x0 + ih,

i = 0, 1, 2, . . . , n, then the value of y, when x = xi , is denoted as yi , where
yi = f (xi ).

Here, the values of x are called arguments and the values of y are known as
entries. The interval h is called the difference interval.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 25 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Forward Difference Operator ∆:

The forward difference or simply difference operator is denoted by ∆ and may be

defined as
∆ f (x) = f (x + h) − f (x)

The symbol ∆ is called the difference operator. At x = xi ,

∆ f (xi ) = f (xi + h) − f (xi )

writing in terms of y, we have

∆ yi = yi − yi−1 i = 0, 1, 2, . . . , (n − 1)

The differences y1 − y0 , y2 − y1 , . . . , yn − yn−1 are called the first differences of the

function y. They are denoted by ∆y0 , ∆y1 , . . ., etc. That is
∆ y0 = y1 − y0
∆ y1 = y2 − y1
..
.
∆ yi = yi+1 − yi
..
.
∆ yn−1 = yn − yn−1

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 26 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Forward Differences contd...

∆ y0 = y1 − y0
∆ y1 = y2 − y1
..
.
∆ yi = yi+1 − yi
..
.
∆ yn−1 = yn − yn−1
The differences of the first differences are called the second differences and they
are denoted by ∆2 y0 , ∆2 y1 , . . . , ∆2 yn−2 . Hence

∆2 y0 = ∆ (y1 − y0 ) = ∆y1 − ∆y0

= (y2 − y1 ) − (y1 − y0 ) .
= y2 − 2y1 + y0

∆2 y1 = ∆ (y2 − y1 ) = ∆y2 − ∆y1

= (y3 − y2 ) − (y2 − y1 )
= y3 − 2y2 + y1

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 27 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Forward Differences contd...:

∆3 y0 = ∆2 (∆y0 ) = ∆2 (y1 − y0 )
= ∆2 y1 − ∆2 y0
= (y3 − 2y2 + y1 ) − (y2 − 2y1 + y0 )
= y3 − 3y2 + 3y1 − y0
∆3 y1 = y4 − 3y3 + 3y2 − y1 , etc.
Generalising, we have

∆n+1 f (x) = ∆n (∆ f (x)) ,

i.e., ∆n+1 yi = ∆n (∆yi ) , n = 0, 1, . . . Also,

∆n+1 f (x) = ∆n (∆f (x))

= ∆n [f (x + h) − f (x)]
= ∆n f (x + h) − ∆n f (x)

and ∆n+1 yi = ∆n yi+1 −∆n yi , n = 0, 1, 2, . . . where ∆0 is the identity operator

i.e., ∆0 f (x) = f (x) and ∆1 = ∆.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 28 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Forward difference table:

x y ∆y ∆2 y ∆3 y ∆4 y ∆5 y ∆6 y
x0 y0
∆y0
x1 y1 ∆2 y0
∆y1 ∆3 y0
2
x2 y2 ∆ y1 ∆4 y 0
∆y2 ∆3 y1 ∆5 y 0
x3 y3 ∆2 y2 ∆4 y 1 ∆6 y 0
3 5
∆y3 ∆ y2 ∆ y1
x4 y4 ∆2 y3 ∆4 y 2
∆y4 ∆3 y3
2
x5 y5 ∆ y4
∆y5
x6 y6

The forward differences for the arguments x0 , x1 , . . . , xn are shown in Table 1. Table 1
is called a diagonal difference table or forward difference table. The first term in Table
1 is y0 and is called the leading term. The differences ∆y0 , ∆2 y0 , ∆3 y0 , . . . ∆y0 , are
called the leading differences. Similarly, the differences with fixed subscript are called
forward differences.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 29 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Backward Difference Operator ∇:

The backward difference operator is denoted by ∇ and it is defined as

∇f (x) = f (x) − f (x − h)

It can can be written as

∇yi = yi − yi−1 ; i = n, n − 1, . . . ., 1.

or, ∇y1 = y1 − y0 , ∇y2 = y2 − y1 , . . . ., ∇yn = yn − yn−1

The differences are called first differences. The second differences are denoted
by
∇2 y2 , ∇2 y3 , . . . ., ∇2 yn .
Hence
∇2 y2 = ∇ (∇y2 )
= ∇ (y2 − y1 )
= ∇y2 − ∇y1
= (y2 − y1 ) − (y1 − y0 )
= y2 − 2y1 + y0 .
V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 30 / 84
Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Backward Differences contd...

Similarly,
∇2 y3 = y3 − 2y2 + y1 ,
∇2 y4 = y4 − 2y3 + y2 ,
Generalising, we have

∇k yi = ∇k−1 yi − ∇k−1 yi−1 ; i = n, n − 1, . . . ., k

where
∇0 yi = yi , ∇1 yi = ∇yi .
The backward differences written in a tabular form is shown in Table 2. In
Table 2, the differences ∇n y with a fixed subscript 0 i0 lie along the diagonal
upward sloping. Here h is the difference interval.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 31 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Backward difference table:

x y ∇y ∇2 y ∇3 y ∇4 y ∇5 y ∇6 y
x0 y0
∇y1
x1 y1 ∇2 y 2
∇y2 ∇3 y 3
2
x2 y2 ∇ y3 ∇4 y 4
3
∇y3 ∇ y4 ∇5 y 5
2 4
x3 y3 ∇ y4 ∇ y5 ∇6 y6
3 5
∇y4 ∇ y5 ∇ y6
x4 y4 ∇2 y 5 ∇4 y 6
∇y5 ∇3 y 6
2
x5 y5 ∇ y6
∇y6
x6 y6

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 32 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Shift Operator E:

The shift operator is defined as

E f (x) = f (x + h)

or E yi = yi+1

Hence, shift operator sifts the function value yi to the next higher value yi+1 .
The second shift operator gives

E 2 f (x) = E [Ef (x)] = E [f (x + h)] = f (x + 2h)

E is linear and obeys the law of indices. Generalising,

E n f (x) = f (x + nh) or E n yi = yi+nh

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 33 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Inverse Shift Operator E −1 :

The inverse shift operator E −1 is defined as

E −1 f (x) = f (x − h)

In a similar manner, second and higher inverse operators are given by

E −2 f (x) = f (x − 2h)

and
E −n f (x) = f (x − nh)
The more general form of E operator is given by

E r f (x) = f (x + rh)

where r is positive as well as negative rationals.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 34 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Differential Operator D:

The differential operator is usually denoted by D, where

d d2
D f (x) = f (x) = f 0 (x) and D2 f (x) = f (x) = f 00 (x)
dx dx2
For linking different operators with differential operator D we consider Taylor’s
formula:
h2
f (x + h) = f (x) + h f 0 (x) + f 00 (x) + . . .
2!
In operator notation, we can write it as:

h2
 
E f (x) = 1 + h D + D2 + . . . f (x) = ehD f (x)
2!

This series in brackets is the expression for the exponential and hence we can
write
E = ehD

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 35 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Properties of the Operator ∆:

1 If c is a constant then ∆ c = 0.
Let f (x) = c, then f (x + h) = c, where h is the interval of differencing.
Hence
∆f (x) = f (x + h) − f (x) = c − c = 0
or ∆c = 0
2 ∆ is distributive, i.e., ∆ [f (x) ± g (x)] = ∆f (x) ± ∆g (x) .

∆ [f (x) + g (x)] = [f (x + h) + g (x + h)] − [f (x) + g (x)]

= f (x + h) − f (x) + g (x + h) − g (x)
= ∆f (x) + ∆g (x) .
Similarly, we have ∆ [f (x) − g (x)] = ∆f (x) − ∆g (x) .
3 If c is a constant then ∆ [cf (x)] = c∆f (x) . From properties 2 and 3
above, it is observed that ∆ is a linear operator.

∆ [c f (x)] = c f (x + h) − c f (x) = c [f (x + h) − f (x)] = c∆f (x)

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 36 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Properties of the Operator ∆ contd...:

4 If m and n are positive integers then ∆m ∆n f (x) = ∆m+n f (x) .

∆m ∆n f (x) = (∆ × ∆ × ∆...m times)(∆ × ∆ . . . n times)f (x)

= (∆∆∆ . . . (m + n) times) f (x)
= ∆m+n f (x) .

In a similar manner, we can prove the following properties:

5 ∆[f1 (x) + f2 (x) + . . . + fn (x)] = ∆f1 (x) + ∆f2 (x) + . . . + ∆fn (x) .
6 ∆ [f (x) g (x)] = f (x) ∆g (x) + g (x) ∆f (x) .

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 37 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Properties of the Operator ∆ contd...:

h i
f (x) g(x)∆f (x)−f (x)∆g(x)
7 ∆ g(x) = g(x)g(x+h) .
 
f (x)
∆
g (x)
f (x + h) f (x)
= −
g (x + h) g (x)
g (x) f (x + h) − f (x) g (x + h)
=
g (x) g (x + h)
g (x) f (x + h) − g (x) f (x) − f (x) g (x + h) + f (x) g (x)
=
g (x) g (x + h)
g (x) [f (x + h) − f (x)] − f (x) [g (x + h) − g (x)]
=
g (x) g (x + h)
g (x) ∆f (x) − f (x) ∆g (x)
= .
g (x) g (x + h)

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 38 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Relation between ∆, E and ∇:

From the definition of ∆, we know that

∆f (x) = f (x + h) − f (x)

where h is the interval of differencing. Using the operator E we can write

∆ f (x) = E f (x) − f (x)

∆ f (x) = (E − 1) f (x)
The above relation can be expressed as an identity

∆ =E−1 or E =1+∆
∇f (x) = f (x) − f (x − h) = f (x) − E −1 f (x)
∇ = 1 − E −1

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 39 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Relation between ∆, E and ∇:

1 (1 + ∆) (1 − ∇) = 1
(1 + ∆) (1 − ∇) f (x) = E E −1 f (x) = E 0 f (x) = 1
2 E ∆ = ∆E.
E ∆ f (x) = E (f (x + h) − f (x))
= E f (x + h) − E f (x)
= f (x + 2h) − f (x + h)
= ∆ f (x + h)
= ∆E f (x)
Hence, E ∆ = ∆E.
3 ∆ = E∇ = ∇E
E ∇ [f (x)] = E [f (x) − f (x − h)]
= E [f (x)] − E [f (x − h)]
= f (x + h) − f (x)
= ∆f (x)

Hence E∇=∆
V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 40 / 84
Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Relation between ∆, E and ∇:

Again, ∇E [f (x)] = ∇f (x + h) = f (x + h) − f (x) = ∆f (x)

Hence, ∇E = ∆

4 ∆∇ = ∇∆ 6 ∇ = E −1 ∆
E −1 ∆f (x)
∇∆ f (x) = E −1 [f (x + h) − f (x)]
= ∇ [f (x + h) − f (x)] = E −1 f (x + h) − E −1 f (x)
= ∇f (x + h) − ∇f (x) = f (x) − f (x − h)
= f (x + h) − f (x) − [f (x) − f (x − h)] = ∇f (x)
= ∆ [f (x) − f (x − h)]
= ∆∇f (x)

5 ∆∇ = ∆ − ∇

∆∇f (x)
= ∆ [f (x) − f (x − h)]
= ∆f (x) − ∆f (x − h)
= ∆f (x) − ∇f (x)
= (∆ − ∇) f (x)
V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 41 / 84
Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Example:
∆2 x E ex
Prove that ex = E e . ∆2 e x , the interval of differencing being h.

Solution: We know that Ef (x) = f (x + h). Hence

E ex = ex+h ,

Again

∆ ex = ex+h − ex = ex (eh − 1) ⇒ ∆2 ex = ex .(eh − 1)2

Hence
∆2
 
ex = (∆2 E −1 )ex = ∆2 ex−h = e−h (∆2 ex ) = e−h ex (eh − 1)2
E

Therefore, the right hand side

ex+h
= e−h ex (eh − 1). = ex
ex (eh − 1)
V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 42 / 84
Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Example:

 
∆2
Evaluate E x3

Solution: Let h be the interval of differencing

 2
∆
E x3
= ∆2 E −1 x3


2
= (E − 1) E −1 x
3
= E − 2E + 1
E −1 x3
2

= E − 2 + E −1 x3
= Ex3 − 2x3 + E −1 x3
3 3
= (x + h) − 2x3 + (x − h)
= 6xh
 
∆2
Note: If h = 1, then E x3 = 6x

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 43 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Introduction Operators Difference Table

Example:
 
∆f (x)
Show that ∆ logf (x) = log 1 + f (x)

Solution: Let h be the interval of differencing

f (x + h) = Ef (x) = (∆ + 1) f (x) = ∆f (x) + f (x)

f (x + h) ∆f (x)
= +1
f (x) f (x)
Taking logarithms on both sides we get
   
f (x + h) ∆f (x)
log = log 1 +
f (x) f (x)
 
∆f (x)
log f (x + h) − log f (x) = log 1 +
f (x)
 
∆f (x)
∆ log f (x) = log 1 +
f (x)
V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 44 / 84
Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Interpolation with Equal Intervals:

Missing Values
Let a function y = f (x) is given for equally spaced values x0 , x1 , . . . , xn of
the argument and y0 , y1 , y2 , . . . , yn denote the corresponding values of the
function. If one or more values of y = f (x) are missing, we can determine the
missing values by employing the relationship between the operators E and ∆.

Newton’s Binomial Expansion Formula

Suppose y0 , y1 , y2 , . . . , yn denote the values of the function y = f (x)
corresponding to the values

x0 , x0 + h, x0 + 2h, . . . ., x0 + nh

of x. Let one of the values of y is missing since n values of the functions are
known. Therefore, we have
n
∆n y0 = 0 or (E − 1) y0 = 0

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 45 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Missing values contd...

or
n
(E − 1) y0 = 0
Expanding, we have
 n n n
E − C1 E n−1 + n C2 E n−2 + . . . + (−1) y0 = 0

n (n − 1) n−2 n
or E n y0 − n E n−1 y0 + E y0 + . . . + (−1) y0 = 0
2
n (n − 1) n
or yn − n yn−1 + yn−2 + . . . + (−1) y0 = 0
2
It is quite useful in determining the missing values without actually
constructing the difference table.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 46 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Example : One Missing value

Example Find the missing entry in the following table.

x 0 1 2 3 4 5
y = f (x) 1 3 11 ? 189 491
Solution: Here, we are given y0 = 1, y1 = 3, y2 = 11, y4 = 189 and y5 =
491. Since five values are given, we assume that y is a polynomial of degree
4. Hence
∆5 y 0 = 0
5
(E − 1) y0 = 0
(1)
E − 5E 4 + 10E 3 − 10E 2 + 5E − 1 y0 = 0
5

y5 − 5y4 + 10y3 − 10y2 + 5y1 − y0 = 0

Substituting the given values fory0 , y1 , . . . ., y5 in Eq.(1), we get

491 − 5 (189) + 10y3 − 10 (11) + 5 (3) − 1 = 0

10y3 = 550 or y3 = 55

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 47 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Example : Two Missing values

Example: Find the missing entries in the following table.

x 0 1 2 3 4 5
y = f (x) 1 ? 11 28 ? 116
Solution: Here, we are given y0 = 1, y2 = 11, y3 = 28, and y5 = 116.
Since three values are known, we assume y = f (x)as a polynomial of
degree three.

∆4 y 0 = 0
4
(E − 1) y0 = 0

E 4 − 4E 3 + 6E 2 − 4E + 1 y0 = 0
(1)
y4 − 4y3 + 6y2 − 4y1 + y0 = 0
y 4 − 4 (28) + 6 (11) − 4y1 + 1 = 0
y4 − 4y1 = 45

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 48 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Example : Missing values contd...

and
∆5 y 0 = 0
5
(E − 1) y0 = 0

E 5 − 5E 4 + 10E 3 − 10E 2 + 5E − 1 y0 = 0
(2)
y5 − 5y4 + 10y3 − 10y2 + 5y1 − y0 = 0
116 − 5y4 + 10 (28) − 10 (11) + 5y1 − 1 = 0
−5y4 + 5y1 = −285
Solving Eqs.(1) and (2), we obtain

y1 = 4 and y4 = 61.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 49 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Newton’s Forward Interpolation Formula:

Let y = f (x), which takes the values y0 , y1 , y2 , . . . , yn , that is the set of

(n + 1) functional values y0 , y1 , y2 , . . . , yn are given corresponding to the set
of (n + 1) equally spaced values of the independent variable, xi = x0 + i h,
i = 0, 1, 2, . . . , n where h is the spacing.

Let φ (x) be a polynomial of the nth degree in x taking the same values as
y corresponding to x = x0 , x1 , . . . , xn . Then, φ (x) represents the continuous
function y = f (x) such that

f (xi ) = φ (xi ) ; for i = 0, 1, 2, . . . , n

f (x) = φ (x) + R (x) ; at all other points

where R (x) is called the error term (remainder term) of the interpolation
formula.
Interpolation function: a function that passes exactly through a set
of data points.
V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 50 / 84
Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Newton’s Forward contd...

Let
φ (x) = a0 + a1 (x − x0 ) + a2 (x − x0 ) (x − x1 )
+a3 (x − x0 ) (x − x1 ) (x − x2 ) (1)
+ . . . + an (x − x0 ) (x − x1 ) . . . (x − xn−1 )
and φ (xi ) = yi ; i = 0, 1, 2, . . . , n (2)
The constants a0 , a1 , a2 , . . . , an can be determined as follows:
Substituting x = x0 , x1 , . . . , xn successively in Eq.(1) and using Eq.(2), we
get
a 0 = y0
y1 = a0 + a1 (x1 − x0 )
or, y1 = y0 + a1 (x1 − x0 )
or, a1 h = y1 − y0 = ∆y0
∆y0
a1 =
h

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 51 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Newton’s Forward contd...

y2 = a0 + a1 (x2 − x0 ) + a2 (x2 − x0 ) (x2 − x1 )

∆y0
= y0 + 2h + a2 . 2h. h
h
2h2 a2 = y2 − y0 − 2 (y1 − y0 )
= y2 − y0 − 2y1 + 2y0
= y2 − 2y1 + y0
2
= (E − 1) y0
∆2 y0
a2 =
2!h2
Similarly, we obtain
∆n y 0
an =
n!hn
Hence, from Eq.(1), we have

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 52 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Newton’s Forward contd...

∆y0 ∆2 y 0
(x − x0 ) +
φ (x) = y0 + (x − x0 ) (x − x1 )
h 2!h2
∆3 y 0 (3)
+ (x − x0 ) (x − x1 ) (x − x2 ) + . . .
3!h3
n
∆ y0
... + (x − x0 ) (x − x1 ) . . . (x − xn−1 )
n!hn
Let x = x0 + u h, we have

x − x0 = u h
x − x1 = (x − x0 ) − (x1 − x0 ) = u h − h = (u − 1) h and (4)
x − x2 = (x − x0 ) − (x2 − x0 ) = uh − 2h = (u − 2) h, etc.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 53 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Newton’s Forward contd...

Using the values from Eq.(4), Eq.(3) reduces to

u (u − 1) 2
φ (x) = y0 + u∆y0 + ∆ y0
2!
u (u − 1) (u − 2) 3
+ ∆ y0 (5)
3!
u (u − 1) . . . (u − n + 1) n
+... + ∆ y0
n!
φ (x) = y0 + u C1 ∆y0 + u C2 ∆2 y0 + . . . + u Cn ∆n y0
where u = (x − x0 ) /h

The formula given in Eq.(5) is called the Newton’s forward interpolation

formula. This formula is used to interpolate the values of y near the beginning
of a set of equally spaced tabular values. This formula can also be used for
extrapolating the values of y a little backward of y0 .

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 54 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Newton’s Backward Interpolation Formula:

Newton’s forward interpolation formula is not suitable for interpolation values

of y near the end of a table of values.

Let y = f (x) be a function which takes the values y0 , y1 , y2 , . . . , yn ,

corresponding to the values x0 , x1 , x2 , . . . , xn of the independent variable x.

Let the values of x be equally spaced with h as the interval of differencing.

That is xi = x0 + i h, i = 0, 1, 2, . . . , n

Let φ (x) be a polynomial of the nth degree in x taking the same values of y
corresponding to x = x0 , x1 , . . . , xn . That is, φ (x) represents y = f (x) such
that f (xi ) = φ (xi ), i = 0, 1, 2, . . . , n.

Hence, we can write φ (x) as φ (xi ) = yi , i = n, n − 1, . . . , 1, 0 and xn−i =

xn − i h, i = 1, 2, . . . , n.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 55 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Newton’s backward contd...

Let
φ (x) = a0 + a1 (x − xn ) + a2 (x − xn ) (x − xn−1 )
+a3 (x − xn ) (x − xn−1 ) (x − xn−2 ) (1)
+ . . . + an (x − xn ) (x − xn−1 ) . . . (x − x1 )
The constants a0 , a1 , a2 , . . . , an can be determined as follows:
Substituting x = xn , xn−1 , . . . , x1 successively in Eq.(1), we get

a0 = yn

yn−1 = a0 + a1 (xn−1 − xn ) or yn−1 = yn + a1 (xn−1 − xn )

or, a1 h = yn − yn−1 = ∆yn−1
∆yn−1
a1 =
h

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 56 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Newton’s backward contd....

Similarly, we obtain

yn−2 = a0 + a1 (xn−2 − xn ) + a2 (xn−2 − xn ) (xn−2 − xn−1 )

1
= yn + ∆yn−1 (−2h) + a2 (−2h) (−h)
h
= yn − 2∆yn−1 + 2h2 a2

2h2 a2 = yn−2 + 2∆yn−1 − yn

= yn−2 + 2 (yn − yn−1 ) − yn
2
= yn − 2yn−1 + yn−2 = (E − 1) yn−2
= ∆2 yn−2
∆2 yn−2
or, a2 =
2h2

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 57 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Newton’s backward contd...

Similarly, we obtain
∆n y 0
an =
n! hn
Substituting the values from a0 , a1 , a2 , . . . , an in Eq.(1), we get

∆yn−1 ∆2 yn−2
φ (x) = yn + (x − xn ) + (x − xn ) (x − xn−1 )
h 2!h2
3
∆ yn−3 (2)
+ (x − xn ) (x − xn−1 ) (x − xn−2 )
3!h3
n
∆ y0
+... + (x − xn ) (x − xn−1 ) . . . (x − x1 )
n!hn
Let x = xn + u h

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 58 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Newton’s backward contd...

x − xn = u h
x − xn−1 = (x − xn ) + (xn − xn−1 ) = u h + h = (u + 1) h (3)
x − xn−2 = (x − xn ) − (xn − xn−2 ) = uh + 2h = (u + 2) h, etc.

Using the values from Eq.(3), Eq.(2) reduces to

u (u + 1) 2 u (u + 1) (u + 2) 3
φ (x) = yn + u∆yn−1 + ∆ yn−2 + ∆ yn−3
2! 3! (4)
u (u + 1) . . . (u + n − 1) n
+... + ∆ y0
n!
where u = x−x h
n
The formula given in Eq.(4) is called the Newton’s backward
interpolation formula. This formula is used for interpolating values of y near
the end of the tabulated values and also used for extrapolating values of y a
little backward of yn .

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 59 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Example:
If y = f (x) is known at the following Newton forward formula, with u = (x − x0 )/h, is
data points: given by

x 0 1 2 3 4 u(u−1)
y (x) = y0 + u∆y0 + 2!
∆2 y0
y 1 7 23 55 109 u(u−1)(u−2) 3
+ 3!
∆ y0
Find y(0.5) using Newton’s
forward difference formula Here x = 0.5, x0 = 0, h = 1 u = 0.5, thus
y(3.5) using Newton’s backward
difference formula 0.5(0.5−1)
y (0.5) = 1 + 0.5 × 6 + 2!
10
Difference Table 0.5(0.5−1)(0.5−2)
+ 3!
6
= 1 + 3 + 2.5 × (−0.5) + (−0.25) (−1.5)
x y ∆y ∆y 2 ∆y 3 ∆y 4 = 3.125
0 1
6 Newton backward formula, with u = (x − xn )/h, is
1 7 10 given by
16 6
2 23 16 0 u(u+1)
32 6 y (x) = yn + u∆yn−1 + 2!
∆2 yn−2
u(u+1)(u+2) 3
3 55 22 + 3!
∆ yn−3
54
4 109 Here x = 3.5, xn = 4.0, h = 1, u = −0.5,
yn = 109, ∆yn−1 = 54, ∆2 yn−2 = 22, ∆3 yn−3 = 6. Thus, we have y(3.5) =?.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 60 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Lagrange’s Formula for Unequal Intervals:

Let y = f (x) be a real valued continuous function defined in an interval [a, b].

Let x0 , x1 , x2 , . . . , xn be (n + 1) distinct points which are not necessarily

equally spaced and the corresponding values of the function are
y0 , y1 , y2 , . . . , yn .

Since (n + 1) values of the function are given corresponding to the (n + 1)

values of the independent variable x, we can represent the function y = f (x)
is a polynomial in x of degree n.
Let the polynomial is represented by

f (x) = a0 (x − x1 ) (x − x2 ) . . . (x − xn )
+a1 (x − x0 ) (x − x2 ) . . . (x − xn ) + . . .
(1)
+ai (x − x0 ) ... (x − xi−1 ) (x − xi+1 ) . . . (x − xn−1 ) + . . .
. . . + an (x − x0 ) (x − x1 ) . . . (x − xn−1 )

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 61 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Lagrange’s contd...

Each term in Eq.(1) being a product of n factors in x of degree n, putting

x = x0 in Eq.(1) we obtain

f (x0 ) = a0 (x0 − x1 ) (x0 − x2 ) . . . (x0 − xn )

f (x0 )
or, a0 =
(x0 − x1 ) (x0 − x2 ) . . . (x0 − xn )
Putting x = x1 in Eq.(1) we obtain

f (x1 )
a1 =
(x1 − x0 ) (x1 − x2 ) . . . (x1 − xn )

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 62 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Lagrange’s contd...

Similarly putting x = xi in Eq.(1) we obtain

f (xi )
ai =
(xi − x0 ) .... (xi − xi−1 ) (xi − xi+1 ) . . . (xi − xn )

and
f (xn )
an =
(xn − x0 ) (xn − x1 ) . . . (xn − xn−1 )

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 63 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Lagrange contd...

Substituting the values of a0 , a1 , ..., an in Eq.(1) we get

(x − x1 ) (x − x2 ) . . . (x − xn )
f (x) = f (x0 )
(x0 − x1 ) (x0 − x2 ) . . . (x0 − xn )

(x − x0 ) (x − x2 ) . . . (x − xn )
+ f (x1 ) + . . .
(x1 − x0 ) (x1 − x2 ) . . . (x1 − xn )

(x − x0 ) ... (x − xi−1 ) (x − xi+1 ) . . . (x − xn−1 )

+ f (xi ) + . . .
(xi − x0 ) ... (xi − xi−1 ) (xi − xi+1 ) . . . (xi − xn−1 )

(x − x0 ) (x − x1 ) . . . (x − xn−1 )
.... + f (xn )
(xn − x0 ) (xn − x1 ) . . . (xn − xn−1 )
(2)
The formula given by Eq.(2) is known as the Lagrange’s interpolation formula.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 64 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Divided Differences

Let the function y = f (x) be given at the point x0 , x1 , x2 , . . . , xn (which need

not be equally spaced) f (x0 ) , f (x1 ) , f (x2 ) , . . . , f (xn ) , denote the (n + 1)
values the function at the points x0 , x1 , x2 , . . . , xn .

Then the first divided differences of f (x) for the arguments x0 , x1 is defined
as
f (x1 ) − f (x0 )
f (x0 , x1 ) =
x1 − x0
It is denoted by f (x0 , x1 )or by [x0 , x1 ]
Similarly, the second divided difference for the arguments x0 , x1 , x2 , is given by

f (x1 , x2 ) − f (x0 , x1 )
f (x0 , x1 , x2 ) =
x2 − x0

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 65 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Divided Differences contd...

In general, the nth divided difference is given by

f (x1 , x2 , . . . , xn ) − f (x0 , x1 , . . . , xn−1 )

f (x0 , x1 , . . . , xn ) =
xn − x0
The first divided differences are called the divided differences of order one, the
second divided differences are called of order two and so on.

Properties of Divided Difference:

(i) The divided differences are symmetric with respect to the arguments.,

f (x1 ) − f (x0 ) − [f (x0 ) − f (x1 )]

f (x0 , x1 ) = = = f (x1 , x0 )
x1 − x0 − (x0 − x1 )

Also, f (x0 , x1 , x2 ) = f (x2 , x1 , x0 )

(ii) The nth divided difference of a polynomial of degree n is constant.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 66 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Newton’s Divided Difference Interpolation Formula:

Let the function y = f (x) be given at the point x0 , x1 , x2 , . . . , xn which are

not equally spaced. Suppose that f (x0 ) , f (x1 ) , f (x2 ) , . . . , f (xn ), denote
the (n + 1) values the function at the points x0 , x1 , x2 , . . . , xn
We know that the first divided differences of f (x) is given by

f (x0 ) − f (x)
f (x, x0 ) =
(x0 − x)

f (x) = f (x0 ) + (x − x0 ) f (x, x0 )

Also, the the second divided difference is given by

f (x0 , x1 ) − f (x, x0 )
f (x, x0 , x1 ) =
(x1 − x)
f (x, x0 ) = f (x0 , x1 ) + (x − x1 ) f (x, x0 , x1 )

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 67 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Newton’s divided contd...

f (x) − f (x0 )
= f (x0 , x1 ) + (x − x1 ) f (x, x0 , x1 )
(x − x0 )
f (x) = f (x0 ) + (x − x0 ) f (x0 , x1 ) + (x − x0 ) (x − x1 ) f (x, x0 , x1 )

f (x) = f (x0 ) + (x − x0 ) f (x0 , x1 ) + (x − x0 ) (x − x1 ) f (x, x0 , x1 )

+ (x − x0 ) (x − x1 ) (x − x2 ) f (x, x0 , x1 , x2 )

Proceeding in the same way, we get

f (x) = f (x0 ) + (x − x0 ) f (x0 , x1 ) + (x − x0 ) (x − x1 ) f (x0 , x1 , x2 )

+ (x − x0 ) (x − x1 ) (x − x2 ) f (x0 , x1 , x2 , x3 ) + · · ·
· · · + (x − x0 ) (x − x1 ) . . . (x − xn−1 ) f (x0 , x1 , . . . , xn ) + · · ·
· · · + (x − x0 ) (x − x1 ) . . . (x − xn ) f (x, x0 , x1 , . . . , xn )

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 68 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Missing Newton Forward Backward Lagrange Divided

Newton’s divided contd...

Since the function f (x) is a polynomial of degree n, therefore

f (x, x0 , x1 , . . . , xn ) = 0

Hence we may write f (x) in terms of divided differences as follows:

f (x) = f (x0 ) + (x − x0 ) f (x0 , x1 ) + (x − x0 ) (x − x1 ) f (x0 , x1 , x2 )

+ (x − x0 ) (x − x1 ) (x − x2 ) f (x0 , x1 , x2 , x3 )
+ (x − x0 ) (x − x1 ) (x − x2 ) (x − x3 ) f (x0 , x1 , x2 , x3 , x4 ) + · · ·

This is known as Newton’s divided difference formula. Newton’s divided

difference formula reduces to Newton’s forward formula if the values of the
arguments are equally spaced.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 69 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Numerical Integration:

If F (x) is a differentiable function whose derivative is f (x), then we can

evaluate the definite integral I as

Zb
I= f (x) dx = F (b) − F (a) ; F 0 (x) = f (x) (1)
a

Eq. (1) is known as the fundamental theorem of calculus. In many cases a

mathematical expression for f (x) is unknown and in some cases even if f (x) is
known its complex form makes it difficult to perform the integration.

Zb
I= f (x)dx = area under the curve f (x) between x = a to x = b.
a

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 70 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Quadrature contd...

y
f(b)
f(x)

f(a)

x0 = a x0 = b x
V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 71 / 84
Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Newton-Cotes closed Quadrature Formula

The general form of the problem of numerical integration may be stated as

follows: Given a set of data points (xi , yi ), i = 0, 1, . . . , n of a function
y = f (x), where f (x) is not explicitly known. It is required to evaluate the
definite integral
Zb Zb
I = f (x) dx = y dx (2)
a a

We replace y = f (x) by an interpolating polynomial φ (x) in order to obtain

an approximate value of the definite integral of Eq.(2).

We derive a general formula for numerical integration by using Newton’s forward

interpolation formula. The interval (a, b) is divided into n-equal subintervals
such that h = (b − a) /n

a = x0 < x1 < x2 < x3 . . . . < xn = b

with xn = x0 + nh, where

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 72 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Newton-Cotes contd...

h is the interval size

n is the number of subintervals
a and b are the limits of integration with b > a.
Hence, the integral in Eq.(2) can be written as

Zxn Zb Zb
I= y dx = f (x) dx ∼
= φ (x) dx
x0 a a

Using Newton’s forward interpolation formula, we have

Zxn 
u (u − 1) 2 u (u − 1) (u − 2) 3
I= y0 + u ∆y0 + ∆ y0 + ∆ y0 + ... dx
2! 3!
x0

where x = x0 + u h, so dx = h du

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 73 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Newton-Cotes contd...

Zn "

#
u2 − u 2 u3 − 3u2 + 2u 3
I=h y0 + u ∆y0 + ∆ y0 + ∆ y0 + ... du
2! 3!
0
n
1 u3 u2 1 u4 3u3 2u2
    
2 3
= h y0 + u ∆y0 + − ∆ y0 + − + ∆ y0 + ...
2 3 2 6 4 3 2 0
Hence, after simplification, we get
Zxn " #
n n (2n − 3) 2 n(n − 2)2 3
I= y dx = n h y0 + ∆y0 + ∆ y0 + ∆ y0 + ... (3)
2 12 24
x0

The formula given by Eq.(3) is known as Newton-Cotes closed quadrature

formula. From the general formula (Eq.(3)), we can derive or deduce different
integration formulae by substituting n = 1, 2, 3, . . . etc.
n = 1 : Trapezoidal rule ; n = 2 : Simpson’s 1/3 rule
n = 3 : Simpson’s 3/8 rule ; n = 6 : Weddle’s rule.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 74 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Trapezoidal Rule:

Substituting n = 1 in Eq.(3) and considering the curve y = f (x) through the

points (x0 , y0 ) and (x1 , y1 ) as a straight line (a polynomial of first degree) so
that the differences of order higher than first become zero), we get
Zx1  
1
I1 = y dx = h y0 + ∆y0
2
x0
 
1 h
= h y0 + (y1 − y0 ) = (y0 + y1 )
2 2

Similarly, we have
Zx2
h
I2 = y dx = (y1 + y2 )
2
x1

Zx3
h
I3 = y dx = (y2 + y3 )
2
x2

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 75 / 84

 Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Trapezoidal contd...
In this method, the known function
values are joined by straight lines. The
area enclosed by these lines between
the given end points is computed to
approximate the integral as shown in Fig. y
y1=f(x1)

Each subinterval with the line y=f(x)

approximation for the function forms
a trapezoid as shown in Fig. The y0=f(x0)
area of each trapezoid is computed by
multiplying the interval size h by the
average value of the function value in
that subinterval.
x0 x1 x
After the individual trapezoidal areas are
obtained, they are all added to obtain the
overall approximation to the integral.
V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 76 / 84
Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Trapezoidal contd...

In general, we have and so on. (see Fig.)

Zxn
h
In = y dx = (yn−1 + yn )
2
xn−1

Adding all the integrals Ii ; i = 1, . . . , n and using the interval additive property
of the definite integrals, we obtain

n Zxn
X h
I= Ii = y dx = [y0 + 2 (y1 + y2 + · · · + yn−1 ) + yn ] (4)
i=1
2
x0

Equation (4) is known as the trapezoidal rule. Summarising, the trapezoidal

rule signifies that the curve y = f (x) is replaced by n-straight lines joining the
points (xn , yn ), i = 0, 1, 2, 3 . . . , n. The area bounded by the curve y = f (x),
the ordinates x = x0 , x = xn and the x-axis is then approximately equivalent
to the sum of the areas of the n-trapezoids so obtained.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 77 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Simpson’s 1/3 Rule:

In Simpson’s rule, the function is approximated by a second degree polynomial

between successive points. Substituting n = 2 in Eq. (3) and taking the curve
through the points (x0 , y0 ) , (x1 , y1 ) and (x2 , y2 ) as a polynomial of second
degree (parabola) so that the differences of order higher than two vanish, we
obtain
Zx2  
1
y dx = 2h y0 + ∆y0 + ∆2 y0
6
x0
 
1
= 2h y0 + (y1 − y0 ) + (y2 − 2y1 + y0 )
6
2h
= (y2 − 2y1 + y0 + 6y1 )
6
h
= (y0 + 4y1 + y2 )
3

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 78 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Simpson’s 1/3 contd...

f(x) yn

Simpson’s 1/3 rule

)
y=f(x
y2
y1
y0

a=x0 x1 x2 xi xi+1 xi+2 b=xn x

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 79 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Simpson’s 1/3 contd...

Similarly,
Zx4
h
y dx = (y2 + 4y3 + y4 )
3
x2
xZn−2
h
y dx = (yn−4 + 4yn−3 + yn−2 )
3
xn−4

In general, we can write

Zxn
h
y dx = (yn−2 + 4yn−1 + yn )
3
xn−2
x
Rn h
In general, we can write y dx = 3 (yn−2 + 4yn−1 + yn )
xn−2

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 80 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Simpson’s 1/3 contd...

Summing up all the above integrals, we obtain

Zxn
h
y dx = [y0 + 2 (y2 + y4 + . . . + yn−2 ) + 4 (y1 + y3 + . . . + yn−1 ) + yn ]
3
x0
(5)
Zxn
h
y dx = (X + 2E + 4O)
3
x0

X : sum of end ordinates

O : sum of odd ordinates
E : sum of even ordinates
Equation (5) is known as Simpson’s 1/3 rule. Simpson’s 1/3 rule requires
the whole range (the given interval) must be divided into even number of equal
subintervals. In this rule, the interpolating polynomial is of degree 2. Therefore,
this rule is also known as parabolic rule.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 81 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Simpson’s 3/8 Rule:

Putting n = 3 into general quadrature formula and taking the curve through
the points in (x0 , y0 ) , (x1 , y1 ) , (x2 , y2 ) and (x3 , y3 ) as a polynomial of degree
three so that the differences of order greater than become zero. We get
Rx3 h i
y dx = 3h y0 + 32 ∆y0 + 3(3) 2 3 3
12 ∆ y0 + 24 ∆ y0
x0
= 3h y0 + 32 (y1 − y0 ) + 14 (y2 − 2y1 + y0 ) + 1
 
8 (y3 − 3y2 + 3y1 − y0 )
= 3h
8 (y0 + 3y1 + 3y2 + y3 )

Similarly in the intervals [x3 , x6 ] , [x6 , x9 ] , . . . , [xn−3 , xn ], we get

Zx6
3h
y dx = (y3 + 3y4 + 3y5 + y6 )
8
x3

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 82 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Simpson’s 3/8 contd...

f(x) yn

Simpson’s 3/8 rule

y2 y3 y=f(x)
y1
y0

a=x0 x1 x2 x3 xi xi+1 xi+2 xi+3 b=xn x

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 83 / 84

Numerical Analysis Least Squares Interpolation Formulae Integration Newton-Cotes Trapezoidal Simpson 1/3 Simpson 3/8

Simpson’s 3/8 contd...

Rx9 3h
y dx = 8 (y6 + 3y7 + 3y8 + y9 )
x6
··· ··· ··· ···
x
Rn
y dx = 3h
8 (yn−3 + 3yn−2 + 3yn−1 + yn )
xn−3

Summing up all the above integrals, we obtain

x
Rn 3h
y dx = 8 [y0 + 3 (y1 + y2 + y4 + y5 + . . . + yn−1 )
x0 (6)
+2 (y3 + y6 + . . . + yn−3 ) + yn ]

Equation (6) is known as Simpson’s 3/8 rule. It requires the whole range (the
given interval) must be divided into equal subintervals. Here, the number of
subintervals should be taken as multiples of 3. In this rule, the interpolating
polynomial is of degree 3.

V. Kumar, DDU Gorakhpur University B.Sc.-I : Numerical Methods Sept. 2018 84 / 84