Math 215/255
Lecture 7 revision: Numerical Integration
refresher
Sections: From your 1st year courses, e.g. Math 101
𝑏𝑏
Interpretations and applications of ∫𝑎𝑎 𝑓𝑓 𝑥𝑥 𝑑𝑑𝑑𝑑
› Area under a curve
› Mean value
› Arc length
› Volume
› Centre of mass
› Work
› Solution of separable DE’s
› Here are some you can do:
1
a) ∫0 𝑠𝑠𝑠𝑠𝑠𝑠 𝑥𝑥 𝑑𝑑𝑑𝑑
b) ∫ 𝑒𝑒 −𝑥𝑥 𝑑𝑑𝑑𝑑
c) ∫ 𝑥𝑥 𝑐𝑐𝑐𝑐𝑐𝑐 𝑥𝑥 𝑑𝑑𝑑𝑑
2
› How about this one: ∫ 𝑒𝑒 −𝑥𝑥 𝑑𝑑𝑑𝑑?
– But the area under the curve is well defined
– Can we approximate the integral between 𝑎𝑎
and 𝑏𝑏?
𝑏𝑏
Riemann sum approximation: 𝐼𝐼 = ∫𝑎𝑎 𝑓𝑓 𝑥𝑥 𝑑𝑑𝑑𝑑
› Idea, divide [𝑎𝑎, 𝑏𝑏] into 𝑁𝑁 intervals of width: ℎ = (𝑏𝑏 − 𝑎𝑎)/𝑁𝑁
› Approximate 𝐼𝐼 with the (left) Riemann sum:
𝐿𝐿𝑁𝑁 = ℎ𝑓𝑓 𝑎𝑎 + ℎ𝑓𝑓 𝑎𝑎 + ℎ + ℎ𝑓𝑓 𝑎𝑎 + 2ℎ +. . . +ℎ𝑓𝑓 𝑏𝑏 − 2ℎ + ℎ𝑓𝑓 𝑏𝑏 − ℎ
𝑅𝑅𝑁𝑁 = ℎ𝑓𝑓 𝑎𝑎 + ℎ + ℎ𝑓𝑓 𝑎𝑎 + 2ℎ +. . . +ℎ𝑓𝑓 𝑏𝑏 − 2ℎ + ℎ𝑓𝑓 𝑏𝑏 − ℎ + ℎ𝑓𝑓 𝑏𝑏
› Notes and questions:
– We know that 𝐿𝐿𝑁𝑁 → 𝐼𝐼 as 𝑁𝑁 → ∞, (this was the definition of the integral)
– For a given 𝑓𝑓 𝑥𝑥 , how large does 𝑁𝑁 need to be for 𝐿𝐿𝑁𝑁 to approximate 𝐼𝐼 to a tolerance ε?
– When is numerical integration needed?
Errors for the Riemann sum
› A rigorous estimate is based on the Example: Suppose that we know that
expression: 𝑚𝑚𝑚𝑚𝑚𝑚 = 𝑓𝑓 ′ 𝑥𝑥 < 5
𝑓𝑓 ′ 𝜁𝜁 𝑥𝑥∈ 0,3
𝐼𝐼 − 𝐿𝐿𝑁𝑁 = ℎ 𝑏𝑏 − 𝑎𝑎
2 How many intervals are needed to evaluate
3
where ζ ∈ [𝑎𝑎, 𝑏𝑏] is some unknown value 𝐼𝐼 = � 𝑓𝑓 𝑥𝑥 𝑑𝑑𝑑𝑑
› We can see the following: 0
– As 𝑁𝑁 → ∞, error goes to zero since ℎ = (𝑏𝑏 − to within two decimal places?
𝑎𝑎)/𝑁𝑁 → 0
– Error grows if the interval is long or if the
function has large derivatives
– If 𝑓𝑓(𝑥𝑥) is increasing, then the left Riemann
sum underestimates the integral (and the right
Riemann sum overestimates).
› Suppose we know that:
𝑚𝑚𝑚𝑚𝑚𝑚 = 𝑓𝑓 ′ 𝑥𝑥 < 𝑀𝑀
𝑥𝑥∈ 𝑎𝑎,𝑏𝑏
› Then the error satisfies:
2
ℎ 𝑀𝑀 𝑏𝑏 − 𝑎𝑎
𝐼𝐼 − 𝐿𝐿𝑁𝑁 = 𝑓𝑓 ′ 𝜁𝜁 𝑏𝑏 − 𝑎𝑎 ≤
2 2𝑁𝑁
Trapezoid rule: uses endpoints of the function
on each interval:
h
TN = [ f (a ) + f (a + h )] + h [ f (a + h ) + f (a + 2h )] + ... + h [ f (b − h ) + f (b )]
2 2 2
h h
= f (a ) + hf (a + h ) + hf (a + 2h ) + ... + hf (b − h ) + f (b )
2 2
› Error estimates from the expression:
″
2
𝑓𝑓 𝜁𝜁
𝐼𝐼 − 𝑇𝑇𝑁𝑁 = −ℎ 𝑏𝑏 − 𝑎𝑎
12
› Suppose we can find a bound on 𝑓𝑓 ″ 𝑥𝑥 :
𝑚𝑚𝑚𝑚𝑚𝑚 = 𝑓𝑓 ″ 𝑥𝑥 < 𝑀𝑀
𝑥𝑥∈ 𝑎𝑎,𝑏𝑏
then the error satisfies:
ℎ2 ″ 𝑀𝑀 𝑏𝑏 − 𝑎𝑎 3
𝐼𝐼 − 𝑇𝑇𝑁𝑁 = 𝑓𝑓 𝜁𝜁 𝑏𝑏 − 𝑎𝑎 ≤
12 12𝑁𝑁 2
Simpson’s rule: based on a quadratic
approximation on each interval
h
SN = [ f (a ) + 4 f (a + h ) + f (a + 2h )] + h [ f (a + 2h ) + 4 f (a + 3h ) + f (a + 4h )] + ...
3 3
h h
... + [ f (b − 4h ) + 4 f (b − 3h ) + f (b − 2h )] + [ f (b − 2h ) + 4 f (b − h ) + f (b )]
3 6
h
= [ f (a ) + 4 f (a + h ) + 2 f (a + 2h ) + 4 f (a + 3h ) + ... + 2 f (b − 2h ) + 4 f (b − h ) + f (b )]
3
› Note the repeating pattern of coefficients: 142424…..24241, and we need 𝑁𝑁 to be even.
› We make error estimates from the expression:
𝑓𝑓 𝑖𝑖𝑖𝑖 𝜁𝜁
𝐼𝐼 − 𝑆𝑆𝑁𝑁 = −ℎ4 𝑏𝑏 − 𝑎𝑎
180
› Suppose we can find a bound on 𝑓𝑓 𝑖𝑖𝑖𝑖 𝑥𝑥 : i.e. 𝑚𝑚𝑚𝑚𝑚𝑚 = 𝑓𝑓 𝑖𝑖𝑖𝑖
𝑥𝑥 < 𝑀𝑀.
𝑥𝑥∈ 𝑎𝑎,𝑏𝑏
› Then the error satisfies:
ℎ4 𝑖𝑖𝑖𝑖
𝑀𝑀 𝑏𝑏 − 𝑎𝑎 5
𝐼𝐼 − 𝑆𝑆𝑁𝑁 = 𝑓𝑓 𝜁𝜁 𝑏𝑏 − 𝑎𝑎 ≤
180 180𝑁𝑁 4
› Lets evaluate:
ORDER OF 1
APPROXIMATION: 𝐼𝐼 = � 𝑠𝑠𝑠𝑠𝑠𝑠 𝑥𝑥 𝑑𝑑𝑑𝑑 = 1 − 𝑐𝑐𝑐𝑐𝑐𝑐 1 ≈ 0.4597
0
Theoretical estimates: › Use left Riemann sum, Trapezoidal
and Simpson’s method, e.g. 𝑁𝑁 = 4:
In other words:
› Increase 𝑁𝑁 : what happens to the
error?
When the errors behave in
this way, i.e.
Error
𝑙𝑙𝑙𝑙𝑙𝑙 = 𝐶𝐶𝑝𝑝
ℎ→0 ℎ𝑝𝑝 › Double 𝑁𝑁, the error decreases by a
method is of p-th order factor of 2 for left Riemann sum, by a
factor of 4 for Trapezoidal, by 16 for
Simpson’s rule
Using Matlab for numerical approximation of
integrals
1
› Example: 𝐼𝐼 = ∫0 𝑠𝑠𝑠𝑠𝑠𝑠 𝑥𝑥 𝑑𝑑𝑑𝑑
› Set up vectors and assign
function values (𝑁𝑁 = 8):
› Left Riemann sum
› Trapezoid rule has special
Matlab function:
› Simpson’s rule:
– Define coefficient vector: simpc
– Now multiply and evaluate sum
Polynomial Approximations via Taylor Series
› Taylor series expansion:
2 3
𝑥𝑥 − 𝑎𝑎 𝑥𝑥 − 𝑎𝑎
𝑓𝑓 𝑥𝑥 = 𝑓𝑓 𝑎𝑎 + 𝑥𝑥 − 𝑎𝑎 𝑓𝑓 ′ 𝑎𝑎 + 𝑓𝑓 ″ 𝑎𝑎 + 𝑓𝑓 ‴ 𝑎𝑎 +. . .
2! 3!
› Truncate this series to give an approximation formula, e.g. at
first order the linear approximation (effectively the tangent
line at x=a)
𝑓𝑓 𝑥𝑥 ≈ 𝑓𝑓 𝑎𝑎 + 𝑥𝑥 − 𝑎𝑎 𝑓𝑓 ′ 𝑎𝑎
› Can we make approximation better by including more terms,
e.g. quadratic approximation:
𝑥𝑥 − 𝑎𝑎 2
𝑓𝑓 𝑥𝑥 ≈ 𝑓𝑓 𝑎𝑎 + 𝑥𝑥 − 𝑎𝑎 𝑓𝑓 ′ 𝑎𝑎 + 𝑓𝑓 ″ 𝑎𝑎
2!
Errors? Quadratic
Linear
› Uses the remainder form of Taylor
series
› Example: linear approximation
– For some 𝜁𝜁: 𝑎𝑎 < 𝜁𝜁 < 𝑥𝑥, we have that
′
𝑥𝑥 − 𝑎𝑎 2 ″
𝑓𝑓 𝑥𝑥 − [𝑓𝑓 𝑎𝑎 + 𝑥𝑥 − 𝑎𝑎 𝑓𝑓 𝑎𝑎 ] = 𝑓𝑓 𝜁𝜁
2!
– The expression on the left is the difference
between the function and the linear
approximation.
› Notes:
– For a larger interval from the point 𝑥𝑥 = 𝑎𝑎
the approximation becomes worse
– If the second derivative is large the
approximation is poor
– The sign of the second derivative
determines if the linear approximation
over/under predicts
EXAMPLE
Use linearly truncated
Taylor series at 𝑥𝑥 = 𝑎𝑎 and
𝑥𝑥 = 𝑏𝑏 to estimate the error
in a Trapezoidal
approximation to
𝑏𝑏
� 𝑓𝑓(𝑥𝑥) 𝑑𝑑𝑑𝑑
𝑎𝑎
Can you use this to derive
the previous error estimate
for the Trapezoid rule?
Numerical differentiation
› Definition of the derivative:
′
𝑓𝑓 𝑎𝑎 + ℎ − 𝑓𝑓 𝑎𝑎
𝑓𝑓 𝑎𝑎 = 𝑙𝑙𝑙𝑙𝑙𝑙
ℎ→0 ℎ
› How do we approximate this from numerical or
experimental data?
› Two simple ideas
𝑓𝑓 𝑎𝑎+ℎ −𝑓𝑓 𝑎𝑎
– Forward Euler difference: 𝑓𝑓 ′ 𝑎𝑎 ≈
ℎ
𝑓𝑓 𝑎𝑎 −𝑓𝑓 𝑎𝑎−ℎ
– Backward Euler difference: 𝑓𝑓 ′ 𝑎𝑎 ≈
ℎ
› Both based on the definition of the derivative
EXAMPLE
Approximate the derivative
of sin 𝑥𝑥 at 𝑥𝑥 = 1. Use both
FE and BE, with
1 1 1 1 1
ℎ= , , , ,
10 20 40 80 160
Compare the errors
› Error is approximately linear in ℎ:
first order method
› Errors of FE and BE are of similar
magnitude but opposite sign.
Central difference method:
› Obtained by averaging FE and BE
′
𝑓𝑓 𝑎𝑎 + ℎ − 𝑓𝑓 𝑎𝑎 − ℎ
𝑓𝑓 𝑎𝑎 ≈
2ℎ
› Error is obtained from a Taylor expansion
Deriving finite difference formulas using Taylor
series
Find an approximation to 𝑓𝑓𝑓𝑓(𝑎𝑎) if 𝑓𝑓(𝑎𝑎), 𝑓𝑓(𝑎𝑎 − ℎ), 𝑓𝑓(𝑎𝑎 + ℎ),
are known. What is the order of your approximation?