Department of Software Engineering

Lakehead University

ESOF 3558: Numerical Methods and Modeling

Due Date: March 17, 2022

Solve the following problems. Each problem will require some analysis, as well as some computations.
Show your work for the computations. You may verify your solution using MATLAB or your choice of
other programming languages.

 Question 1 (Interpolation)(30 points)

Ohm’s law states that the voltage drop V across an ideal resistor is linearly proportional to the current i
flowing through the resistor as in V = iR, where R is the resistance. However, real resistors may not obey
Ohm’s law. Suppose that you performed some very precise experiments to measure the voltage drop and
corresponding current for a resistor. The following results suggest a curvilinear relationship rather than
the straight line represented by Ohm’s law:

i -1 -0.5 -0.25 0.25 0.5 1

a) Use a fifth-order Newton interpolating polynomial to fit the data and compute V for i = 0.10.
b) Use a fifth-order Lagrange interpolating polynomial to fit the data and compute V for i = 0.10.
Compare your solution with part a).
c) Fit a fourth-order Newton interpolating polynomial to the first 5 points. Use it to compute the
voltage at i = 1 by extrapolation and compare your prediction with the actual result shown in the
table.

Solutions to Q1:

General formula for (n-1)th order Newton interpolating polynomial:

f n1 ( x)  b1  b2 ( x  x1 )  ...  bn ( x  x1 )( x  x 2 )...(x  xn1 )

where,
𝑏 𝑓 𝑥
𝑏 𝑓 𝑥 ,𝑥
𝑏 𝑓 𝑥 ,𝑥 ,𝑥
.
.
𝑏 𝑓 𝑥 ,𝑥 ,...,𝑥 ,𝑥
1st finite difference: 𝑓 𝑥 , 𝑥

, ,
2nd finite difference: 𝑓 𝑥 , 𝑥 , 𝑥
…

, ,… , , ,… ,
(n-1)th finite difference: 𝑓 𝑥 , 𝑥 ,…,𝑥 ,𝑥

Solution to 1.a

5th order Newton’s interpolating polynomial:

f 5 ( x)  b1  b2 ( x  x1 )  b3 ( x  x1 )( x  x 2 )  ...  b6 ( x  x1 )( x  x 2 )( x  x3 )( x  x 4 )( x  x5 )

1st finite 2nd finite 3rd finite 4th finite 5th finite
i V(i) difference difference difference difference difference
1 ‐637 1081 ‐1036 592 0 0
 0.5 ‐96.5 304 ‐296 592 0
 0.25 ‐20.5 82 296 592
0.25 20.5 304 1036
0.5 96.5 1081
1 637

The values on the first row are coefficients of the Newton Polynomial:
b1 =  637
b2 = 1081
b3 =  1036
b4 = 592
b5 = 0
b6 = 0

5th order Newton interpolating polynomial:

𝑉 𝑖 𝑏 𝑏 𝑖 𝑖 𝑏 𝑖 𝑖 𝑖 𝑖 ... 𝑏 𝑖 𝑖 𝑖 𝑖 𝑖 𝑖 𝑖 𝑖 𝑖 𝑖

where i1 = 1, i2 =  0.5, i3 =  0.25, i4 = 0.25, i5 = 0.5, i6 = 1.

𝑉 𝑖 637 1081 𝑖 1 1036 𝑖 1 𝑖 0.5 592 𝑖 1 𝑖

Hence:

𝑉 𝑖 637 1081 𝑖 1 1036 𝑖 1 𝑖 0.5 592 𝑖 1 𝑖 0.5 𝑖 0.25

Solution to 1.b

Lagrange (n-1)th order Polynomial

Hence:

f 5 ( x)  L1 ( x) f ( x1 )  L2 ( x) f ( x 2 )  L3 ( x) f ( x3 )  L4 ( x) f ( x4 )  L5 ( x) f ( x5 )  L6 ( x) f ( x6 )

where

n x  xj
Li ( x )  
j 1 xi  x j
j i

( x  x 2 )( x  x 3 )( x  x 4 )( x  x5 )( x  x 6 )
L1 ( x) 
( x1  x 2 )( x1  x 3 )( x1  x 4 )( x1  x5 )( x1  x 6 )

( x  x1 )( x  x3 )( x  x 4 )( x  x5 )( x  x 6 )
L2 ( x) 
( x 2  x1 )( x 2  x3 )( x 2  x 4 )( x 2  x5 )( x 2  x 6 )
...

( x  x1 )( x  x 2 )( x  x3 )( x  x 4 )( x  x 5 )
L6 ( x ) 
( x 6  x1 )( x 6  x 2 )( x 6  x 3 )( x 6  x 4 )( x 6  x 4 )

xi f(xi)
1  637
 0.5  96.5
 0.25  20.5
0.25 20.5
0.5 96.5
1 637

Hence:
L1(0.1) = 0.008064
 L2(0.1) =  0.14784
L3(0.1) = 0.405504
L4(0.1) = 0.946176
L5(0.1) =  0.22176
L6(0.1) = 0.009856

V (0.1)  f 5 (0.1)  L1 (0.1) f ( x1 )  L2 (0.1) f ( x 2 )  L3 (0.1) f ( x 3 )  L4 (0.1) f ( x 4 )  L5 (0.1) f ( x 5 )  L6 (0.1) f ( x 6 )

V (0.1)  L1 (0.1)(637)  L 2 (0.1)(96.5)  L3 (0.1)(20.5)  L 4 (0.1)(20.5)  L5 (0.1)(96.5)  L6 (0.1)(637)

Solution to 1.c

Using the results of part (a), the 4th order Newton interpolating polynomial:

𝑓 𝑥 𝑏 𝑏 𝑥 𝑥 𝑏 𝑥 𝑥 𝑥 𝑥 ... 𝑏 𝑥 𝑥 𝑥 𝑥 𝑥 𝑥 𝑥 𝑥

where

b1 =  637
b2 = 1081
b3 =  1036
b4 = 592
b5 = 0
and
x1 =  1, x2 =  0.5, x3 =  0.25, x4 = 0.25, x5 = 0.5.

Hence:
V 4 (i )  637  1081(i  (1))  (1036)(i  (1))(i  (0.5))  592(i  (1))(i  (0.5))(i  (0.25))  0

V 4 (i )  637  1081(i  1)  1036(i  1)(i  0.5)  592(i  1)(i  0.5)(i  0.25)

V4(1) = 637

In this case the prediction by extrapolation yields the exact/true value of V(1).
 Question 2 (Numerical Integration) (30 points)

A wind force distributed against the side of a skyscraper is measured as:

Height l, m 0 30 60 90 120 150 180 210 240

1 b
b  a a
The net force is defined as F (l )dl . Compute the net force using:
a) Composite trapezoidal rule
b) Composite Simpson’s 1/3 rule
c) Combination of Simpson’s 3/8 rule and Simpson’s 1/3 rule.

Solutions to Q2:

Part 2.(a) Trapezoidal rule:

ℎ
𝐼 𝑓 𝑥 2∑ 𝑓 𝑥 𝑓 𝑥 h = 30

I = (30/2) [ 0 + 2(340 + 1200 + 1550 + 2700 + 3100 + 3200 + 3500) + 3750 ] = 523950

Net Force = I / (b – a) = 523950 / (240 – 0) = 2183.125

Part 2.(b) Composite Simpson’s 1/3 rule

ℎ
𝐼 𝑓 𝑥 4∑ 𝑓 𝑥 2∑ 𝑓 𝑥 𝑓 𝑥 , h = 30
, odd , even

I = (30/3) [ 0 + 4(340 + 1550 + 3100 + 3500) + 2(1200 + 2700 + 3200) + 3750 ] = 519100.

Net Force = I / (b – a) = 519100 / (240 – 0) = 2162.917

Part (c) Since we have 9 points (hence odd-points event-segments), we can begin by applying 3/8 rule
two times followed by one application of the 1/3 rule at the end.
h
Single 1/3 rule: I  f x0   4 f x1   f x 2 
3
3h
Single 3/8 rule: I   f  x0   3 f x1   3 f x2   f  x3 
8
h = 30
I = (30)(3/8) [ 0 + 3(340 + 1200) + 1550 ]
+ (30)(3/8) [ 1550 +3(2700 + 3100) + 3200 ]
+ (30)(1/3) [ 3200 + 4(3500) + 3750 ]
= 528100.

Net Force = I / (b – a) = 528100/ (240 – 0) = 2200.417.

Another alternative solution is to apply the 1/3 rule one time followed by two applications of the 3/8 rule.