现代数值计算（第3版）课件 第6、7章 数值积分与数值微分、非线性方程求解

第六章AdvancedNumericalComputing数值积分与数值微分现代数值计算同济大学数学科学学院目录/Contents第六章数值积分与数值微分第一节几个常用积分公式及其复合公式第二节变步长方法与外推加速技术第三节高斯公式第四节多重积分的计算第五节数值微分引言Newton-Leibniz公式:数值(近似)积分的必要性：无解析原函数的被积函数：

，

，等等利用数据表给出的函数的积分计算机程序提供的函数的积分被积函数表达式太复杂设在 (不妨先设

为有限数)上，

，

为某个较“简单”的函数,则有误差为:因此只要

，就有误差估计引言6.1几个常用积分公式及其复合公式6.1.1几个常用积分公式6.1几个常用积分公式及其复合公式6.1.1几个常用积分公式6.1几个常用积分公式及其复合公式6.1.1几个常用积分公式6.1几个常用积分公式及其复合公式6.1.1几个常用积分公式6.1几个常用积分公式及其复合公式6.1.1几个常用积分公式6.1几个常用积分公式及其复合公式6.1.1几个常用积分公式设如果对函数

，用

上的函数值近似代替，即得中点公式推导误差：设

，由Taylor公式,两边积分，即得6.1几个常用积分公式及其复合公式6.1.1几个常用积分公式下面我们通过插值节点

和

作线性插值函数

，利用

得梯形公式：上面求积公式的右端值可看成是由线段

，过点

的直线以及

轴围成的梯形面积.如果

，则由线性插值函数的误差公式(见第3章)以及积分中值定理得：6.1几个常用积分公式及其复合公式6.1.1几个常用积分公式6.1几个常用积分公式及其复合公式6.1.1几个常用积分公式若

用通过节点

的二次插值多项式

代替，6.1几个常用积分公式及其复合公式6.1.1几个常用积分公式则得抛物型公式(或称Simpson公式)可以证明：若

，则有误差公式:6.1几个常用积分公式及其复合公式6.1.1几个常用积分公式求积公式

、节点

、系数

、余项

.6.1几个常用积分公式及其复合公式6.1.1几个常用积分公式【定义6.1】如果对于所有次数

的多项式

，等式

精确成立，但对于某一次数为

的多项式不精确成立，则称此求积公式的代数精度为

次。6.1几个常用积分公式及其复合公式6.1.2代数精度【例6.1】试确定系数

，

，使得求积公式

有尽可能高的代数精度，并求出此求积公式的代数精度。解：分别令

, , ,代入使积分公式精确成立，得到线性方程组6.1几个常用积分公式及其复合公式6.1.2代数精度解得

，

这样求积公式为该公式对

精确成立,但

时不精确成立。因此具有三次代数精度。6.1几个常用积分公式及其复合公式6.1.2代数精度设在节点上的函数值为作

次Lagrange插值多项式，有插值型求积公式代数精度有几次?6.1几个常用积分公式及其复合公式【定理6.3】求积公式

至少具有

次代数精度的充分必要条件是它是插值型的。

次数

时，

，代数精度至少

次。反之，若代数精度至少

次，则必定是插值型的：用

代入,6.1几个常用积分公式及其复合公式Newton-Cotes公式即是等距节点的插值型求积公式。将区间

作

等分,步长

，等距节点

，

.

作变量代换

，代入

中，有(NC系数)6.1几个常用积分公式及其复合公式常见的Newton-Cotes公式梯形公式(一次)抛物线(Simpson)公式(三次)Cotes公式(五次)6.1几个常用积分公式及其复合公式有负系数6.1几个常用积分公式及其复合公式

n12345678【定理6.4】当

为偶数时，牛顿-柯特斯公式的代数精度至少为 .证明:只需要证明

时，积分余项为 .作变量替换

.对上述积分再做变量代换

，得6.1几个常用积分公式及其复合公式稳定性即数据

有误差，只能有近似值

，是否有若

，则当

时，

，且

稳定.而当

时，由于求积系数

有正有负，

一般无界，因此稳定性不成立，收敛性也不成立.6.1几个常用积分公式及其复合公式设

，由误差公式知当

很小时，中点公式，梯形公式及抛物型公式的误差为

，

，

通常情况下，积分区间

的长度

不是非常小，此时误差就比较大。为了确保计算精度，提出了复合求积公式。即将积分区间分为若干份，在每一个“小区间”上用低阶求积公式如中点公式，梯形公式及抛物型公式进行计算，再将计算值相加即得原积分的近似值.6.1几个常用积分公式及其复合公式6.1.2代数精度复合中点公式等分求积区间,记小区间

中点为

，记步长

，每个小区间上运用中点求积公式设

，则(待证)存在

6.1几个常用积分公式及其复合公式6.1.3积分公式的复合functionI=fmid(fun,a,b,n)

h=(b-a)/n;

x=linspace(a+h/2,b-h/2,n);

y=feval(fun,x);

I=h*sum(y);6.1几个常用积分公式及其复合公式6.1.3积分公式的复合复合梯形公式等分求积区间，节点为，每个小区间上运用梯形求积公式设

，则存在

及

6.1几个常用积分公式及其复合公式6.1.3积分公式的复合functionI=ftrapz(fun,a,b,n)

h=(b-a)/n;

x=linspace(a,b,n+1);

y=feval(fun,x);

I=h*(0.5*y(1)+sum(y(2:n))+0.5*y(n+1));6.1几个常用积分公式及其复合公式6.1.3积分公式的复合复合Simpson公式等分求积区间，记小区间

中点为

，在每个小区间上Simpson求积公式,设

，则6.1几个常用积分公式及其复合公式6.1.3积分公式的复合functionI=fsimpson(fun,a,b,n)

h=(b-a)/n;

x=linspace(a,b,2*n+1);

y=feval(fun,x);

I=(h/6)*(y(1)+2*sum(y(3:2:2*n-1))+4*sum(y(2:2:2*n))+y(2*n+1));6.1几个常用积分公式及其复合公式6.1.3积分公式的复合【定理6.1】设

，则存在

，使得证：设

，

其中

.

则有6.1几个常用积分公式及其复合公式6.1.3积分公式的复合利用连续函数的中值定理，存在

满足：所以6.1几个常用积分公式及其复合公式6.1.3积分公式的复合下面设，推导复合中点公式(5.13)的误差.首先由(5.12),(5.13)及误差公式(5.3),利用定理5.2.3，便有6.1几个常用积分公式及其复合公式6.1.3积分公式的复合同理可推得复合梯形公式、复合Simpson公式的误差估计：6.1几个常用积分公式及其复合公式6.1.3积分公式的复合【例6.2】设

在9个节点处的值由下表给出，分别用复合梯形公式和复合Simpson公式计算积分6.1几个常用积分公式及其复合公式6.1.3积分公式的复合01/81/43/81/210.99739780.98961580.97672670.95385105/83/47/810.93615560.90885160.87719250.8414709解：将积分区间

八等分，由复合梯形公式(此时 )计算得:如果将

四等分，由复合Simpson公式( )得：6.1几个常用积分公式及其复合公式6.1.3积分公式的复合functiony=f(x)

x=x+(x==0)*eps;

y=sin(x)./x;

调用:ftrapz(@f,0,1,8)或fsimpson(@f,0,1,4)与精确值I=0.9460831相比较：复合梯形公式两位有效数字,复合Simpson公式有六位有效数字!6.1几个常用积分公式及其复合公式6.1.3积分公式的复合【例6.3】利用n=5的复合Simpson公式计算积分将

五等分，由复合Simpson公式( )得:6.1几个常用积分公式及其复合公式6.1.3积分公式的复合【例6.4】使用三种不同的计算公式：复合中点公式,复合梯形公式,复合Simpson公式计算下列积分值,并分别取

，比较三种不同算法的收敛速度.积分的精确值为6.1几个常用积分公式及其复合公式6.1.3积分公式的复合解：三种复合求积公式的计算结果分别用

表示，则计算结果如下表.6.1几个常用积分公式及其复合公式6.1.3积分公式的复合

n10.97510.15890.703021.03700.56700.502140.12220.23483.139×10-382.980×10-25.635×10-21.085×10-3166.748×10-31.327×10-27.381×10-5321.639×10-33.263×10-34.682×10-6644.066×10-48.123×10-42.936×10-71281.014×10-42.028×10-41.836×10-82562.535×10-55.070×10-51.148×10-9三者的收敛速度比较图,其中横坐标为

，纵坐标为

，从图中可看出三种计算公式的误差分别为:与理论误差估计相当吻合.6.1几个常用积分公式及其复合公式6.1.3积分公式的复合对于节点

及函数值

，定义三个分段函数

则在

有:6.1几个常用积分公式及其复合公式6.1.3积分公式的复合则三个函数都可看成是

的某种逼近6.1几个常用积分公式及其复合公式6.1.3积分公式的复合复合中点公式,复合梯形公式和复合Simpson公式可看成

用代替后的近似积分值易证当

为连续函数时， .从而有6.1几个常用积分公式及其复合公式6.1.3积分公式的复合复合求积公式的稳定性:设

在节点

处的精确值为

，实际值为

.

节点

处的误差为：记由数值

及

计算所得的梯形公式的值为

和

，则设

，则6.1几个常用积分公式及其复合公式6.1.3积分公式的复合目录/Contents第六章数值积分与数值微分第一节几个常用积分公式及其复合公式第二节变步长方法与外推加速技术第三节高斯公式第四节多重积分的计算第五节数值微分引言：(1)复合求积公式是有效的求积方法，当步长

越小,计算精度越高.(2)实际运用中，需选取一个合适的步长

：若步长太大，计算精度就难以保证，若步长太小，

则会增加不必要的计算开销；(3)在给定计算精度的情形下，往往通过不断调整步长的方式进行计算：如采用让步长逐次折半的方式，反复使用复合求积公式直至相邻两次计算结果之差的绝对值小于给定的计算精度为止。这种方法即称为变步长算法。6.2变步长方法与外推加速技术6.2.1变步长梯形法下面以变步长的梯形公式加以说明：由求积公式误差估计，若

，则有6.2变步长方法与外推加速技术6.2.1变步长梯形法只要

和

充分接近,就能保证

的误差很小。算法(区间折半法)：取初始步长h=b-a；

计算

；

取步长

，计算出相应的积分值

；若条件

满足，则取

为最后积分计算的近似值,否则让步长折半，回到第(3)步.6.2变步长方法与外推加速技术6.2.1变步长梯形法变步长梯形积分区间折半,新增区间中点，

6.2变步长方法与外推加速技术6.2.1变步长梯形法【例6.5】用变步长梯形公式计算积分

的近似值,要求计

算精度

将积分区间等分

份复合梯形公式才满足计算精度。6.2变步长方法与外推加速技术6.2.1变步长梯形法kk00.920735560.946076910.939793370.946081520.944513580.946082730.945690990.946083040.9459850100.946083150.9460596

例如上例中，

， (各具有2,3个有效数字) 新的近似值具有6位有效数字。实际上，6.2变步长方法与外推加速技术6.2.2Romberg算法例如，加上修正部分使精度从梯形的

提高到了抛物形的 .6.2变步长方法与外推加速技术6.2.2Romberg算法实际上，从抛物形积分组合得到Cotes积分：从Cotes积分组合得到Romberg积分：6.2变步长方法与外推加速技术6.2.2Romberg算法【例6.6】用Romberg算法求积分

，要求精度

.解：按照公式及不同的步长(初始步长 )，计算结果如下表：注：

表示积分区间

的等分次数，节点数为

.Romberg方法最后采用的步长为

，而上例中普通变步长法为

.两者的计算精度大致相同,因此Romberg加速收敛的效果是非常明显的.6.2变步长方法与外推加速技术6.2.2Romberg算法k

00.920735510.93979330.946145920.94451350.94608690.946083030.94569090.94608340.94608310.9460831Richardson外推加速收敛技术，其实质是利用不同步长的复合梯形公式以及外推技术加速收敛速度，通常也称为龙贝格Romberg求积方法，为推导更一般的结论，我们给出【定理6.2】欧拉-麦克劳林(Euler-MacLaurin)公式:设 (k为非负整数)，

，则存在伯努利数

，使得下面关系成立：6.2变步长方法与外推加速技术6.2.2外推加速技术与Romberg算法由上面定理知复合梯形公式的误差为两边同乘以

相减，有以

代替

6.2变步长方法与外推加速技术6.2.2Romberg算法由此可见，

的计算精度为

，而

的计算精度为

一般记

，有这个计算过程即称为龙贝格（Romberg）方法.的计算精度为

6.2变步长方法与外推加速技术6.2.2Romberg算法目录/Contents第六章数值积分与数值微分第一节几个常用积分公式及其复合公式第二节变步长方法与外推加速技术第三节高斯公式第四节多重积分的计算第五节数值微分若上面求积公式对任何次数小于等于

的多项式

等式成立，但对于某一个次数为

的多项式不成立，则称积分公式的代数精度为

次。6.4高斯(Gauss)公式6.4.1高斯公式的定义及性质考虑带权函数的求积公式其中

称为权函数.显然

时上述定义与第二节一致，因此可以将这里的代数精度的概念其概念的推广.

为一些特殊的函数。如

，

等.而

一般比较光滑的函数。问：对于等分节点的Newton-Cotes公式，代数精度一般为

或

，而如果采用非等分节点，即节点

与求积系数

都待定的话，能否提高代数精度呢?要使求积公式具有

次代数精度，应对

，，，成立下列等式(为叙述简单，以

为例)：关于未知量

与

的非线性代数方程组,如何求解？6.4高斯(Gauss)公式6.4.1高斯公式的定义及性质如何选择

和

使代数精度达到最高?6.4高斯(Gauss)公式6.4.1高斯公式的定义及性质最简单的例子

：试确定

和

使下面的求积公式有尽量高的代数精度假设

，

.令

, , ,,代入因此，6.4高斯(Gauss)公式6.4.1高斯公式的定义及性质【定义6.2】若对于节点

及求积系数

，求积公式(5.26)的代数精度为

则称节点

为高斯点，

为高斯系数，相应的求积公式称为带权的高斯公式。一旦确定了高斯点，由于它的代数精度为

，因此仍然可用待定系数法的方法来确定系数

特别地，若在高斯公式中令

，则有因此高斯公式仍可看成是一种插值型求积公式。6.4高斯(Gauss)公式6.4.1高斯公式的定义及性质【定理6.5】

是求积公式的高斯点的充分必要条件是多项式

与任意次数不超过

的多项式

关于权函数

正交：6.4高斯(Gauss)公式6.4.1高斯公式的定义及性质必要性：因为

次数不超过

，其中

，

次数不超过 .充分性：对于任意给定的次数不超过

次的多项式

，令

.

6.4高斯(Gauss)公式6.4.1高斯公式的定义及性质6.4高斯(Gauss)公式6.4.1高斯公式的定义及性质【定理6.6】插值型求积公式的代数精度最高不超过 .证：令

，则

是

次的多项式。该积分公式的代数精度不可能达到 .6.4高斯(Gauss)公式6.4.1高斯公式的定义及性质求积公式的稳定性令

为

次多项式.由于Gauss积分精确成立，所以若在高斯公式中取

，则有利用

的性质及稳定性的讨论知高斯公式对任意

都是稳定的.6.4高斯(Gauss)公式6.4.1高斯公式的定义及性质【例6.7】求高斯型求积公式的系数

和节点

。6.4高斯(Gauss)公式6.4.1高斯公式的定义及性质【定理6.7】设

为有限数，则对任意的函数

，当

时，高斯求积公式均收敛，即高斯型公式有三次代数精度，令 , , ,,代入因此，6.4高斯(Gauss)公式6.4.1高斯公式的定义及性质解之得几个常见的Gauss-Legendre公式1.高斯--勒让德(Gauss-Legendre)求积公式

取

，知高斯点

即为

次勒让德正交多项式的零点，

.6.4高斯(Gauss)公式6.4.2常用Gauss求积公式6.4高斯(Gauss)公式6.4.2常用Gauss求积公式然后再用相应的高斯--勒让德公式来计算对于一般区间

上的积分,我们可以做变量代换将积分区间化为

，其中

即为高斯点，

为高斯系数.6.4高斯(Gauss)公式6.4.2常用Gauss求积公式2.高斯--切比雪夫(Gauss-Chebyshev)公式取 .知此时的高斯点为区间

关于权函数

的正交多项式--切比雪夫正交多项式

的零点.高斯点高斯系数6.4高斯(Gauss)公式6.4.2常用Gauss求积公式3.高斯--拉盖尔(Gauss-Laguerre)求积公式：取

，高斯点应为区间

上关于权函数

正交的多项式零点-- 次拉盖尔(Laguerre)正交多项式

的零点.高斯系数6.4高斯(Gauss)公式6.4.2常用Gauss求积公式4.高斯--埃尔米特(Gauss-Hermite)求积公式取

，高斯点为

上关于权函数

正交的多项式-- 次埃尔米特正交多项式

的零点.高斯系数

的计算公式为：求积公式为：6.4高斯(Gauss)公式6.4.2常用Gauss求积公式【例6.8】用高斯--勒让德公式计算积分 (精确值 )解：先作变换，则积分区间[0,1化为[-1,1]，且用两点高斯-勒让德公式得(四位有效数字):6.4高斯(Gauss)公式6.4.3高斯公式的应用【例6.8】用高斯--勒让德公式计算积分 (精确值 )解：用三点高斯-勒让德公式(7位有效数字)：复合梯形公式要求

个点上的值才能达到同样的精度.6.4高斯(Gauss)公式6.4.3高斯公式的应用【例6.9】用高斯--勒让德公式计算下列积分分别取

，观察误差随高斯点数的变化规律.解：当

时，

在

处分别为连续、一阶连续可导以及二阶连续可导的.令误差用

表示横坐标和纵坐标，则计算结果如下图.6.4高斯(Gauss)公式6.4.3高斯公式的应用6.4高斯(Gauss)公式6.4.3高斯公式的应用Figure1：被积函数的光滑性与高斯公式收敛速度之间的关系从上图可看到无论取0,1或2，随着高斯点的增加，误差也相应地减小，即高斯公式是收敛的.当函数为无穷次可微的解析函数时，可以证明此时误差为指数收敛（或称为无穷次收敛速度）

为正常数.6.4高斯(Gauss)公式6.4.3高斯公式的应用当

较大时，误差

有近似表达式

(为与

无关的正常数).刻划收敛速度的指标

随着被积函数

的光

滑性的提高而变大.从而可知高斯公式的收敛速度与被积函数的光滑性有关.目录/Contents第六章数值积分与数值微分第一节几个常用积分公式及其复合公式第二节变步长方法与外推加速技术第三节高斯公式第四节多重积分的计算第五节数值微分【例6.10】用复合Simpson公式,取

计算积分(真值=1)：解：设

，则其中，

.同理，再次利用Simpson公式得6.5多重积分的计算6.5.1二重积分的计算使用同样思路，可得到矩形区域

上的复合Simpson公式.将

方向和

方向分成

份和

份，记

则计算二重积分的复合Simpson公式为：6.5多重积分的计算6.5.1二重积分的计算6.5多重积分的计算6.5.1二重积分的计算上述复合Simpson公式原则上可推广至任何有界区域上的多重积分.推导如下：取

维长方体将被积函数

延拓为

上的函数 ,则

化为

维长方体上的积分：由于

为长方体，可以使用上述复合Simpson公式计算

的值.6.5多重积分的计算6.5.1二重积分的计算维数

的增大，节点数

急剧增多.这导致计算量大大增加，而计算的误差

的减少却很缓慢，即维数灾难.

维复合Simpson公式的误差R与节点数

的关系为则每个方向的计算误差约为

，而总的求积误差 (及

均为正常数).多重积分的计算误差与计算节点数之间的关系以复合Simpson公式为例，设函数

充分光滑.复合Simpson公式计算时每个方向取

个点由于区域是

维，因此总的节点数为 .计算高维

数值积分的蒙特卡罗(MonteCarlo)方法.6.5多重积分的计算6.5.1二重积分的计算考虑区间[0,1]上的积分，可直接推广到多维的情形.设随机变量

服从[0,1]上的均匀分布，即 . 为任意可积函数，则随机变量

的函数

的数学期望函数

在区间[0,1]上的积分.6.5多重积分的计算6.5.2蒙特卡罗(MonteCarlo)模拟求积法简介记则

是

的一个无偏估计量：由概率论中的大数定理，当

时，

依概率收敛于

，即6.5多重积分的计算6.5.2蒙特卡罗(MonteCarlo)模拟求积法简介估计误差： .由于

仍为随机变量，我们估计其方差(一般认为均方差即是误差).由于

，所以由于

相互独立，可推出当

时，因此6.5多重积分的计算6.5.2蒙特卡罗(MonteCarlo)模拟求积法简介高维区域积分可取

次独立取样值的算术平均方差仍为

，即收敛速度不受积分区域维数

的影响！正因为蒙特卡罗方法具有这个特性，才能成为计算高维积分的有效数值方法.另外我们可以看出，蒙特卡罗方法的收敛速度为

，与被积函数的光滑性无关.这一性质与插值型积分公式大不相同!作为积分

的无偏估计量.6.5多重积分的计算6.5.2蒙特卡罗(MonteCarlo)模拟求积法简介【例6.11】用蒙特卡罗方法计算积分

.解：首先令

，将积分化为标准区间[0,1]上的积分：其次取

为

次相互独立的均匀分布

的样本，则分别取

，计算出积分的近似值 .6.5多重积分的计算6.5.2蒙特卡罗(MonteCarlo)模拟求积法简介程序如下：x=[102050100200

500

100020005000100002000050000];I=quad(‘sin(x)’./x’,0,2,1.e-16);fori=1:length(x)z=rand(x(i),1);In=sin(2*z)./z;In=mean(In);y(i)=abs(In-I);endx=log(x);y=log(y);plot(x,y,’.-’);xlabel(‘logN’);ylabel(‘logR_N’);6.5多重积分的计算6.5.2蒙特卡罗(MonteCarlo)模拟求积法简介6.5多重积分的计算6.5.2蒙特卡罗(MonteCarlo)模拟求积法简介Figure2：误差

随节点的变化规律目录/Contents第六章数值积分与数值微分第一节几个常用积分公式及其复合公式第二节变步长方法与外推加速技术第三节高斯公式第四节多重积分的计算第五节数值微分插值余项令6.6数值微分6.6.1基于拉格朗日插值多项式的求导方法插值余项两点公式同理6.6数值微分6.6.1基于拉格朗日插值多项式的求导方法三点公式6.6数值微分6.6.1基于拉格朗日插值多项式的求导方法【例6.12】给定函数的下列数据表，试利用二点、三点微分公式计算处的一阶导数值。解：两点公式6.6数值微分6.6.1基于拉格朗日插值多项式的求导方法2.52.62.72.82.912.182513.463714.879716.444618.1741解：三点公式6.6数值微分6.6.1基于拉格朗日插值多项式的求导方法【例6.12】给定函数的下列数据表，试利用二点、三点微分公式计算处的一阶导数值。2.52.62.72.82.912.182513.463714.879716.444618.17416.6数值微分6.6.1基于拉格朗日插值多项式的求导方法解：【例6.12】给定函数的下列数据表，试利用二点、三点微分公式计算处的一阶导数值。2.52.62.72.82.912.182513.463714.879716.444618.17416.6数值微分6.6.1基于拉格朗日插值多项式的求导方法【例6.13】给定下列数据表,试利用三点微分公式计算

处的一阶和二阶导数值。解：1.11.21.30.48600.86161.5975

越小一般精度越高.但在实际计算中，数据

有误差，并不是

越小计算效果越好.理由说明如下.设

，并令

，则可得6.6数值微分6.6.2基于样条函数的求导方法由此可见系数

有正数，也有负数.并且当数据

有一定误差以及

很小时，由于

的表达式中会出现分母为“小数”

，因此会造成较大的数值误差，即算法不是一个稳定的算法.因此实际应用中步长

不要取得太小.基于样条函数的求导方法我们将用样条函数

代替拉格朗日插值多项式

作为函数

的近似.若

，则有估计其中

.即此时不但

与

的函数值很“接近”，它们的导数值也很“接近”.如何求

？三弯矩方法、三转角方法。计算三次样条函数的三弯矩方法已在第3章介绍，下面我们将推导直接以节点处的导数值

作为未知量的三转角方程组.6.6数值微分6.6.2基于样条函数的求导方法设 .假定样条函数

在

处的导数值

，则在区间

上，

为三次多项式.利用Hermite插值公式，可得求出

的值，需要

在

处的连续性条件:以及在区间

端点

及

处的附加条件.6.6数值微分6.6.2基于样条函数的求导方法1.设

为已知.使用与第3章类似的方法可知

满足下列方程组6.6数值微分6.6.2基于样条函数的求导方法2.若

已知，则方程组为6.6数值微分6.6.2基于样条函数的求导方法3.若满足周期性条件：

，则方程组变为6.6数值微分6.6.2基于样条函数的求导方法AdvancedNumericalComputing现代数值计算同济大学数学科学学院学海无涯，祝你成功！第七章AdvancedNumericalComputing非线性方程求解现代数值计算同济大学数学科学学院目录/Contents第七章非线性方程求解第一节非线性方程求解的基本问题第二节非线性方程基本迭代方法第三节牛顿法和割线法第四节非线性方程组简介第五节非线性最小二乘问题第六节大范围求解方法7.1非线性方程求解的基本问题在科学工程计算中常常遇到非线性方程和方程组的问题，例：高次代数方程超越方程【定义7.1】求解非线性方程组

若 ,称

是

的根或

的零点。【定义7.2】代数方程:重零点：且非线性方程求根的基本问题：根的存在性，个数和重数有根区间的确定求出足够精度的近似根【定理7.1】（代数基本定理）

次多项式函数 ,其中 ,在复数域中恰有

个根,重根按其重数计算。【定理7.2】（中值定理）若连续函数

在某两个点

上满足

和

，

则在区间

（或

）上至少存在函数

的一个零点。7.1非线性方程求解的基本问题【例7.1】考察非线性方程

的根的个数。7.1非线性方程求解的基本问题7.1非线性方程求解的基本问题【例7.1】考察非线性方程

的根的个数。7.1非线性方程求解的基本问题【例7.1】考察非线性方程

的根的个数。7.1非线性方程求解的基本问题【例7.1】考察非线性方程

的根的个数。7.1非线性方程求解的基本问题【例7.1】考察非线性方程

的根的个数。【例7.2】求解方程

，

是一给定参数。7.1非线性方程求解的基本问题【定义7.3】收敛速度设

为第

个迭代点的误差,若称

为

阶收敛。线性收敛( )、超线性收敛( )、平方收敛( )【例7.3】考察如下几个序列的收敛速度:其中，

是Fibonacci数列，即

，7.1非线性方程求解的基本问题7.1非线性方程求解的基本问题目录/Contents第七章非线性方程求解第一节非线性方程求解的基本问题第二节非线性方程基本迭代方法第三节牛顿法和割线法第四节非线性方程组简介第五节非线性最小二乘问题第六节大范围求解方法7.2.1二分法【定理7.3】（中值定理）若连续函数

在某两个点

上满足

和

，则在区间 (或 )上至少存在函数

的一个零点。若

在

连续，

与

异号，中点 .若

与

同号,则有有根区间 ;若

与

同号,则有有根区间

；反复上面的过程。几何意义:7.2.1二分法【算法7.1】（二分法）(1)给定初始区间 ,满足 ,以及计算精度 ;(2)令 ;(3)若

或者 ,停止算法;(4)若 ,令

,

否则 ;转步骤(2).Matlabprogram

(二分法)7.2.1二分法functionx=bisect(f,a,b,tol)

ifnargin<4,tol=1e-12;

end

fa=feval(f,a);fb=feval(f,b);

whileabs(a-b)>tol,

x=(a+b)/2;

fx=feval(f,x);7.2.1二分法

ifsign(fx)==sign(fa),

a=x;fa=fx;

elseif

sign(fx)==sign(fb),b=x;

fb=fx;

elsereturn;

end

end【例7.4】用二分法求解非线性方程 .函数

满足 ,方程的近似根为 .7.2.1二分法用erfen.m7.2.1二分法【例7.4】用二分法求解非线性方程 .7.2.1二分法【例7.4】用二分法求解非线性方程 .1.00000-1.000002.000005.000001.00000-1.000001.500000.875001.25000-0.296881.500000.875001.25000-0.296881.375000.224611.31250-0.051511.375000.224611.31250-0.051511.343750.082611.31250-0.051511.328130.014581.32031-0.018711.328130.014581.32422-0.002131.328130.014581.32422-0.002131.326170.006211.32422-0.002131.325200.002047.2.2不动点迭代方法类似于线性方程组的迭代解法，考虑非线性方程等价迭代方程迭代算法【定义7.4】若把函数

看成是一个映射，求解

满足

相当于求解

，即求在映射

下不动的点。因此，方程

称为不动点方程，

称为函数

的不动点,该问题也称为不动点问题。【例7.5】求解非线性方程

.

取 .

，对应迭代方法 ;

，对应迭代方法 ;

，对应迭代方法

；

，对应迭代方法7.2.2不动点迭代方法7.2.2不动点迭代方法1.5001.5000001.5000001.500000000000002.3751.3572090.8000001.3478260869565212.3961.330861-2.7777781.325200398950911904.0031.3258840.1488971.324718173999051.324939-1.0226731.324717957244791.32476021.8054621.324717957244751.3247260.0021081.32471795724475四个不同的迭代方法（neex5.m)【定理7.4】不动点迭代设一元函数

在区间

上一阶连续可导,且(1) 对一切

成立;(2)存在常数

，

，使得

对一切

成立。则成立如下结论:(a)]对任何

产生的迭代序列

必收敛于

在区间

上的唯一不动点,即

；(b)]序列

的收敛速度有估计7.2.2不动点迭代方法证：令

，则 .因此，存在

， .

即

若

，则由数学归纳法可得

，7.2.2不动点迭代方法因此7.2.2不动点迭代方法证：全局收敛给出

，使得

。局部收敛性非线性方法的收敛与初始点有关,收敛速度与不动点有关若

(线性收敛)，7.2.2不动点迭代方法(1)若

，其中

是不动点,则存在

的某个邻域 ,使得对于任意

,

有 .(2)原因:无法说明

.几何含义7.2.2不动点迭代方法【例7.6】求解方程

的根。方程等价于有根区间 .7.2.2不动点迭代方法7.2.2不动点迭代方法1.7501.750000000000001.750000000000001.750000000000001.7871.770794952435151.756890790371001.763254784462251.7231.758951903005441.760194735991501.763222834536611.8391.765652618481771.761775690112281.763222834351901.6421.761847231831661.762531452898811.76322283435190-27.4621.763222834351901.76322283435190

7.2.3迭代加速假设有不动点迭代，设不动点为，若变化不大，得用两步迭代近似值得到的某种平均是更好的近似；另外，可以看成是一个新的迭代方法，【例7.7】加速例3.2中的函数的迭代.导数变化不大令7.2.3迭代加速7.2.3迭代加速加速方法1.750000000000001.750000000000001.756890790371001.763379158584341.760194735991501.763220601443681.761775690112281.763222866181401.762531452898811.763222833898161.763222834358361.763222834351901.763222834351901.763222834351901.76322834351907.2.3Aitken加速若仅知

而不知其具体数值假若

变化不大，注意：

都是靠近

的数,，公式出现了分子分母都很靠近

的情况。采用双精度计算方式。7.2.3迭代加速另外的方式:7.2.3Aitken加速【定理7.5】设不动点迭代

的迭代函数

在

的某个邻域内具有二阶连续导数，且

，

，则相应的艾特肯(Aitken)迭代加速是二阶收敛的，迭代序列的极限仍为 .7.2.3迭代加速{【例7.8】对例3.3中的迭代函数做Aitken加速.Aitken加速为该方法每一步的计算量约为原迭代每一步计算量的两倍.7.2.3迭代加速7.2.3Aitken加速艾特肯加速1.750000000000001.750000000000001.770794952435151.763249280656661.758951903005441.763222834456391.765652618481771.763222834351901.76184723183166NaN1.7632283435190目录/Contents第七章非线性方程求根第一节非线性方程求根的基本问题第二节非线性方程基本迭代方法第三节牛顿法和割线法第四节非线性方程组简介第五节非线性最小二乘问题第六节大范围求解方法7.3牛顿法和割线法考虑非线性方程

，近似根为

，在近似根处的Taylor展开为取其线性部分有Newton迭代格式为注：可看作不动点迭代考虑有若迭代格式

，取

，则其加速为特殊的加速法7.3牛顿法和割线法几何含义7.3牛顿法和割线法【定理7.6】牛顿法设

是方程

的根，

在某个包含

为内点的区间内足够光滑，且

那么存在

的一个邻域

，使得对于任意

牛顿法产生的迭代序列以不低于二阶的收敛速度收敛于解

.证：牛顿法是对应于函数

的不动点迭代。若

，则有

.局部收敛7.3牛顿法和割线法证：因此【定理7.6】牛顿法设

是方程

的根，

在某个包含

为内点的区间内足够光滑，且

那么存在

的一个邻域

，使得对于任意

牛顿法产生的迭代序列以不低于二阶的收敛速度收敛于解

.7.3牛顿法和割线法证：

定号，根是唯一的。

定号，可得四种情形：

且

；

且

；

且

；

且

；【定理7.7】给定非线性函数

，若它在区间

上二阶连续可微,满足

，并且对于所有

有

， .若选定初始点

满足

，则牛顿迭代法收敛于方程的唯一解 .7.3牛顿法和割线法

且

；

且

，函数单调递增。由于

，

，因此 .归纳证明：【定理7.7】给定非线性函数

，若它在区间

上二阶连续可微,满足

，并且对于所有

有

， .若选定初始点

满足

，则牛顿迭代法收敛于方程的唯一解 .证：

定号，根是唯一的。

定号，可得四种情形：7.3牛顿法和割线法

且

；

单调下降且有下界

.两边取极限。证：【定理7.7】给定非线性函数

，若它在区间

上二阶连续可微,满足

，并且对于所有

有

， .若选定初始点

满足

，则牛顿迭代法收敛于方程的唯一解 .7.3牛顿法和割线法【定理7.8】设在区间

上有二阶连续可微函数

，

，并且对于所有

，有

，

.若

两点满足则对于任何

，牛顿法迭代收敛于方程的唯一根 .【定理7.7】给定非线性函数

，若它在区间

上二阶连续可微,满足

，并且对于所有

有

， .若选定初始点

满足

，则牛顿迭代法收敛于方程的唯一解 .7.3牛顿法和割线法计算

相当于求解方程

的正根。可证【例7.9】设计一个算法计算

，其中 .7.3牛顿法和割线法【例7.9】设计一个算法计算

，其中 .7.3牛顿法和割线法有效位02.00000000000000011.50000000000000121.41666666666667331.41421568627451641.414213562374691251.414213562373091561.4142135623730915

有重根的情形：设

，其中，

且 .单根7.3牛顿法和割线法7.3牛顿法和割线法牛顿下山法的迭代公式参数

依次取

使

成立.7.3牛顿下山法【算法7.2】牛顿下山法(1)给定初始值

，精度

，

；(2)若

，近似解为

，停止迭代；(3)令

，

；(4)若

，则, 转(5)；否则，

，重复步骤(4)；(5) ，转(2).7.3牛顿法和割线法几何含义7.3牛顿下山法7.3牛顿法和割线法function[x,it,convg]=newton(x0,f,g,maxit,tol)%findthezerooffunctionf,withgradientgprovided%Usage:[x,it,convg]=newton(x0,f,g,maxit,tol)

ifnargin<5,

tol=1e-10;

ifnargin<4,maxit=100;

end;end

x=x0;

fx=feval(f,x);

convg=0;

it=1;

while~convg,

it=it+1;

ifnorm(fx)<=tol,

fprintf('NewtonIterationsuccesses!!\n');

convg=1;

return;

end7.3牛顿下山法7.3牛顿法和割线法

d=-feval(g,x)\fx;

lambda=1;

lsdone=0;

while~lsdone,

xn=x+lambda*d;

fn=feval(f,xn);

ifabs(fn)<abs(fx),

lsdone=1;

else

lambda=1/2*lambda;

iflambda<=eps,

convg=-1;

error('linesearchfails!!');

end

end

end7.3牛顿下山法

x=xn;

fx=fn;

ifit>maxit,

convg=0;

error('Newtonmethodneedsmoreiterations.!!');

end

end7.3牛顿法和割线法例：计算方程的根。编写两个Matlab文件:functionv=f(x)

v=x.^2+\sin(10*x)-1;和functionv=g(x)

v=2*x+10*\cos(10*x);7.3牛顿下山法调用>>[x,it,convg]=newton(30,'f','g’)NewtonIterationsuccesses!!

x=

-0.412101013664971it=

12convg=

1令代入牛顿法几何含义:用割线代替牛顿法中的切线7.3牛顿法和割线法【例6.10】用割线法求方程的解。有根区间:.7.3牛顿法和割线法【例6.10】用割线法求方程的解。有根区间:.7.3牛顿法和割线法牛顿法割线法01.0000000000000000.00000000000000010.8000000000000001.00000000000000020.7568181818181820.50000000000000030.7548814744397500.69230769230769240.7548776662613990.77560339204174850.7548776662466930.75352325251062460.75484958576524170.75487770485289880.75487766624559390.754877666246693【定理7.9】割线法给定非线性方程 .若函数

在其解

的某个邻域内二阶连续可导，且

，则存在

的一个邻域

使得对于任意

，双点割线法产生的序列收敛于解

，且收敛阶为

.单点割线法7.3牛顿法和割线法目录/Contents第七章非线性方程求根第一节非线性方程求根的基本问题第二节非线性方程基本迭代方法第三节牛顿法和割线法第四节非线性方程组简介第五节非线性最小二乘问题第六节大范围求解方法7.4非线性方程组简介记：

定义:若对于某个

，存在

的一个邻域

，则称

是

的内点。称

在

点处可微,若【算法7.3】非线性方程组的高斯--赛德尔迭代方法(1)给定

，

，控制精度

；(2)对

，以

为初始值求解如下问题：

并把其解记为

；(3)若

，则迭代收敛，方程组的近似解为

；

否则，

，转步骤(2).7.4非线性方程组简介【定理7.10】压缩映像定理设

在某区域上满足

，并且存在压缩因子

，使得对于任意

成立 .那么下面的结论成立:在

上存在唯一不动点，即

；对于任意的

，不动点迭代

收敛到该唯一不动点

，且收敛速度有估计7.4非线性方程组简介牛顿法设函数

连续可微，有近似值

，向量形式7.4非线性方程组简介【例7.11】求解非线性方程组不动点迭代高斯---赛德尔迭代7.4非线性方程组简介牛顿法7.4非线性方程组简介不动点高斯-赛德尔牛顿法01(0.45969769,0.15852902)(0.45969769,1.44367720)(0.49023685,1.03629258)2(0.01253943,1.44367720)(0.87322296,1.76640321)(0.24236719,0.74784106)3(0.87322296,1.01253910)(1.19436188,1.92998127)(0.26208909,0.74085322)4(0.47029118,1.76640321)(1.35151131,1.97605323)(0.26211952,0.74087174)5(1.19436188,1.45314588)(1.39425486,198445699)(0.26211953,0.74087174)6(0.88262077,1.92998127)(1.40196392,1.98578163)(0.26211953,0.740