1 of 24

UNIT-1

ROOT FINDING

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

1

2 of 24

Computational METHOD

3RD SEMESTER

ES-201

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

2

Lecture Notes

3 of 24

1-D Root Finding

  • Given some function, find location where f(x)=0
  • Need:
    • Starting position x0, hopefully close to solution
    • Ideally, points that bracket the root

f(x+) > 0

f(x–) < 0

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

3

4 of 24

1-D Root Finding

  • Given some function, find location where f(x)=0
  • Need:
    • Starting position x0, hopefully close to solution
    • Ideally, points that bracket the root
    • Well-behaved function

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

4

5 of 24

What Goes Wrong?

Tangent point:

very difficult�to find

Singularity:

brackets don’t�surround root

Pathological case:

infinite number of�roots – e.g. sin(1/x)

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

5

6 of 24

Example: Press et al., Numerical Recipes in C�Equation (3.0.1), p. 105

2

2 ln((Pi - x) )

f := x -> 3 x + ------------- + 1

4

Pi

> evalf(f(Pi-1e-10));

30.1360472472915692...

​

> evalf(f(Pi-1e-100));

25.8811536623525653...

​

> evalf(f(Pi-1e-600));

2.24285595777501258...

​

> evalf(f(Pi-1e-700));

-2.4848035831404979...

​

> evalf(f(Pi+1e-700));

-2.4848035831404979...

​

> evalf(f(Pi+1e-600));

2.24285595777501258...

​

> evalf(f(Pi+1e-100));

25.8811536623525653...

​

> evalf(f(Pi+1e-10));

30.1360472510614803...

​

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

6

7 of 24

Bisection Method

  • Given points x+ and x– that bracket a root, find� xhalf = ½ (x++ x–)�and evaluate f(xhalf)
  • If positive, x+ ← xhalf else x– ← xhalf
  • Stop when x+ and x– close enough
  • If function is continuous, this will succeed�in finding some root

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

7

8 of 24

Bisection

  • Very robust method
  • Convergence rate:
    • Error bounded by size of [x+… x–] interval
    • Interval shrinks in half at each iteration
    • Therefore, error cut in half at each iteration:� |εn+1| = ½ |εn|
    • This is called “linear convergence”
    • One extra bit of accuracy in x at each iteration

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

8

9 of 24

Faster Root-Finding

  • Fancier methods get super-linear convergence
    • Typical approach: model function locally by something whose root you can find exactly
    • Model didn’t match function exactly, so iterate
    • In many cases, these are less safe than bisection

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

9

10 of 24

Secant Method

  • Simple extension to bisection: interpolate or extrapolate through two most recent points

1

2

3

4

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

10

11 of 24

Secant Method

  • Faster than bisection:� |εn+1| = const. |εn|1.6
  • Faster than linear: number of correct bits multiplied by 1.6
  • Drawback: the above only true if sufficiently close to a root of a sufficiently smooth function
    • Does not guarantee that root remains bracketed

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

11

12 of 24

False Position Method

  • Similar to secant, but guarantee bracketing

​

​

​

​

​

  • Stable, but linear in bad cases

1

2

3

4

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

12

13 of 24

Other Interpolation Strategies

  • Ridders’s method: fit exponential to�f(x+), f(x–), and f(xhalf)
  • Van Wijngaarden-Dekker-Brent method:�inverse quadratic fit to 3 most recent points�if within bracket, else bisection
  • Both of these safe if function is nasty, but�fast (super-linear) if function is nice

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

13

14 of 24

Newton-Raphson

  • Best-known algorithm for getting quadratic�convergence when derivative is easy to evaluate
  • Another variant on the extrapolation theme

1

2

3

4

Slope = derivative at 1

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

14

15 of 24

Newton-Raphson

  • Begin with Taylor series

​

  • Divide by derivative (can’t be zero!)

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

15

16 of 24

Newton-Raphson

  • Method fragile: can easily get confused

​

​

​

​

  • Good starting point critical
    • Newton popular for “polishing off” a root found approximately using a more robust method

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

16

17 of 24

Newton-Raphson Convergence

  • Can talk about “basin of convergence”:�range of x0 for which method finds a root
  • Can be extremely complex:�here’s an example�in 2-D with 4 roots

Yale site

D.W. Hyatt

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

17

18 of 24

Popular Example of Newton: Square Root

  • Let f(x) = x2 – a: zero of this is square root of a
  • f’(x) = 2x, so Newton iteration is

​

​

  • “Divide and average” method

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

18

19 of 24

Reciprocal via Newton

  • Division is slowest of basic operations
  • On some computers, hardware divide not available (!): simulate in software

​

​

​

​

  • Need only subtract and multiply

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

19

20 of 24

Rootfinding in >1D

  • Behavior can be complex: e.g. in 2D

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

20

21 of 24

Rootfinding in >1D

  • Can’t bracket and bisect

​

​

​

​

​

​

  • Result: few general methods

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

21

22 of 24

Newton in Higher Dimensions

  • Start with

​

​

  • Write as vector-valued function

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

22

23 of 24

Newton in Higher Dimensions

  • Expand in terms of Taylor series

​

​

  • f’ is a Jacobian

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

23

24 of 24

Newton in Higher Dimensions

  • Solve for δ

​

​

  • Requires matrix inversion (we’ll see this later)
  • Often fragile, must be careful
    • Keep track of whether error decreases
    • If not, try a smaller step in direction δ

DEPARTMENT OF APPLIED SCIENCES, BVCOE NEW DELHI

24