常微分方程数值解

第八章 常微分方程数值解 姓名 学号 班级 习题主要考察点：欧拉方法的构造，单步法的收敛性和稳定性的讨论，线性多步法中亚当姆斯方法的构造和讨论。 1 用改进的欧拉公式，求以下微分方程

]1,0[1)0(2xyy

xyy

的数值解（取步长2.0h），并与精确解作比较。（改进的尤拉公式的应用） 解：原方程可转化为 xyyy22，令22yz，有xzdxdz22

解此一阶线性微分方程，可得 12xy。 利用以下公式

)4,3,2,1,0()(21)2(2.0)2(2.01iyyyyxyyyyxyyy

cpipipic

iiiip

求在节点)5,4,3,2,1(2.0iixi处的数值解iy，其中，初值为1,000yx。 MATLAB程序如下： x(1)=0;%初值节点 y(1)=1;%初值 fprintf('x(%d)=%f,y(%d)=%f,yy(%d)=%f\n',1,x(1),1,y(1),1,y(1)); for i=1:5 yp=y(i)+0.2*(y(i)-2*x(i)/y(i));%预报值 yc=y(i)+0.2*(yp-2*x(i)/yp);%校正值 y(i+1)=(yp+yc)/2;%改进值 x(i+1)=x(i)+0.2;%节点值 yy(i+1)=sqrt(2*x(i+1)+1);%精确解 fprintf('x(%d)=%f,y(%d)=%f,yy(%d)=%f\n',i+1,x(i+1),i+1,y(i+1),i+1,yy(i+1)); end 程序运行的结果如下： x(1)=0.000000, y(1)=1.000000, yy(1)=1.000000 x(2)=0.200000, y(2)=1.220000, yy(2)=1.183216 x(3)=0.400000, y(3)=1.420452, yy(3)=1.341641 x(4)=0.600000, y(4)=1.615113, yy(4)=1.483240 x(5)=0.800000, y(5)=1.814224, yy(5)=1.612452 x(6)=1.000000, y(6)=2.027550, yy(6)=1.732051 2用四阶龙格－库塔法求解初值问题0)0(1yyy，取2.0h, 求4.0,2.0x时的数值解. 要求写出由nnyxh,,直接计算1ny的迭代公式，计算过程保留3位小数。（龙格－库塔方法的应用） 解：四阶龙格-库塔经典公式为

11234(22)6nnhyykkkk

1(,)nnkfxy

2111(,)22nnkfxhyhk

3211(,)22nnkfxhyhk

43(,)nnkfxhyhk 由于yyxf1),(，在各点的斜率预报值分别为： nyk11 )21)(1()1(21)2(112hyyhykhyknnnn

)]21(21)[1()21)(1(21)2(123hhyhyhykhyknnnn ))]21(21(1)[1()]21(21)[1(1)(134hhhyhhyhyhkyknnnn 四阶经典公式可改写成以下直接的形式：

)436)(1(6321hhhyhyynnn

在2.01xx处，有

1813.0)4)2.0()2.0(2.036)(01(62.00321y 在4.02xx处，有 3297.0)4)2.0()2.0(2.036)(1813.01(62.01813.0322y 注：这两个近似值与精确解xey1在这两点的精确值十分接近。 3 用梯形方法解初值问题 1)0(0yyy

证明其近似解为 nnhhy2

2

并证明当0h时，它收敛于原初值问题的准确解xey。 解：显然，xey是原初值问题的准确解。 求解一般微分方程初值问题的梯形公式的形式为 )],(),([2111nnnnnnyxfyxfhyy 对于该初值问题，其梯形公式的具体形式为 )(211nnnnyyhyy，nnyhyh)21()21(1，nnyhhy))22(1

于是： 101121222222))22(nnnnnhhyhhyhhyh

hy

亦即：nnhhy22 注意到：nhnhxn0，hxnn，令hht22，2111th有 22)1()1()1(221nnnnnxtxxtxh

x

nttthhy

从而 nnnxxttxtnhetty2000)1(lim)1(limlim 即：当0h时，ny收敛于原初值问题的准确解nxnexy)(。 4对于初值问题1)0(10yyy，证明当2.0h时，欧拉公式绝对稳定。（显式和隐式欧拉公式的稳定性讨论） 证明：显式的欧拉公式为nnnnnyhyxhfyy)101(),(1 从而nnehe)101(1，由于2.00h，11011h，nnee1 因此，显式欧拉公式绝对稳定。 隐式的欧拉公式为111110),(nnnnnnhyyyxhfyy

hyynn1011，heenn1011

由于h0，110110h，nnee1 因此，隐式的欧拉公式也是绝对稳定的。 5证明：梯形公式)],(),([2111nnnnnnyxfyxfhyy无条件稳定。（梯形公式的稳定性讨论） 解：对于微分方程初值问题

)0(1)0(



yyy

其隐式的梯形公式的具体形式可表示为 ][211nnnnyyhyy，nnyhyh)21()21(1，nnyhhy)22(1

从而nnehhe)22(1

由0h，0可知，nnneehhe)22(1，故隐式的梯形公式无条件稳定。

6设有常微分方程的初值问题00)(),(yxyyxfy，试用泰勒展开法，构造线性两步法数值计算公式)()(11011nnnnnffhyyy，使其具有二阶精度，并推导其局部截断误差主项。（局部截断误差和主项的计算） 解：假设)(nnxyy，)(11nnxyy，利用泰勒展式，有

32116)(2)()()()(hxyhxyhxyxyxyynnnnnn

)())(,(),(nnnnnnxyxyxfyxff

21111112)()()()())(,(),(hxyhxyxyxyxyxfyxffnnnnnnnnn

3121101)()26()()2()()()(2hxyhxyhxyxyynnnnn



 又 321)(61)(21)()()(hxyhxyhxyxyxynnnnn 欲使其具有尽可能高的局部截断误差，必须 12，110，2121

从而 21，470，411 于是数值计算公式为 )4147()(21111nnnnnffhyyy。 该数值计算公式的局部截断误差的主项为

33111)(245)()2661()(hxyhxyyxynnnn



7已知初值问题







01.0)1.0(0)0(2yyxy

取步长1.0h，利用阿当姆斯公式)3(211nnnnffhyy，求此微分方程在[0，10]上的数值解，求此公式的局部截断误差的首项。（阿当姆斯公式的应用） 解：假设)(nnxyy，)(11nnxyy，利用泰勒展开，有

)(nnxyy，)(nnxyf，2112)()()()(hxyhxyxyxyfnnnnn

3214)()(21)()(hxyhxyhxyxyynnnnn

而321)(61)(21)()()(hxyhxyhxyxyxynnnnn 3311)(125)()4161()(hxyhxyyxynnnn

该阿当姆斯两步公式具有2阶精度，其局部截断误差的主项为3)(125hxyn。 取步长1.0h，节点nxn1.0（100,,2,1,0n），注意到xyxf2),(，其计算公式可改写为

01.002.0)26(21.011nyxxyynnnnn

仅需取一个初值00y，可实现这一公式的实际计算。 其MATLAB下的程序如下： x0=0;%初值节点 y0=0;%初值