不动点迭代求一元函数f(x)的近似根 --北太天元学习11

通过北太天元学习1-10, 我们已经掌握了北太天元的很多编程知识，在这一节，
应用前面学习的技巧，我们讲求解形式为f（x）=0的方程,
另外学习匿名函数以及函数句柄的知识。

使得f(x) = 0成立的x称为f(x)的零点或者函数f(x)的根。当f是特殊的函数时，
我们曾经学习了如何手工求解，例如 ax = b 这样的一元一次方程, ax^2+bx+c=0
这样的一元二次方程。 简单的函数，我们手工通过有限步的初等运算
就找到了解析解。然而，许多方程并不那么容易手工求解, 例如，
x^7+x−1=0，或sin(x^3+1)=2x.  

此时，我们可以退而求其次，想办法找到这些函数的解的近似值, 称为数值解或者
近似解。例如，我们可以通过绘制 f(x) 的图像，然后观察 f(x) 与 y = 0 直线的
交点大致在哪个区间之内。 但是此时得到的数值解往往是非常粗糙的(或者说误差是
多么大)，我们这一节学习如何编写一个计算机程序来找到误差更小的数值解。

不动点迭代:

对于函数 f(x), 如果 x 使得 f(x) = x 成立，那么称 x 是函数 f 的不动点。
对于一些问题，我们可以通过不动点迭代来求解函数的根。
例: 求  x^3 -5x -2 = 0 ;

北太天元画上面函数图像的代码
x = -1: 0.1: 1;
fh = @(t) t.^3 - 5*t -2;  % 这里定义了一个匿名函数
plot(x, fh(x), 'LineWidth', 5, 'Color', 'Red')

从图中，我们可以观察到在 -0.6, -0.2 之间有一个根，而且我们可以计算
fh(0.1)  得 f在 -0.6 的函数值是 -0.78, 在 -0.2 处的函数值是 1.01

我们可以把上面的求 f(x) = x^3 - 5*x -2 转成一个求某个函数的不动点的问题:
  x^3 -5x -2 = 0
  5x = x^3 - 2
    x = (x^3 - 2 )/5;
于是构造不动点迭代
  x_{n+1} = ( x_{n}^3- 2)/ 5;
可以取 x_{0} =  -0.6, 然后用上面的递推关系式计算 x_{n}, n = 1, 2, 3, ...
如果数列 { x_{n} } 极限存在，\lim_{n\rightarrow \infty}x_{n} = xi
那么
   xi = ( xi^3 - 2 ) /5,
xi 就是 g(x) = (x^3 - 2)/5 的不动点，也是 f(x) = x^3 - 5x-2 的一个根。

值得花一点时间来考虑如何得出最初的“猜测”。猜测初始值可能会导致不同的结果，
特别是如果方程 g(x) −x=0有多个根。有时需要一点试错才能确定一个可行的初步猜测。
我们可能知道关于函数的一些东西，或者能够使用数学理论来帮助缩小其中可能出现不动点。

迭代的停机准则(stopping critrion)的选择也是很多选择，我们可以选择两步迭代之间
的差小于某个给定的忍量(tolerance), 例如
| x_{n+1} - x_{n} | < 1e-4.

这里给出了一个求不动点的基本北太天元脚本:

x= -0.6；%初始化x
gh = @(t) (t.^3 -2)/5; %用匿名函数给出函数定义
tol=1e-4；%设定忍量
n=1；%设置索引

x(n+1) = gh( x(n) );

while abs(x(n+1)-x(n))> 1e-4；%应用不动点迭代, 知道两步迭代之间的更新小于给定忍量

    n = n + 1; ；%更新索引
     x(n+1) = gh( x(n) ); %应用不动点迭代, 得到下一个点
end
x(end) %输出迭代的最后一个x点
%画出迭代过程
hold on
plot(0:length(x)-1,x, 'r-', 'LineWidth', 3);
sh = scatter(0:length(x)-1, x, 'filled');
set(sh, 'SizeData', 500)
title(['不动点迭代 x(end) = ', num2str(x(end))]);
hold off