非线性方程(一)

此为全书第一章，主matlab入门——经过学习各类插值法：反线性插值、牛顿法之类。

 

一、diff

>> diff('x^2')

ans =

   -26   -44


>> diff(x^2)
 
ans =
 
2*x

　书上用的是加引号那个版本的，以致于我困了很久。

 

二、牛顿法

对  f(x) = (x*exp(x))/2 - 2*x^2 进行牛顿插值

解

>> x = 0: 0.05: 3;
>> y = 0.5*x.*exp(x) - 2*x.^2;
>> plot(x,y);grid

 

牛顿插值实现

1 >> f =inline('(x*exp(x))/2 - 2*x^2');
2 >> fp = inline('exp(x)/2 - 4*x + (x*exp(x))/2');  %f的导数
3 >> x0 = 0.4;
4 >> x0 = x0 - f(x0)/fp(x0)
5 
6 x0 =
7 
8     0.3611

其中line 4为公式。

 

三、偏差初探

关于有限精度的做用能够用解析方法进行研究，也是必须讨论的内容，但这里只研究其对割线法和有限差分牛顿法的影响。

考虑后一种算法。

example：

　　首先将其写成拟牛顿法的形式，x_k+1 = x_k - f(x_k)/s_k , 其中  _sk  = (f(x _k +h) - f(x   _k ))/h；

如今用牛顿法近似f(x)在 x=1的导数。

>> x=1; h=0.1;
>> (log(x+h) - log(x))/h

ans =

    0.9531

对于0.1的步长，这个结果还算合理【书上原话= =

>> h=10.^(-(1:20))
>> (log(x+h) - log(x))./h

ans =

  Columns 1 through 5

    0.9531    0.9950    0.9995    1.0000    1.0000

  Columns 6 through 10

    1.0000    1.0000    1.0000    1.0000    1.0000

  Columns 11 through 15

    1.0000    1.0001    0.9992    0.9992    1.1102

  Columns 16 through 20

         0         0         0         0         0

>> plot(h, abs(ans-1), 'r'), grid

 

有图能够看到，偏差有一个时期是降低的，从大约h = 1E-12开始，偏差又开始增长。

这里主要指出，当h很小时，上面的计算中有两个因素在做怪：一个是分子中两个接近相等的数之间的减法，另外一个是用很小的数做除数。

产生这些缘由在于计算机只能表示数字的有限多位。

附几段关于计算机精度代码

>> inf*2
ans =    Inf

>> inf -inf
ans =    NaN

>> eps
ans =    2.2204e-16


>> 1+eps
ans =    1.0000

>> 1+eps == 1
ans =    0



>> 1+eps/2
ans =      1

>> ans == 1
ans =     1



>> (1+eps/2) + eps/2
ans =    1

>> ans == 1
ans =    1



>> 1 + (eps/2 + eps/2)  %计算机加法不知足结合律
ans =    1.0000

>> ans == 1
ans =     0

>> x1= 0.1; x2= pi/2;
>> abs(sin(x1+err)- sin(x1))

ans =  9.9500e-11

>> abs(sin(x2+err)- sin(x2))

ans =  0

%在x接近(pi/2)【此时导数为0】时sin(x)的求值要比x接近于零时的求值更加精确，由于零点的导数为1.
%能够看到f的导数(变化速度)越大，偏差也越大

>> err = 1e-6

err =    1.0000e-06

>> abs(sin(x1+err)- sin(x1))

ans =    9.9500e-07

>> abs(sin(x2+err)- sin(x2))

ans =    5.0004e-13