第3章-常微分方程的差分方法 1VIP

1,第3章常微分方程的数值解法,1引言2欧拉方法3龙格-库塔方法4阿达姆斯方法5算法的稳定性及收敛性6方程组及高阶方程的数值解法7边值问题的数值解法,2,1引言在工程和科学技术的实际问题中，常需要解常微分方程。但常微分方程组中往往只有少数较简单和典型的常微分方程（例如线性常系数常微分方程等）可求出其解析解。对于变系数常微分方程的解析求解就比较困难，而一般的非线性常微分方程就更不用说了。在大多数情况下，常微分方程只能用近似法求解。这种近似解法可分为两大类：一类是近似解析法，如级数解法、逐次逼近法等；另一类则是数值解法，它给出方程在一些离散点上的近似解。在具体求解微分方程时，需要具备某种定解条件，微分方程和定解条件合在一起组成定解问题。定解条,3,件有两种：一种是给出积分曲线在初始点的状态，称为初始条件，相应的定解问题称为初值问题；另一种是给出积分曲线首尾两端的状态，称为边界条件，相应的定解问题则称为边值问题。例如，弹簧-质量系统的振动问题（图7-1），作一定的简化后，可用一个二阶常微分方程,来描述。式中，x是质量m离平衡位置（0点）的距离；t是时间；c是弹簧常数。,4,当弹簧在振动开始时刻t=t0时的初始位置x(t0)=x0和初速度,确定时，弹簧的振动规律x(t)也就唯一确定。这就是一个常微分方程的初值问题，可写成：,5,图7-1,本章先从一阶常微分方程的初值问题：,（1.1）,出发进行讨论。,6,由常微分方程的理论知，只要上式中的函数f(x,y)在区域内连续，且关于y满足李普希兹(Lipschitz)条件，即存在常数L(它与x,y无关)使,对内任意两个y1和y2都成立，则方程的解必定存在且唯一。下面的分析均假定满足上述条件。初值问题（1.1）的数值解法，常采用差分方法，即把一个连续的初值问题离散化为一个差分方程来求解。即将（1.1）离散化后，求找其解y=y(x)在一系列离散节点,7,上的近似值y0,y1,yn。两相邻节点间的距离,称为步长。当（常值）时称为等步长，有,或,a=x0x1xi0都是绝对稳定的。,稳定区域：能保证算法稳定的的取值范围，称为稳定区域。对于（5.1）有,102,显然，稳定区域越大，意味着该算法的稳定性越好。由（5.2）可知，显式欧拉方法的稳定区域为：,由（5.3）可知，隐式欧拉方法的稳定区域为：,103,可见隐式欧拉方法的稳定性比显式欧拉方法的好。可用与上述类似的方法来分析k步阿达姆斯显式和隐式方法的稳定性。对于典型方程，它们的稳定区域如分别记为和，则其数值如表7-9所示，由表79可见：（1）阿达姆斯隐式方法的稳定区域都比同阶的显式方法的大，这是隐式方法最大的优点；（2）k越大（步数、阶数越大）时稳定区域就越小，只有隐式欧拉方法是绝对稳定的。,104,105,5.2收敛性微分方程在离散为差分方程来求解时，存在着差分方程的解yi当步长h0时能否收敛到微分方程的准确解y(xi)的问题，对于任意节点xi=x0+ih，当h0（同时i）时，yi都能趋向于准确解y(xi)的算法，称为是收敛的，否则称为不收敛。这里简单地以欧拉方法为例，来分析其收敛性。在欧拉格式,（5.4）,106,中，如取yi=y(xi)，所得yi+1为,（5.5）,由（2.3）可知，其局部截断误差为：,必定存在常数C，使,（5.6）,107,局部截断误差的积累误差即为整体截断误差。下面考虑整体截断误差。记,其中y(xi)为（1.1）的精确解，yi为由欧拉公式得到的近似解。,（5.4）（5.5）有,（5.7）,108,将（5.7）和（5.6）代入上式，可得,式中L是f关于y的李普希兹(Lipschitz)常数。因,将上式反复递推后可得,109,（5.8）,或,设xi-x0=ihT(T为常数)1+hLehL，(1+hL)ieihLeTL把上式代入（5.8），得,110,如不计初值误差，即0=0，则有,（5.9）,上式说明，当h0时，i，即yiy(xi)，因此，欧拉方法是收敛的，且其收敛速度为O(h)，即具有一阶收敛速度。（5.9）还说明欧拉方法的整体截断误差为O(h)，即是一阶的，因此算法的精度为一阶。由（2.3）有,111,即欧拉方法的局部截断误差是二阶的，可见其整体截断误差（一阶）要比局部截断误差低一阶，这个结论对其他算法也适用，是带有一般性的。,112,6方程组及高阶方程的数值解法*6.1一阶方程组对于一阶微分方程组的初值问题：,（6.1）,可以将单个方程中的f和y看作向量来处理，这样就可把前面介绍的各种差分算法推广到一阶方程组中应用。,113,设xi=x0+ih(i=1,2,3,),yi，zi为节点xi上的近似解，则有改进的欧拉格式为：,（6.2）,114,再如经典R-K格式为：,（6.3）,其中,(6.4),115,将节点xi上的yi和zi值代入（6.4），依次算出K1，K2，L2，K3，L3，K4和L4，然后代入（6.3）算出节点xi+1上的yi+1和zi+1值。对于具有三个或三个以上方程的方程组的初值问题，也可用类似方法处理.。此外，多步方法也同样可以应用于求解方程组的初值问题。,(6.4),116,6.2高阶方程对于高阶微分方程（或方程组）的初值问题，可以化为一阶方程组来求解。例如，有二阶微分方程的初值问题,（6.5）,（6.6）,令新变量，代入（6.5）得,117,（6.6）为一个一阶微分方程组的初值问题，对此可用6.1节中的方法来求解。例如用经典R-K方法（6.3）：,(6.7),由（6.4）可知，式中,118,(6.8),119,将上式中的K1，K2，K3，和K4代入（6.7），可消去K1，K2，K3，和K4，化简成,（6.9）,式中,（6.10）,120,上述方法同样可以用来处理三阶或更高阶的微分方程（或方程组）的初值问题。例5试求解节下列二阶微分方程的初值问题：,解先作变换：令，代入上式得一阶方程组：,121,用经典R-K方法求解，按（6.7）及（6.8）进行计算：取h.，,i=0时：,122,123,124,125,然后计算i=1时的K1，L1，K2，L2，K3，L3，K4，L4，y2和z2;再计算i=2时K1，L1，y3和z3；依此类推，直到i=9时的y10和z10，即可得到数值解：,y1，y2，y10,126,7边值问题的数值解法*以二阶常微分方程为例来讨论边值问题：,（2）第二边值条件：,（1）第一边值条件：,其边值条件可分为三类：,127,式中,（3）第三边值条件：,常微分方程的边值问题有时可化为初值问题来求解，从而可采用前面已讲的各种求解初值问题的单步法和多步法。一般常采用差分方法，即将边值问题离散为差分方程组来求解。现以下列二阶线性常微分方程为例来讨论边值问题的差分解法：,128,(7.1),7.1差分方程组的建立首先将区间均分为n等份，即步长,各节点为：,129,将在上求解的问题化成在这些节点上求的近似值的问题。首先将（7.1）在节点上离散化，对于用一阶和二阶差商来分别近似代替一阶和二阶导数：,130,并设,即得（7.1）的下列近似方程，即差分方程：,（7.2）,131,的联立线性代数方程组，而其方程的个数只有n-1个，要使方程组有唯一解成在，还需补充两个方程，这可由不同的边值条件来补足。对于第一边值条件，直接就可得到另外两个方程：,（7.2）是具有n+1个未知数,于是可得第一边值问题的差分方程组：,132,（7.3）,对于第二边值条件,可用最简单的一阶差商来近似代替，,即,133,如要求误差达到O(h2)，则可用过三点的一阶微分公式来近似代替，即,134,从而可得第二边值问题差分方程组：,(7.4),同样，可得第三边值条件的差分方程组：,135,综上所述，通过离散化过程，可将求微分方程的近似解化成差分方程组的求解问题。7.2差分方程组的求解以由第一边值问题建立起来的差分方程组（7.3）为例，来讨论它的求解。在（7.3）中，可设法消去已,(7.5),136,(7.6),知的和，经整理后可得到关于的下列方程组：,137,显然，方程组（7.6）是对角占优的对角线型方程组，其系数矩阵是严格对角占优的三对角型矩阵，故用追赶法来求解特别有效。,138,小结本章着重讨论常微分方程初值问题的数值解法，主要有欧拉方法，龙格库塔方法和阿达姆斯方法等。它们都是基于把一个连续的定解问题离散化为一个差分方程来求解，是一种步进式的方法。前两种是单步方法；后一种则是多步方法，可以获得较高的精度。通过对算法的数值稳定性和收敛性概念的阐述和分析，了解到隐式结构的稳定性比显示结构的好，而预测校正方法则是把显示和隐式结合起来，使隐式结构既能通过显式来计算，又能保持隐式结构的稳定性好等优点。,139,在实际应用中，如何选择算法并非一件容易的事，这既要考虑到算法的精确能否满足实际问题的要求，也要考虑算法的简易程度和计算工作量的大小，更要注意步长的选择、算法的收敛性和稳定性等。否则会导致计算失败。一般说来，龙格库塔法比较为常用，它特别适用于多步方法（如阿达姆斯方法）中作初值计算；也可独立应用于精度要求不高，函数f较为简单的场合，因为它需要计算函数值的次数多。当函数较复杂时，可用显示阿达姆斯方法或阿达姆斯预测校正方法，后者的稳定性较好。变步长处理和利用李查逊外推法来提高算法的,140,计算精度，在某些场合中也是值得推荐的。对于一阶常微分方程组的初值问题，单个微分方程初值问题的上述各种数值解法都能推广应用，只要把方程组中的y和f都理解为向量即可。至于高阶常微分方程或方程组的初值问题，原则上总是可以把它化成一阶方程组来求解。本章最后也简要阐述了常微分方程边值问题的数值解法。介绍了离散后差分方程组的建立和求解的方法。,141,习题七试分别用欧拉方法和改进的欧拉方法求解下列边值问题，比较它们所得结果，并与精确解作比较：（1）（2）（取h=0.2）,（取h=0.1）,（已知精确解为）,142,(4)2（1）试用欧拉方法求解初值问题：,(取h=0.5),(取h=0.2),143,并证明其截断误差为（2）试证明改进的欧拉方法能准确地解上述初值问题。,3试用梯形格式求解初值问题（取h=0.1）：并证明：（1）其近似解为（2）当h0时，近似解,144,试分别用欧拉方法（取h=0.025），改进的欧拉方法（取h=0.05）和经典R-K方法（取h=0.1）求下列初值问题，并比较它们计算结果的精度：(1)(2)(3),145,试分别用显式和隐式二阶阿达姆斯方法求解下列初值问题：取h=0.2,y(0.2)=0.181来