Hello, dear friend, you can consult us at any time if you have any questions, add WeChat: daixieit

MATH0033: Numerical Methods

Final 7-day coursework, 2019-2020

Answer all questions

1. Let f(x) = x2 − x− 6 and denote by ↵+  and ↵respectively the positive and negative roots of f(x) = 0.

(a)  Suppose that the bisection method is used to nd ↵+  with initial interval [1, 4].

Compute the rst two approximations x0  and x1  and estimate the number of iterations required for the absolute error to be smaller than 106 .

(b) For each of the three functions

'1 (x) = x2 6,    '2 (x) = ^x +6,    '3 (x) = 1 + 6/x  (for x  0),

determine whether the corresponding xed point iteration xk+1   =  '(xk ) can converge to ↵+  or to ↵, given a sufficiently good initial guess.   Justify your answers carefully.

(c) For those combinations of function and root for which convergence is possible, determine a suitable interval [a,b] such that the xed point iteration is guaranteed to converge to the root in question for every initial guess x0  2 [a,b].

Solution

(a)  First note that the roots of f(x) = 0 are +  = 3 and −  = −2.

With  initial  interval  [1 , 4]  the  bisection  method  produces  the  approximations x0  = (4 + 1)/2 = 2 .5 and x1  = (4 + 2 .5)/2 = 3 .25.

The error after k iterations is guaranteed to satisfy

|xk  − ↵+ | <   =  ,

so to ensure that  |xk  + | < 10 6  it suces to take 2k+1  > 3 ⇥ 106 , i .e .

log 3 + 6log 10

(b)  Consider rst '1 (x) = x2 − 6 .  The equation '1 (x) = x is globally consistent with

f(x) = 0 so both roots ↵± are xed points of '1 .  Furthermore, '1 is continuously

di↵erentiable .   However,  '(x)  =  2x,  so  that  |'1 (↵+ )|  =  |'1 (3)|  =  6  > 1  and

 

to either root .

Next consider '2 (x) = − ^x + 6 .  Global consistency now fails, since the positive

root  ↵+   is  not  a  xed  point  of  '2 ,  so  ↵+   certainly  cannot  be  found  by  the

 

−  > 6 .  Moreover, since '(x) = −1/(2^x +6) we have |'2 (↵)| =  |'2 (−2)| =

1/(2^4) = 1/4 < 1 so the xed point iteration converges to −  for a suciently

good initial guess .

Finally consider '3 (x) = 1 + 6/x.  The equation '3 (x) = x is globally consistent

with  f(x)  =  0  so  both  roots  ↵±   are  xed  points  of  '3 .    Furthermore,  '3   is

continuously di↵erentiable for x  0 with '(x) = −6/x2 .  Hence the xed point

 

3/2  > 1 .   But  it  can  converge  to  ↵+  for  a  sufficiently  good  initial  guess  since

|'(↵+ )| = |'(3)| = 6/32  = 2/3 < 1.

(c)  By the CMT the iteration for '2  converges to ↵−  = −2 for all x0  2 [−3 , −1] .  To justify this, note that  (i) '2  : [−3 , −1] ! [−3 , −1] since '2  is strictly decreasing on  [−3 , −1]  and  '2 (−3)  =  − ^3  2  [−3 , −1]  and  '2 (−1)  =  − ^5  2  [−3 ,  −1];

and  (ii)  '2  is  a  strict  contraction  on  [−3 , −1]  since  |'(x)|  =  |1/(2^x +6)|  < 1/(2^5) < 1 for all x 2 [−3 , −1] .

By the  CMT the  iteration for  '3  converges to  ↵+  =  3 for  all x0  2  [2 .5 , 4] .   To justify this, note that  (i) '3  :  [2 .5 , 4] !  [2 .5 , 4] since '3  is strictly decreasing on [2 .5 , 4] and '3 (2 .5) =  1 + 6/2 .5 =  17/5  2  [2 .5 , 4] and '3 (4) =  1 + 6/4 = 5/2  2 [2 .5 , 4];  and  (ii)  '3  is  a  strict  contraction  on  [2 .5 , 4]  since  |'(x)|  =  |6/x2 |  < 6/(2 .5)2  = 24/25 < 1 for all x 2 [2 .5 , 4] .


Solution

(a) Taylor-expanding f(↵) about x = xk  gives

0 = f(↵) = f(xk )+ f\ (xk )(↵ − xk )+ f\\ (⇠k )(↵ − xk )2 ,

for some k  between and xk .

Combining this with the definition of the Newton iteration gives

 


xk+1 =


 

f\\ (k )

2f\ (xk )

 


(xk  )2 .


(†)


Given 6 > 0 define I6  := [↵ − 6, ↵ + 6].

If xk  2 I6  then also ⇠k  2 I6 , which implies that

|xk+1 − ↵| <  ·  · 6 · |xk  − ↵| =  |xk  − ↵|.

Hence if x0  2 I6  then by induction it follows that

|xk  | < 2k |x0 |,  k = 0, 1, . . . ,

so that xk  (and hence also ⇠k ) converges to ↵ as k ! 1.               Furthermore, from (†) we have by the continuity of f\  and f\\  that

 


xk+1     f\\ (k )        f\\ ()

(xk  − ↵)2         2f\ (xk )      2f\ (↵)


as k ! 1.


(b) From (*) we know that for every C > |c| there exists k0  2 N such that

k ≥ k0   implies   l l < C,  i.e.  |xk+1 − ↵| < C|xk  − ↵|2 ,

so the convergence is quadratic.

(c) We argue by contrapositive.  Suppose that the convergence is order p for some p > 2. Then there exist k0  2 N and C > 0 such that

k ≥ k0   implies   < C.

But then for k ≥ k0

l l =  · |xk  |p 2  < C|xk  |p 2 ,

which tends to zero as k ! 1. Hence c = 0.


(d) For f(x) = x/(1 + x2 ) we have

1 x2                                              2x(3 x2 )

(1 + x2 )2                                             (1 + x2 )3   .

Hence, for |x|, |y| < 1/4,

 


 

|f\ (y)| ≥

1 (1/4)2    (1 + (1/4)2 )2

=

15/16   (17/16)2

 


 

(Please turn over)


3.  Consider the solution of the Cauchy problem

{ y\ (t) = f(t,y(t)),  t > 0,

y(0) = y0 ,

using the following numerical method:

u0  = y0 ,

 p1  = f(tn ,un )

>

For n = 1, 2, . . . ,  < p2  = f(tn + h,un + a1 hp1 )

>

 un+1  = un + h(p1 + a2p2 ),

where a1 ,a2  2 R.

(a) Write the method in one-step form un+1  = un  + h (tn ,un ,un+1;h), specifying the increment function    . Is the method explicit or implicit?

(b) Determine conditions on the constants a1 ,a2   under which the method  (*) is

second order consistent as h ! 0.   State clearly the smoothness assumptions on y and f under which your analysis is valid.

(c) Determine for which h the method (*) is absolutely stable, in the special case where a1  = 2/3 and a2  = 3/4.


Solution

(a) The method can be written as un+1  = un + h (tn ,un ;h) with

(tn ,un ;h) = f(tn ,un )+ a2 f tn + h,un + a1 hf(tn ,un ).

The method is explicit, since     is independent of un+1 .

(This is a two-stage Runge-Kutta method.)

(b) The truncation error of a one-step method un+1  = un + h (tn ,un ;h) is

yn+1 yn

 

where yn  = y(tn ) denotes the exact solution of the Cauchy problem at time tn .  Assuming sufficient smoothness (see later for specific assumptions), Taylor ex- pansion gives (with y denoting y\ (tn ) etc.)

h2               h3

 

 

 

 

 

 

so that

Tn  =(y − g(0)) + 2 (y − 2g\ (0)) +  6  (y\\\ (tn + ⇠) − 3g\\ (⌘)) .    (**)

Now,

g(0) =  + a2 f(tn ,yn ),

and f(tn ,yn ) = y (since y\  = f), so the O(1) term in (**) vanishes provided

 

a2  = 3

Also, by the chain rule,

g\ (h) =    (t,y)+ a1 f(tn ,yn )  (t,y),

where t⇤  = tn + h and y⇤  = yn + a1 hf(tn ,yn ), so that

g\ (0) =    (tn ,yn )+ a1 f(tn ,yn )  (tn ,yn ),

while applying the chain rule to y\  = f gives

?f                  ?f


so that the O(h) term in (**) vanishes provided

 

a1  = 2

Applying the chain rule again we nd that

g\\ (h) =   (t,y)+ f(tn ,yn )  (t,y)+ f(tn ,yn )2  (t,y).

Hence, from (**) we see that if a1  =  and a2  =  , and if y\\\ , f,  ,  and  all exist and are bounded, then there exists a constant C > 0 such that

|Tn | W Ch2 ,

so the method is of second order.

(c) A method is absolutely stable if, when applied to the model problem

for λ < 0, the numerical approximations un  satisfy un  t 0 as n t 8.            Applying the proposed method to this problem (with f(t,y(t)) = λy(t)) gives

un+1  = un + h f(tn ,un )+ f(tn + h,un + hf(tn ,un )

= un + h λun + λ(un + hλun )

1+ hλ + un .

Hence, by induction,

un  = 1+ hλ + n u0 ,    n = 1, 2, . . . .

Thus the method is absolutely stable if and only if

()2

2

The left-hand inequality is always satisfied since the quadratic 2 + σ + σ2 /2 is positive for all real σ (it is positive at infinity and has no real roots). The right- hand inequality holds if and only if

hλ +  = hλ 1+ < 0,

which, since hλ < 0, holds if and only if 1 + hλ/2 > 0, i.e.

h <  2


4. Let A 2 Rn n  be an invertible matrix, b 2 Rn , and let x ⇤  2 Rn  be the unique solution of the equation Ax = b. Consider the two-step iteration

xk+1  = B0xk + B1xk 1 + c,  k = 1, 2, . . . ,

where B0 ,B1  2 Rn n  and c 2 Rn . We say the iteration is consistent if Ax = b , x = B0x + B1x + c.

(a) Rewrite (#) as an equivalent one-step iteration

uk+1  = Buk + d,  k = 1, 2, . . . ,

for uk  = xk(x)1  2 R2n, specifying B 2 R2n⇥2n  and d 2 R2n .

(b) Assume that (#) is consistent. Prove that if the sequence (uk ) dened in part (a)

converges for every initial guess u0  = x(x)0(1)  2 R2n then the limit is independent of u0  and equals  x(x)⇤(⇤)   .

Hint: show that if uk  ! u ⇤   then u ⇤   is a xed point of the map u 7! Bu + d .

(c) Now let B0  = ↵I − βA and B1  = (1 − ↵)I for some ↵ , β 2 R with β  0. Determine the choice of c ensuring that (#) is consistent in this case.

(d) Assuming the choices of B0 , B1  and c from part (c), show that (#) converges to x  for all initial guesses x0 ,x1  2 Rn  if and only if

 


 ± s 2 +1  < 1

                                                      


 

Aλ 2 σ(A).


Solution

(a) The iteration can be written in the prescribed one-step form with

B =  ,    d = 0(c) .

(b) We rst show that if uk  ! u*  then u*  is a xed point of the map u 7! Bu + d. For this we use the definition of the iteration and the fact that u 7! Bu + d is continuous to see that

u*  =  lim uk+1  =  lim (Buk  + d) = B (lim uk )+ d = Bu* + d.

k!1                     k!1                                            k!1

This implies by the denition of B that if u*  =  y(x) then y = x and x =

B0x + B1y + c, which combine to give that x = B0x + B1x + c. By the assumed consistency of (#) we conclude that x = y = x* .

(c)  Consistency for (#) is equivalent to requiring that I − B0 − B1  is invertible and

c = (I − B0 − B1 )x* . With the definitions of B0 , B1  given, we have I − B0 − B1  = I − (↵I − βA) − (1 − ↵)I = βA.

Hence (#) is consistent if and only if c = βAx* , i.e. c = βb.  Note that invert- ibility of I − B0 − B1  follows from that of A because β  0.

(d) The iteration converges for all initial guesses u0  2 R2n  if and only if p(B) < 1.   The constant µ 2 C is an eigenvalue of B if and only if there exists 0  u =  y(x) 2 R2n  such that Bu = µu. Since

B =  ,

this is equivalent to saying that

(↵I − βA)x +(1 − ↵)y = µx     and    x = µy. By the second equation, the rst equation is equivalent to

µβAy = (↵µ +1 − ↵ − µ )y2 .

Since β  0, this implies that for µ  0 (we note that the case µ = 0 is not relevant for calculating the spectral radius)

↵µ +1 µ2

 

for some λ 2 σ(A). Then


which has solutions

µ± (λ) =  ± s ✓2 +1 ↵ .

Hence p(B) < 1 if and only if |µ± (λ)| < 1 for all λ 2 σ(A).