CHAPTER 2

Solutions of Equations of One Variable

Theorem (Intermediate Value theorem (IVT)). If f 2 C[a, b] and K is any number between f (a) and f (b), then there exists c 2 (a, b) such that f (c) = K.

In particular, if f (a)f (b) < 0, then there exists c in (a, b) such that f (c) = 0. Thus, one can use the IVT as a root locator. 2. The Bisection Method Problem: find the roots of f (x) = 0 where f is continuous on an open interval containing the roots. Solution: first graph to locate roots by IVT. For each root p, find a closed interval [a, b] where p 2 (a, b) and f (a) · f (b) < 0. Make sure these intervals do not overlap. Note. It is often better to use sgn(f (a)) · sgn(f (b)) < 0 rather than f (a) · f (b) < 0. 8 > < 1 x :1 x>0 16

2. THE BISECTION METHOD

17

Graphical idea:

f(a1)f(b1) 1, to a root of f (x) = 0 is computed from the approximations pn and pn 1 using the equation f (pn)(pn pn 1) pn+1 = pn . f (pn) f (pn 1) Example. Find a root of ex + 2 x + 2 cos x 6 = 0 on [1, 2] to within 10 6. Use the following keystrokes to do this on the TI-89 (where ) means STO.): e^(x)+2^(-x)+2cos(x)-6)y1(x) 1.0)xx 2.0)yy xx)zz:yy)xx:xx-(y1(xx)(xx-zz))/(y1(xx)-y1(zz)))yy [keep repeating this last line by hitting ENTER] n 0 1 2 3 4 5 6 7 8

pn 1.0 2.0 1.67830848477 1.80810287702 1.83229846352 1.82933117293 1.82938347398 1.82938360194 1.82938360193

Thus p ⇡ 1.829384 and f (p) ⇡ 1.6 ⇥ 10 6. ⇤ Note that p is not between p2 and p3. What if we want to keep our root between successive approximations, i.e., want |pn+1

p| < |pn+1

pn|,

providing an error bound? We use a hybrid of the secant and bisection methods.

26

2. SOLUTIONS OF EQUATIONS OF ONE VARIABLE

Method of False Position (Regula Falsi) As previously, choose p0 and p1 such that p0 < p1 and f (p0) · f (p1) < 0. Then an interval [an+1, bn+1], for n > 1, containing an approximation to a root of f (x) = 0 is found from an interval [an, bn] containing the root by first computing f (an)(bn an) pn+1 = an . f (bn) f (an) Then set an+1 = an and bn+1 = pn+1 if f (an)f (pn+1) < 0, and an+1 = pn+1 and bn+1 = bn otherwise. Example. Find a root of ex + 2 x + 2 cos x 6 = 0 on [1, 2] to within 10 6. From the previous example, we know the root is in [1.8, 1.9], so we start there. n 0 1 2 3 4 5 6

an 1.8 1.8276767517 1.8292856191 1.8293779810 1.8293832795 1.8293835834 1.8293836009

bn 1.9 1.9 1.9 1.9 1.9 1.9 1.9

pn+1 1.8276767517 1.8292856191 1.8293779810 1.8293832795 1.8293835834 1.8293836009 1.8293836019

f (an)f (pn+1) + + + + + +

Since f (1.8293836) · f (1.8293837) < 0, |p

p6| < |1.8293837

1.8293836| = .00000001 = 1 ⇥ 10 8.

Thus we approximate p by p6 = 1.8293836009 and f (p6) ⇡ 4.2 ⇥ 10 9. Although the Method of False Position provides us with error bounds, the method converges more slowly than the secant method and possibly even more slowly than the bisection method. ⇤

3. THE SECANT METHOD

27

The pseudocode algorithms for the Secant and False Position Methods follow.

Maple. See secantalg.mw or secantalg.pdf for the implementation of this pseudocode as a Maple procedure. To use the algorithm in problem solving, see secant.mw or secant.pdf.

28

2. SOLUTIONS OF EQUATIONS OF ONE VARIABLE

Maple. See falsepositionalg.mw or falsepositionalg.pdf for the implementation of this pseudocode as a Maple procedure. To use the algorithm in problem solving, see falseposition.mw or falseposition.pdf.

4. NEWTON’S METHOD

29

4. Newton’s Method The line that best approximates the graph of a function at a point is the tangent line at the point. Using this line instead of the secant line gives us Newton’s Method (also known as the Newton-Raphson Method). Suppose p is a root of f (x) = 0, and f 0(x) exists in an interval containing all approximations to p. Let p0 be a first approximation. The equation of the tangent line at (p0, f (p0)) is y

f (p0) = f 0(p0)(x

p0)

and crosses the x-axis at p1 if f 0(p0) 6= 0. Then 0

f (p0) = f 0(p0)(p1

p0) =) p1 = p0

f (p0) . f 0(p0)

Continuing, a sequence {pn} of successive approximations is generated.

30

2. SOLUTIONS OF EQUATIONS OF ONE VARIABLE

Newton’s Method The approximation pn+1 to a root of f (x) = 0 is computed from the approximation pn using the equation pn+1 = pn

f (pn) . f 0(pn)

We will usually stop our sequence of approximations when |pn+1 for some prescribed tolerance ✏. Example. Find a root of 2 + cos(ex

pn| < ✏

ex on [.5, 1.5].

2)

We have f (x) = 2 + cos(ex

2)

ex,

f 0(x) =

2)

ex,

so ex sin(ex

giving pn+1 = pn

2 + cos(epn epn sin(epn

2) 2)

epn . epn

Choose p0 = 1.0. You can use the following keystrokes to do this on the TI-89 (where ) means STO.): 2+cos(e^(x)-2)-e^(x))y1(x) d(y1(x),x) ENTER ans(1))y2(x) 1.0)xx xx-y1(xx)/y2(xx))xx [keep repeating this last line by hitting ENTER]

4. NEWTON’S METHOD

n 0 1 2 3 4

31

pn 1.0 1.00768903626 4 place accuracy 1.00762397633 8 place accuracy 1.00762397166 11 place accuracy 1.00762397166

Thus p ⇡ 1.00762397166 and f (p) ⇡

1 ⇥ 10

11

.

Compare this with the bisection method, where we had three decimal place accuracy after 10 steps. ⇤ Note. 1. We need p0 close to p, or {pn} may not converge or converge to a di↵erent root. 2. We want f 0 bounded away from 0 near p 3. If points 1 and 2 are met and f 2 C 2[a, b], then {pn} ! p.

Example. f (x) = x1/3 has a single root at x = 0, but Newton’s method 1 will not work. Since f 0(x) = x 2/3, 3 1/3 pn pn+1 = pn = pn 3pn = 2pn, 2/3 1 3 pn so each approximation is twice as far from the origin as the previous one and in the opposite direction. ⇤

32

2. SOLUTIONS OF EQUATIONS OF ONE VARIABLE

The pseudocode for the Newton-Raphson Method follows.

Maple. See newtonalg.mw or newtonalg.pdf for the implementation of this pseudocode as a Maple procedure. To use the algorithm in problem solving, see newton.mw or newton.pdf.

5. ERROR ANALYSIS AND ACCELERATING CONVERGENCE

33

5. Error Analysis and Accelerating Convergence Definition. Suppose{pn}1 n=0 ! p, pn 6= p for all n. If there exist ↵ > 0, and M > 0 with |pn+1 p| = n!1 |pn p|↵ lim

> 0,

 M,

then {pn} converges to p of order ↵ with asymptotic error constant . An iterative technique of the form pn+1 = g(pn) is said to be of order ↵ if the sequence {pn} converges to the solution p = g(p) of order ↵. ↵ = 1 is linear convergence: |pn+1

p|  M |pn

p| for large n.

↵ = 2 is quadratic convergence: |pn+1

p|  M |pn

p|2 for large n.

The order of convergence is most important here. Problem (Page 56 #5(b)(i-ii)). Show the sequence {pn} ! 0 linearly for {pn} = n12 , n 1. lim

1 (n+1)2

n!1

1 n2

0 0

1

⇣ n ⌘2 = lim = 12 = 1 n!1 n + 1

Thus we have linear convergence with error constant 1. What if we check for quadratic convergence? lim

n!1

1 (n+1)2 1 n2

0 0

2

n4 = lim =1 n!1 (n + 1)2

So a limit of infinity means we need to check lower. In the same way, if the limit is 0, we need to check higher.

34

2. SOLUTIONS OF EQUATIONS OF ONE VARIABLE

What value of n gives |pn p|  5 ⇥ 10 2? 1 2 2 0  5 ⇥ 10 = .05 () n 20 (= n n2 Example (Linear vs. Quadratic Convergence).

5.

Suppose {pn} ! p = 0 linearly and {ˆ pn} ! p = 0 quadratically with the constant M = 0.5 in each case. Linear Case |pn

0| = |pn|  .5|pn 1|  .52|pn 2|  · · ·  .5n|p0|

Quadratic Case |ˆ pn

0| = |ˆ pn|  .5|ˆ pn 1|2

⇥ ⇤2  .5 .5|ˆ pn 2|2 = .53|ˆ pn 2|4 ⇥ ⇤4  .53 .5|ˆ pn 3|2 = .57|ˆ pn 3|8

Now suppose |p0| = |ˆ p0| = 1. n 1 2 3 4 5 6 7

pn = .5n 0.5 0.25 0.125 0.0625 0.03125 0.015625 0.0078125

n

 · · ·  .52

1

n

|ˆ p0|2

n

pˆn = .52 1 0.5 0.125 0.0078125 3.0518 ⇥ 10 4.6566 ⇥ 10 1.0842 ⇥ 10 5.8775 ⇥ 10

5 10 19 39

Note. While bisection provides only linear convergence, Newton’s method provides quadratic convergence.

5. ERROR ANALYSIS AND ACCELERATING CONVERGENCE

Aitkin’s

2

35

method (for accelerating convergence)

Linearly convergent methods tend to be more common than quadratically convergent methods. Aitkin’s 2 method can be used to accelerate the convergence of a linearly convergent sequence. Theorem (Aitkin’s

then {qn}1 n=0

2

method). If {pn}1 n=0 ! p linearly, and if

(pn+1 pn)2 qn = pn , pn+2 2pn+1 + pn ! p also and, in general, more rapidly.

Problem (Page 54 #1c). Generate the first 5 terms of the sequence {qn} for pn = 3 pn 1 (n 1) if p0 = 0.5. n pn qn 0 .500000 .548101 (3) 1 .577350 (1) .547915 (5) 2 .530315 (2) .547847 (7) 3 .558439 (4) .547823 (9) 4 .541448 (6) .547814 (11) 5 .551650 (8) .547811 (13) 6 .545502 (10) 7 .549199 (12) We achieve in 13 steps what would take about 20 steps without the use of Aitkin’s 2 method. To 6 decimal places, the solution is .547809.

36

2. SOLUTIONS OF EQUATIONS OF ONE VARIABLE

Why is this called a

2

method? For {pn}1 n=0 , define the

forward di↵erence k

pn =

(

k 1

pn = pn+1

pn), k

pn for n

0;

2.

Then 2

pn =

( pn) =

(pn+1 pn) = pn+1 pn = (pn+2 pn+1) (pn+1 pn) = pn+2

so qn = pn

(pn+1 pn)2 = pn pn+2 2pn+1 + pn

2pn+1 + pn,

( pn)2 . 2p n

1 The sequence {qn}1 n=0 ! p more rapidly than {pn }n=0 ! p in the following sense:

Theorem (Aitken’s 2 Convergence). If {pn}1 n=0 ! p linearly and (pn p)(pn+1 p) > 0 for large values of n, and qn = pn then

qn n!1 pn lim

( pn)2 , 2p n p = 0. p

5. ERROR ANALYSIS AND ACCELERATING CONVERGENCE

Note the term construction order for Aitkin’s

2

37

method:

p0, p1, p2, q0, p3, q1, p4, q2, . . . Ste↵ensen’s Method readjusts this term construction Method: (0)

p0 ,

(0)

p1 = g(p0),

(0)

p2 = g(p1),

(1)

p0 = {

2

(0)

}(p0 )

(1)

is the same as Aitkin’s method, but the next step uses q0 = p0 instead of (0) p2 = p2 . Continuing, every third term uses Aitkin. (1)

(1)

p1 = g(p0 ),

(1)

(1)

p2 = g(p1 ),

Problem (Page 54 #1c (continued)). k 0 1 2

(k)

p0 .500000 .548101 .547809

(2)

p0 = {

(k)

p1 .577350 .547633 .547809

2

(1)

}(p0 ),

...

(k)

p2 .530315 .547915

We now achieve the solution in 6 steps, with a 7th giving verification. Note. (1) If any imation. (2) {

2

2

(n 1)

pn = 0, we terminate the sequence and take p2

(1)

(1)

(1)

(1)

}(p0 ) means apply Aitkin to p0 , p1 , p2 , etc., not

as the approx2 (1) p0 .

38

2. SOLUTIONS OF EQUATIONS OF ONE VARIABLE

The pseudocode for the Ste↵ensen’s Method follows.

Maple. See ste↵ensenalg.mw or ste↵ensenalg.pdf for the implementation of this pseudocode as a Maple procedure. To use the algorithm in problem solving, see accellconverge.mw or accellconverge.pdf.

6. MULLER’S METHOD

39

6. Muller’s Method Roots of Polynomials Definition. A solution of f (x) = 0 is a zero of multiplicity m if for x 6= p, f (x) = (x

p)mq(x)

where lim q(x) 6= 0.

x!p

A simple root is a root of mutiplicity one. Theorem. f 2 C 1[a, b] has a simple zero at p 2 (a, b) () f (p) = 0 but f 0(p) 6= 0. Theorem. f 2 C m[a, b] has a zero of multiplicity m at p 2 (a, b) () but f (m)(p) 6= 0.

0 = f (p) = f 0(p) = · · · = f (m

1)

(p)

Note. Zeros of multiplicity greater than 1 may not reach quadratic convergence with Newton’s method. A polynomial of degree n is of the form P (x) = anxn + an 1xn

1

+ · · · + a1x + a0,

an 6= 0.

The polynomial P (x) = 0 has no degree. a Theorem (Rational Root). If 2 Q is a rational root of P (x) = 0 where b a, b 2 Z, then a is a divisor of a0 and b is a divisor of an. Theorem (Fundamental Theorem of Algebra). If P (x) is a polynomial of degree n 1 with real or complex coefficients, then P (x) = 0 has at least one (possibly) complex root.

40

2. SOLUTIONS OF EQUATIONS OF ONE VARIABLE

Corollary. Allowing for the use of complex coefficients, P (x) = an(x

x1)m1 (x

x2)m2 · · · (x

where

xk )mk

m1 + m2 + · · · + mk = n and the xi are real or complex roots, simple or multiple. For real coefficients, any complex roots come in conjugate pairs a ± bi. Horner’s Method or synthetic substitution This is a method used to divide P (x) by x P (x) = (x

x0 to get

x0)Q(x) + P (x0)

where Q(x) is the quotient and P (x0) is the remainder or to find P (x0). Problem (Page 58 # 2b). Suppose P (x) = x4 2x3 12x2 + 16x 40 n⇥ o ⇤ = (x 2)x 12 x + 16 x 40

Find P ( 2) or divide P (x) by x ( x0 = 2 1 2 2 1 4 Then P (x) = (x + 2) (x3 | Note.

2) = x + 2. 12 16 40 8 8 48 4 24 | 88 = P ( 2) 4x2 4x + 24) |{z} 88 . {z } Q(x)

P ( 2)

(1) x0 is a root of P (x) or P (x0) = 0 if the remainder P (x0) = 0. (2) P 0(x) = Q(x) + (x

x0)Q0(x) =) P 0(x0) = Q(x0).

6. MULLER’S METHOD

x0 =

2

1

2 2 4 2 6

1 1

12 8 4 12 8

16 8 24 16 |8

41

40 48 | 88 = P ( 2) = P 0( 2)

Finding roots by deflation To find roots of P (x): (1) Find xˆ1 by using Newton. Then use Horner to get P (x) = (x

xˆ1)Q1(x) + P (ˆ x ) ⇡ (x | {z1}

xˆ1)Q1(x).

xˆ2)Q2(x) + P (ˆ x ) ⇡ (x | {z2}

xˆ2)Q2(x).

⇡0

(2) Find xˆ2 by using Newton on Q1(x) (deflation). Then use Horner to get Q1(x) = (x

⇡0

Continue until 2 roots remain. (3) P (x) = (x

xˆ1)(x

xˆ2) + · · · + (x

xˆk )Qk (x).

Use the quadratic formula to find the last two roots. (4) Since the approximations have introduced errors, use each xˆi as p0 for Newton with P (x) to improve accuracy. Note. In steps (1) and (2), begin with p0 = complex to find complex roots.

42

2. SOLUTIONS OF EQUATIONS OF ONE VARIABLE

The pseudocode for Horner’s method follows.

Maple. See horneralg.mw or horneralg.pdf and poly.mw or poly.pdf.

6. MULLER’S METHOD

43

Muller’s Method This is an extension of the secant method (and is not just for polynomials).

For the secant method (graph (a) above), we draw the secant line through the points (x0, f (x0)) and (x1, f (x1)) and take our new approximation to be the point where the secant meets the x-axis. For Muller’s method (graph (b) above), we find the parabola through the points (x0, f (x0)), (x1, f (x1)), and (x2, f (x2)) and then take the intercept closest to x2 as our new approximation. Suppose we have 3 initial approximations p0, p1, and p2 for a root p of f (x) = 0. To find p3, consider P (x) = a(x

p2)2 + b(x

p2) + c

that passes through the points (p0, f (p0)), (p1, f (p1)), and (p2, f (p2)). For p0, p1, and p2 real, the graph of P is a parabola. The constants a, b, c of P (x) can be found by substituting each of the points above into P (x). f (p0) = a(p0

p2)2 + b(p0

p2) + c

f (p1) = a(p1

p2)2 + b(p1

p2) + c

f (p2) = a · 02 + b · 0 + c

44

2. SOLUTIONS OF EQUATIONS OF ONE VARIABLE

We solve this system of equations simultaneously for a, b, and c to get c = f (p2) b= a=

(p0 (p1

⇥ ⇤ ⇥ ⇤ p2)2 f (p1) f (p2) (p1 p2)2 f (p0) f (p2) (p0 p2)(p1 p2)(p0 p1) ⇥ ⇤ ⇥ ⇤ 2 p2) f (p0) f (p2) (p0 p2) f (p1) f (p2) (p0 p2)(p1 p2)(p0 p1)

We now solve P (x) = 0 for p3. If P (x) is constant (a = b = 0), the method breaks down. If P (x) is linear, (a = 0, b 6= 0), p3 is the solution of P (x) = 0.

If P (x) is quadratic (a 6= 0), p3 is the solution closest to x2, 2c p p3 = p2 b ± b2 4ac where the sign in the denominator is chosen so as to give the denominator the largest magnitude, which in the real case is 2c p p3 = p2 . 2 b + sgn(b) b 4ac The method continues until we get an approximation of desired accuracy. Note. If b2 4ac < 0, and if we use complex arithmatic, we get complex approximations.

6. MULLER’S METHOD

45

The pseudocode for Muller’s method follows.

Maple. See mulleralg.mw or mulleralg.pdf for the algorithm. See muller.mw or muller.pdf for the implementation of the algorithm, and mapleroots.mw or mapleroots.pdf .