[整理]常微分方程数值解

1、第八章常微分方程数值解姓名 学号 班级习题主要考察点：欧拉方法的构造，单步法的收敛性和稳定性的讨论，线性多步法中亚当姆 斯方法的构造和讨论。1用改进的欧拉公式，求以下微分方程,2xy =y -一，yxwo,ij（0） =1（改进的尤拉公式的应用）.dz,有一- 2z = 2x dx的数值解（取步长 h =0.2）,并与精确解作比较。2解：原方程可转化为yy'= y2 _2x,令z =y-2解此一阶线性微分方程，可得y = J2x +1。利用以下公式'2为、yp = y +0.2 (yi ) y2xUc =yi +0.2 (yp ) (i =0,1,2, 3, 4) ypyf =

2、；(yp +丫。)求在节点xi = 0.2，i (i =1,2,3,4,5)处的数值解yi,其中，初值为x0=0,y0=1。MATLAB程序如下：x(1)=0;%初值节点y(1)=1;涮值fprintf( 'x(%d)=%f,y(%d)=%f,yy(%d)=%fn',1,x(1),1,y(1),1,y(1);for i=1:5yp=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)

3、;研确解fprintf( 'x(%d)=%f,y(%d)=%f,yy(%d)=%fn',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.000000x(2)=0.200000, y(2)=1.220000, yy(2)=1.183216x(3)=0.400000, y(3)=1.420452, yy(3)=1.341641x(4)=0.600000, y(4)=1.615113, yy(4)=1.483240x(5)=0.800000, y(5)=1.814

4、224, yy(5)=1.612452x(6)=1.000000, y(6)=2.027550, yy(6)=1.732051'yy 12用四阶龙格库塔法求解初值问题y y ，取h=0.2，求x = 0.2,0.4时的数值解7(0)=0要求写出由h,xn,yn直接计算yn中的迭代公式，计算过程保留3位小数。(龙格库塔方法的应用)解：四阶龙格-库塔经典公式为hyn 1 = yn(ki2k2 2k3 kJ6k2k3一1.1、f (xn h, yn -hki) 22一1 .1、f (x- h, yn - hk2)22k4f (xn h, yn hk3)由于f (x, y) =1 - y ,在

5、各点的斜率预报值分别为: k1 =1 - ynk2 =1 - (yn hk1)=1 - yn -h(1 - yn) =(1 - yn )(1-7)222hhhhhk3=1-(ynk2)=1-yn 。7口)(1)=(1-丫口)1(1)22222hhhhk4 =1 -(ynhk3) =1 -yn -h(1 -yn)1-(1 -) =(1-yn)1-h(1-(1 -)2222四阶经典公式可改写成以下直接的形式：h2 hyn 1 = yn(1 - yn)(6 -3h h )64在x = x1 =0.2处，有3y1 =0 02(1 -0)(6 -3 0.2 (0.2)2 -(-) =0.181364在x

6、 =x2 = 0.4处，有30.22 (0.2)y2 =0.1813 (1 -0.1813)(6-3 0.2 (0.2) -(-) -0.329764注：这两个近似值与精确解 y =1 -e在这两点的精确值十分接近。3用梯形方法解初值问题y+ y =0：y(0) =1证明其近似解为并证明当hT 0时，它收敛于原初值问题的准确解y = e。解：显然，y =e”是原初值问题的准确解。求解一般微分方程初值问题的梯形公式的形式为 hyn 1 =yn - f (Xn, yn) f (Xn 1, 丫口 .)对于该初值问题，其梯形公式的具体形式为hhh2 - hyn4 =yn +二(yn yn4) , (1

7、 +；)yn 书=(1 二)yn， yn 书=(；TZ7)yn2222 h)2 -h2 h)yn但二h仔2 - h俨<2 + h J 九，西,亦即:、. xc2h 111.注息、至 U: xn = 0+nh=nh, n =，令 t = -, = - Wh2 h h t 2xn2 2h V能e,Xn2Xnyn= 1-=(1+t) t 2 =(1+t) t (1+t) 2、2 +h JXnXn从而 lhimQyn= |imQ(1 t) t lim (1 t) 2 = e%即：当hT 0时，yn收敛于原初值问题的准确解y(Xn)=ei。4对于初值问题；y'= -10yJ(0)=1,证明

8、当h <0.2时，欧拉公式绝对稳定。(显式和隐式欧拉公式的稳定性讨论)证明：显式的欧拉公式为yn 1 = yn，hf (Xn, yn) = (1 - 10h)yn隐式的欧拉公式为从而 en4=(110h)en,由于 0<h<0.2, -1 <1-10h <1, en| <|en 因此，显式欧拉公式绝对稳定。yn 1 =ynhf (Xn 1, Yn 1) = Yn -1。卜丫口 1一ynyn + )1 10henen 1 二1 - 10h因此，隐式的欧拉公式也是绝对稳定的。h .5证明：梯形公式yn =yn+f(xn,yn)+ f(xnMynQ无条件稳定。(梯

9、形公式的稳定性讨论) 解：对于微分方程初值问题r Vr.y = 一九 y/7 (人 A 0 )、y(0) =1其隐式的梯形公式的具体形式可表示为、，2-h、)%' yn +=(二:F)yn2 hhh、yn +=yn +二一九yn 九yn+，(1 + )yn+ =(1 222 - h、从而 en 1- =()en2 h由h>0,九>0可知,c /2 +九 hJc ° 6书|<(2)6=异,故隐式的梯形公式无条件稳定。6设有常微分方程的初值问题；y'= f (x,y)=丫幅)二v。，试用泰勒展开法，构造线性两步法数值计算公式yn由=a(yn + y。)+

10、h(P0fn + P f。)，使其具有二阶精度，并推导其局部截断误差主项。(局部截断误差和主项的计算) 解：假设yn =y(Xn), yn4=y(Xn)，利用泰勒展式，有y (Xn), 2 y (Xn) .3 .yn"y(Xn)= y(Xn) -丫(Xn)h h - 一二一 h 26fn = f(Xn,yn) = f (Xn, y(Xn) = y (Xn) fnL W yn,)二 f (XnJ, y(XnJ 二 /(心二 y(Xn) - y (Xn)h 号 h2 -yn 1, =2：y(Xn) ( F J -1)y (Xn)h - (- - 、)y (Xn)h2 - (- 一)y (

11、Xn)h3 2621 .213又 y(Xn 1) = y(Xn) y (Xn)hy (Xn)hy (Xn)h2 6欲使其具有尽可能高的局部截断误差，必须2豆=1,久 +P _口 =1 , - ?1 =- 22从而 口 =1, P。=7, P1 =-1 244一，171于是数值计算公式为L书=(yn +yn)+h( % 二3)。244该数值计算公式的局部截断误差的主项为15y(Xn1)-Yn1 =(-v)y (Xn)h3:丫6 62247已知初值问题y y = 2xy(0) =oy(0.1) =0.01取步长h =0.1 ,利用阿当姆斯公式ynd1 = yn +-(3fn - fn)，求此微分方

12、程在0 ,102上的数值解，求此公式的局部截断误差的首项。(阿当姆斯公式的应用)解：假设yn =y(xn) , yn=y(xn)，利用泰勒展开，有 yn = y(xn), fn = yH(xn), fn=yxn)=yxn) -y"(xn)h + y 卜)h2 -2yn 1 = y(xn) y (Xn)h 1y (Xn)h2 - y (Xn)h3241 213而 y(xn 1) = y(xn) y (xn)h y (xn)hy (xn)h2 611. Q5. Qy(xn 1) - yn 1 =(- -)y (Xn)hy (xn)h5c该阿当姆斯两步公式具有 2阶精度，其局部截断误差的主

13、项为-5y7xn)h3o取步长h = 0.1 ,节点Xn =0.1n ( n = 0,1 , 2，100),注意到f (X, y) = 2x ,其计算公 式可改写为0.1 - 一 、 yn 1 = yn k(6xn -2xn4) = yn 0.02n 0.012仅需取一个初值y0 = 0 ,可实现这一公式的实际计算。其 MATLAB 下的程序如下：x0=0; %初值节点y0=0; %初值for n=0:99y1=y0+0.02*n+0.01;x1=x0+0.1;,n+1,x1,n+1,y1);fprintf( 'x(%3d)=%10.8f,y(%3d)=%10.8fn'x0=x

14、1;y0=y1;end运行结果如下：x( 1)=0.10000000,y( 1)=0.01000000x( 2)=0.20000000,y( 2)=0.04000000x( 3)=0.30000000,y( 3)=0.09000000x( 4)=0.40000000,y( 4)=0.16000000x( 5)=0.50000000,y( 5)=0.25000000x( 6)=0.60000000,y( 6)=0.36000000x( 7)=0.70000000,y( 7)=0.49000000x( 8)=0.80000000,y( 8)=0.64000000x( 9)=0.90000000,y

15、( 9)=0.81000000x( 10)=1.00000000,y( 10)=1.00000000x( 11)=1.10000000,y( 11)=1.21000000x( 12)=1.20000000,y( 12)=1.44000000x( 13)=1.30000000,y( 13)=1.69000000x( 14)=1.40000000,y( 14)=1.96000000x( 15)=1.50000000,y( 15)=2.25000000x( 16)=1.60000000,y( 16)=2.56000000x( 17)=1.70000000,y( 17)=2.89000000x( 18

16、)=1.80000000,y( 18)=3.24000000x( 19)=1.90000000,y( 19)=3.61000000x( 20)=2.00000000,y( 20)=4.00000000x( 21)=2.10000000,y( 21)=4.41000000x( 22)=2.20000000,y( 22)=4.84000000x( 23)=2.30000000,y( 23)=5.29000000x( 24)=2.40000000,y( 24)=5.76000000x( 25)=2.50000000,y( 25)=6.25000000x( 26)=2.60000000,y( 26)=

17、6.76000000x( 27)=2.70000000,y( 27)=7.29000000x( 28)=2.80000000,y( 28)=7.84000000x( 29)=2.90000000,y( 29)=8.41000000x( 30)=3.00000000,y( 30)=9.00000000x( 31)=3.10000000,y( 31)=9.61000000x( 32)=3.20000000,y( 32)=10.24000000x( 33)=3.30000000,y( 33)=10.89000000x( 34)=3.40000000,y( 34)=11.56000000x( 35)=

18、3.50000000,y( 35)=12.25000000x( 36)=3.60000000,y( 36)=12.96000000x( 37)=3.70000000,y( 37)=13.69000000x( 38)=3.80000000,y( 38)=14.44000000x( 39)=3.90000000,y( 39)=15.21000000x( 40)=4.00000000,y( 40)=16.00000000x( 41)=4.10000000,y( 41)=16.81000000x( 42)=4.20000000,y( 42)=17.64000000x( 43)=4.30000000,y

19、( 43)=18.49000000x( 44)=4.40000000,y( 44)=19.36000000x( 45)=4.50000000,y( 45)=20.25000000x( 46)=4.60000000,y( 46)=21.16000000x( 47)=4.70000000,y( 47)=22.09000000x( 48)=4.80000000,y( 48)=23.04000000x( 49)=4.90000000,y( 49)=24.01000000x( 50)=5.00000000,y( 50)=25.00000000x( 51)=5.10000000,y( 51)=26.010

20、00000x( 52)=5.20000000,y( 52)=27.04000000x( 53)=5.30000000,y( 53)=28.09000000x( 54)=5.40000000,y( 54)=29.16000000x( 55)=5.50000000,y( 55)=30.25000000x( 56)=5.60000000,y( 56)=31.36000000x( 57)=5.70000000,y( 57)=32.49000000x( 58)=5.80000000,y( 58)=33.64000000x( 59)=5.90000000,y( 59)=34.81000000x( 60)=

21、6.00000000,y( 60)=36.00000000x( 61)=6.10000000,y( 61)=37.21000000x( 62)=6.20000000,y( 62)=38.44000000x( 63)=6.30000000,y( 63)=39.69000000x( 64)=6.40000000,y( 64)=40.96000000x( 65)=6.50000000,y( 65)=42.25000000x( 66)=6.60000000,y( 66)=43.56000000x( 67)=6.70000000,y( 67)=44.89000000x( 68)=6.80000000,y

22、( 68)=46.24000000x( 69)=6.90000000,y( 69)=47.61000000x( 70)=7.00000000,y( 70)=49.00000000x( 71)=7.10000000,y( 71)=50.41000000x( 72)=7.20000000,y( 72)=51.84000000x( 73)=7.30000000,y( 73)=53.29000000x( 74)=7.40000000,y( 74)=54.76000000x( 75)=7.50000000,y( 75)=56.25000000x( 76)=7.60000000,y( 76)=57.76000000x( 77)=7.70000000,y( 77)=59.29000000x( 78)=7.80000000,y( 78)=60.84000000x( 79)=7.90000000,y( 79)=62.41000000x( 80)=8.000