提供了编程的基础技术教程

zl程序教程

您现在的位置是：首页 >  其它

当前栏目

非线性方程(一)

方程 非线性
2023-06-13 09:12:58 时间

大家好，又见面了，我是你们的朋友全栈君。

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

1、diff

>> diff('x^2')

ans =

   -26   -44


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

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

2、牛顿法

对  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为公式。

3、误差初探

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

考虑后一种算法。

example：

　　首先将其写成拟牛顿法的形式，xk+1 = xk – f(xk)/sk , 其中 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

发布者：全栈程序员栈长，转载请注明出处：https://javaforall.cn/167251.html原文链接：https://javaforall.cn

相关文章