Читать книгу Fundamentals of Numerical Mathematics for Physicists and Engineers - Alvaro Meseguer - Страница 15
1.2.1 The Bisection Method
ОглавлениеIn order to clarify the concept of tolerance we introduce here the bisection method (also called the interval halving method). Let us assume that is a continuous function within the interval within which the function has a single root . Since the function is continuous, . First, we define , , and as the starting interval whose midpoint is . The main goal of the bisection method consists in identifying that half of in which the change of sign of actually takes place. Use the simple rules:
Bisection Rules:
Finally set and .
Figure 1.1 (a) A simple exploration of reveals a change of sign within the interval (the root has been depicted with a gray bullet). (b) We start the bisection considering the initial interval . The bisection rule provides the new halved interval .
Figure 1.1 shows the result of applying the bisection rule to find roots of the cubic polynomial already studied in Eq. (1.2). As shown in Figure 1.1a, this polynomial takes values of opposite sign, and , at the endpoints of the interval whose midpoint is . Since , the new (halved) interval provided by the bisection rules is . The midpoint of (white triangle in Figure 1.1b) provides a new estimation of the root of the polynomial. We could resume the process by checking whether or in order to obtain the new interval containing the root and so on. After bisections we have determined the interval with midpoint . In general, the interval will be obtained by applying the same rules to the previous interval :
Bisection Method: Given such that , compute and set
This general rule that provides from is an example of what is generally termed as algorithm,5 i.e. a set of mathematical calculations that sometimes involves decision‐making. Since this algorithm must be repeated or iterated, the bisection method described above constitutes an example of what is also termed as an iterative algorithm.
Now we revisit the concept of tolerance seen from the point of view of the bisection process. For , the estimation of the root was , which is the midpoint of the interval , whereas for we obtain a narrower region of existence of such root, as well as its corresponding improved estimation . In other words, for the root lies within , with a tolerance , whereas for the interval containing the root is , with . Overall, after bisections, the tolerance is , becoming halved in the next iteration.
If we evaluate the radicals that appear in (1.3) with a scientific calculator, we can check that Cardano's analytical solution is approximately . Expression (1.3) is mathematically elegant, but not very practical if one needs an estimation of its numerical value. For we already have that estimation nearly within a or relative error. The reader may keep iterating further to provide better approximations of , overall obtaining the sequence . A natural question is whether this sequence has a limit. This leads to the mathematical concept of convergence of a sequence:
Convergence (Exact): The sequence is said to be convergent to if or, equivalently, if , where .
The difference appearing in the previous definition is called the error associated with . However, since may be positive or negative and we are just interested in the absolute discrepancy between and the root , it is common practice to define the absolute error of the th iterate. Checking convergence numerically using the previous definition has two main drawbacks. The first one is that we do not know the value of , which is precisely the goal of root‐finding. The second is that in practice we cannot perform an infinite number of iterations in order to compute the limit as , and therefore we must substitute the convergence condition by a numerically feasible one. For example, looking at the sequence resulting from the bisection method, it is clear that the absolute difference of two consecutive elements decreases when increasing (for example, and ). Since, by construction, the bisection sequence satisfies , it is common practice to consider a sequence as converged when this difference becomes smaller than a prescribed threshold or tolerance :
Convergence (Practical): A sequence is said to be converged to the desired tolerance after iterations if , .
The criterion above only refers to the convergence of a sequence, and not necessarily to the convergence to a root. After an iterative method has apparently converged to some stagnated value , it is always advisable to check the magnitude of in order to confirm if the method has succeeded. The code below is a simple implementation of the bisection method using the tolerance criterion previously described:
% Code 1: Bisection method for solving f(x) = 0 in [a,b] % Input: 1. [a,b]: interval (it assumes that f(a)f(b) < 0) % 2. tol: tolerance so that abs(x_k+1 - x_k) < tol % 3. itmax: maximum number of iterations allowed % 4. fun: function's name % Output: 1. xk: resulting sequence % 2. res: resulting residuals % 3. it: number of required iterations function [xk,res,it] = bisection(a,b,tol,itmax,fun) ak = a; bk = b; xk = []; res = []; it = 0; tolk = abs(bk-ak)/2 ; while it < itmax & tolk > tol ck = (ak + bk)/2; xk = [xk ck]; if it > 0; tolk = abs(xk(end)-xk(end-1)); end fa = feval(fun,ak); fc = feval(fun,ck); res = [res abs(fc)]; if fc*fa < 0; bk = ck; else ak = ck; end it = it + 1 ; end
In the previous code, the iterations stop either when the number of iterations reaches itmax
or when tol
, i.e. the condition for practical convergence. This condition therefore provides what is usually known as a stopping criterion for root‐finding methods such as the bisection and other algorithms. However, two words of caution deserve to be mentioned at this point. First, we are implicitly assuming that since seems to be a Cauchy Sequence6 (that is, decreases with increasing ), its convergence is guaranteed. While a convergent sequence must necessarily be of Cauchy type, the reverse statement is, in general, not true (although here it will be assumed to be). Second, the quantity has disappeared from the convergence criteria. In other words, once for some and beyond, we only know that the sequence has converged (in the practical sense) to some value , that henceforth will play the role of the numerical root within the prescribed tolerance.
Finally, the bisection code also provides the vector res
containing the sequence of absolute quantities , usually called the residuals. Since this quantity should approach zero when converges to a root, the reader may wonder why should not be used as a measure for stopping the iterations. The reason is that the residual is not always a reliable measure (although it is often used in practice to decide convergence). We will address this question later in this chapter, when we introduce the concept of condition number of a root.