二分法和if … else … end 语句

先回顾一下二分法。要求方程\[f(x)=0\]的根。假设$c = f(a) < 0$和$d = f(b) > 0$，如果$f(x)$是连续函数，那么方程的根$x^*$一定位于$a$和$b$之间。然后，我们看一下$a$和$b$中点$x=(a+b)/2$，计算函数值$y=f(x)$，如果函数不为0，比较$c$、$d$和$y$的符号，确定新的二分区间。具体来说，如果$c$和$y$同号，新的二分区间就是$[x,b]$，如果$c$和$y$异号，新的二分区间就是$[a,x]$。如图5.1所示。

图5.1 二分法

在不同的情况下做不同的事情，在程序里，这叫做流程控制。达到这个目的，最通常的方法是通过if … else … end语句来实现，这个语句是我们前面讲过的if … end语句的扩展。

误差界

二分法的一个好处是，我们一直知道方程的真实的根$x^*$一定位于二分区间内，这样我们就可以知道最大误差可以是多少，一定小于区间大小的一半，即

\begin{equation*}
\mid x-x^* \mid \lt \frac{b-a}{2}
\end{equation*}

其中$x=(a+b)/2$，为区间的中点。

二分法程序如下：

function [x e] = mybisect (f,a,b,n)
% function [x e] = mybisect (f,a,b,n)
% Does n iterations of the bisection method for a function f
% Inputs : f -- an inline function
% a,b -- left and right edges of the interval
% n -- the number of bisections to do.
% Outputs : x -- the estimated solution of f(x) = 0
% e -- an upper bound on the error
format long
% evaluate at the ends and make sure there is a sign change
c = f(a); d = f(b);
if c*d > 0.0
error (’Function has same sign at both endpoints .’)
end
disp (’ x y’)
for i = 1:n
% find the middle and evaluate there
x = (a + b )/2;
y = f(x);
disp ([ x y])
if y == 0.0 % solved the equation exactly
e = 0;
break % jumps out of the for loop
end
% decide which half to keep , so that the signs at the ends differ
if c*y < 0
b=x;
else
a=x;
end
end
% set the best estimate for x and the error bound
x = (a + b )/2;
e = (b-a )/2;

此程序不仅给出了方程的近似根，还给出了最大可能误差。

二分法总是可以进行下去，而牛顿法，如果初值$x_0$不足够接近$x^*$的话，程序可能不会收敛。而二分法，每一步的区间都可以缩一半，最终想缩多小就能缩多小。

找根

二分法和牛顿法都可以用来找根的任意近似值，但都需要先给定初值。二分法需要两个初值$a$和$b$，并要求真实的根位于二者之间，牛顿法只需要一个初值$x_0$，要求不能离真实的根太远。如何确定合适的初值？看情况而定。如果你一次只解一个方程，最好先把函数画出来，通过图上的函数零点，可以很容易定出合适的初值。

有些情况下你不是一次就解一个方程，而是多次求解同一个方程，只是系数不同。当你开发某特定用途的软件的时候，常常会遇到这种情况。在这种情况下，首先你要利用问题的自然条件，比如选定的区间要符号实际情况。然后通过区间内某些点的符号，很容易就可以逐渐接近方程的真实的根了。方程的根一定位于符号相异的两点之间。下面的程序就是在一给定的区间$[a_0,b_0]$之间找方程的根：

function [a,b] = myrootfind (f,a0 ,b0)
% function [a,b] = myrootfind (f,a0 ,b0)
% Looks for subintervals where the function changes sign
% Inputs : f -- an inline function
% a0 -- the left edge of the domain
% b0 -- the right edge of the domain
% Outputs : a -- an array , giving the left edges of subintervals
% on which f changes sign
% b -- an array , giving the right edges of the subintervals
n = 1001; % number of test points to use
a = []; % start empty array
b = [];
% split the interval into n -1 intervals and evaluate at the break points
x = linspace (a0 ,b0 ,n);
y = f(x);
% loop through the intervals
for i = 1:(n -1)
if y(i)*y(i +1) < 0 % The sign changed , record it
a = [a x(i)];
b = [b x(i +1)];
end
end
if size (a ,1) == 0
warning (’no roots were found ’)
end

如果没有任何已知信息来找方程的根，事情就比较难办，好在这种情况在工程问题里不常见。

一旦确定出根所在的区间$[a, b]$，$a$和$b$就可以作为二分法和割线法（见下一讲）的初值。对于牛顿法，初值$x_0$一个明显的选法是$x_0 = (a + b)/2$。一个更好的选法是利用割线法选$x_0$。

练习

5.1 修改程序mybisect，在给定公差下解方程。应用你的程序，求解方程$f(x) = 2x^3 + 3x − 1=0$的根，初始区间是$[0,1]$，公差为$10^{-8}$。程序需要运行多少步能达到公差的要求？最终残量有多大？

5.2 用纸和计算器对于函数$f(x) = x^3 − 4$和初始区间$[1,3]$进行3步二分法迭代，计算每步结果的误差和相对误差，并与第3讲练习3的结果比较。