Do
No
t S
ha
re
th
is
do
cu
m
en
t
Problem 1 [3 points] Without pivoting. Multiply the first row by 1/3 = 3.3333× 10−1
and subtract from second and multiply by 2 and subtract from third:
3 4 33.6667 −2.0000
−5 1
x1x2
x3
=
103.6667
−5
Multiply second row by −5/3.6667 = −1.3636 and subtract from third:
3 4 33.6667 −2.0000
−1.7272
x1x2
x3
=
103.6667
−8.7880× 10−5
Then
x3 = −8.7880× 10−5/(−1.7272)
= 5.0880× 10−5
x2 = (b2 − a23x3)/a22
= (3.6667 + (−2.0000)× 5.0880× 10−5)/3.6667
= (3.6667− 1.0176× 10−4)/3.6667
= 3.6666/3.6667
= 9.9997× 10−1
x1 = (b1 − a12x2 − a13x3)/a11
= (10− 4× 9.9997× 10−1 − 3× 5.0880× 10−5)/3
= (10− 3.9999− 1.5264× 10−4)/3
= (6.0001− 1.5264× 10−4)/3
= 5.9999/3
= 2.0000.
With partial pivoting. Multiply third row by 0.5 and subtract from first, and multiply by
1/6 = 1.6667× 10−1 and subtract from second:
2.5 −5.0000× 10−14.5 −2.1667
6 3 7
x1x2
x3
=
2.54.5
15
(1)
Multiply second row by 2.5/4.5 = 5.5556× 10−1 and subtract from first:
7.0373× 10−14.5 −2.1667
6 3 7
x1x2
x3
=
−2.0000× 10−54.5
15
(2)
1
Do
No
t S
ha
re
th
is
do
cu
m
en
t
Then
x3 = −2.0000× 10−5/7.0373× 10−1
= −2.8420× 10−5
x2 = (4.5− 4.5×−2.8420× 10−5)/4.5
= (4.5 + 1.2789× 10−4)/4.5
= 4.5001/4.5
= 1
x1 = (15− 3× 1− 7×−2.8420× 10−5)/6
= (15− 3 + 1.9894× 10−4)/6
= (12 + 1.9894× 10−4)/6
= 1.2000× 101/6
= 2
Problem 2 [3 points] The inverse of A is
A−1 =
1
�2
[
1 −1− �
−1 + � 1
]
.
In the infinity norm, for small �,
cond(A) = ‖A‖∞ · ‖A−1‖∞ = (2 + �)
(2 + �)
�2
=
(2 + �)2
�2
≈
4
�2
.
Hence for � ≈
√
�mach,
cond(A) ≈
4
�mach
.
In double precision the machine epsilon is ≈ 2.2204× 10−16 and and hence
cond(A) ≈ 1.8014× 1016.
The script eps_system.m produces
eps |x1-1| |x2-eps|/eps cond(A)
0.9*sqrt(eps) 1.00e+00 7.46e+07 5.81e+16
1.0*sqrt(eps) 0.00e+00 0.00e+00 1.80e+16
1.5*sqrt(eps) 1.11e-01 4.97e+06 8.01e+15
2.0*sqrt(eps) 0.00e+00 0.00e+00 4.50e+15
If the condition number is ≈ 10k, we expect to lose about k digits of accuracy.
However, how to explain the second and fourth lines? With partial pivoting, we have[
1 1 + �
0 1− (1 + �)(1− �)
] [
x1
x2
]
=
[
1 + (1 + �)�
1−
(
1 + (1 + �)�
)
(1− �)
]
.
2
Do
No
t S
ha
re
th
is
do
cu
m
en
t
Consider the (2,2) entry of the matrix:
1− (1 + �)(1− �) = �2.
When � is small, (1 + �)(1− �) ≈ 1. If in the evaluation of (1 + �)(1− �) there are roundoff
errors, after the cancellations the result can be very inaccurate. However, for � =
√
�mach
and 2
√
�mach, 1 − (1 + �)(1 − �) is the same as �2, so no errors. Try computing for various
values the error
1− (1 + �)(1− �)
�2
.
Problem 3 [12 points] The relative error in the solution is
‖x− x̃‖
‖x̃‖
≤ cond(A)
‖r‖
‖b‖
,
where r = b−Ax̃. Gauss elimination with partial pivoting produces relative residual which
is about the machine epsilon, in double precision of the order of 10−16. For cond(A) ≈ 10k,
the error is at most 10k−16, so we expect about k − 16 accurate digits.
With no pivoting, the residual can be larger, so the larger errors in the third column.
Random matrices of size 2000
n A\b no pivot. pivoting cond(A)
———————————–
1 1.4e-13 1.9e-12 1.2e-13 4.4e+03
2 8.9e-14 2.1e-12 9.5e-14 3.9e+03
3 2.0e-13 1.6e-12 1.5e-13 4.0e+03
4 1.4e-12 9.8e-10 3.3e-13 1.8e+05
5 4.9e-13 9.6e-12 5.3e-13 2.2e+04
Problem 4 [4 points] To have degree n = 5, we need (n + 1) points. Let h = 1/5 = 0.2.
Since ex is bounded by e on [0, 1], the error is
M
4(n+ 1)
hn+1 ≤
e
4(5 + 1)
0.25+1 ≈ 7.2488× 10−6.
Let h = 1/n. Then we want
M
4(n+ 1)
hn+1 ≤
e
4(n+ 1)
(1/n)
n+1 ≤ 10−8.
By trial and error, one can easily find that n ≥ 8, so we need at least 9 equally spaced
points.
3
Do
No
t S
ha
re
th
is
do
cu
m
en
t
Problem 5 [5 points]
a. The script interp_sqrt.m produces
sqrt(0.05) ~ 1.024700: error 4.923e-06
sqrt(0.15) ~ 1.072350: error 3.053e-05
b. We have a polynomial of degree n = 3 and spacing h = 0.1. The error bound is
|
√
x+ 1− p3(x)| ≤
M
4(n+ 1)
hn+1, (3)
where
|f (4)(x)| ≤M for all x ∈ [a, b].
The 4th derivative of f(x) =
√
x+ 1 is
√
x+ 1
(4)
=
−15
16(x+ 1)7/2
On [0, 0.3]
|
√
x+ 1
(4)
| =
15
16(x+ 1)7/2
≤
15
16(0 + 1)7/2
= 0.9375 = M.
The error bound is
|
√
x+ 1− p3(x)| ≤
M
4(n+ 1)
hn+1 =
0.9375
4× 4
0.14 ≈ 5.8594× 10−6.
c. The first error is within the error bound, but not the second. The reason is the y values
are approximations of
√
x+ 1 and accurate up to 4 digits. Replace in this script y by
y = sqrt(x+1) and compare the errors.
Problem 6 [4 points] Run main_interp_absx.m.
Problem 7 [3 points] Run main_interp_sinx.m.
In the bound for the error, we have the (n + 1)st derivative of the function f(x) being
interpolated, and this derivative should be continuous for the error bound to hold. When
f(x) = sin(x), |f (n+1)(x)| ≤ 1 for any x, so the small errors. When f(x) = |x|, the first
derivative if discontinuous.
Problem 8 [4 points] Run main_interp.m
4