常微分方程解法小结

第六章 常微分方程数值解法 ．143． 第六章 常微分方程数值解 一、基本提要 1 . Euler Euler 方法的基本思 路是用差 商近似导 数，即在 等距节点 nhxxn += 0 ),,2,1,0( Nn L= 上，用向前差商 h xyxy nn )()( 1 −+ 代替 )( nxy′ ，然后代入微分方程 ⎪⎩ ⎪⎨ ⎧ = = )( ),( 00 yxy yxf dx dy 则得 ))(,()()( 1 nnnn xyxfh xyxy ≈−+ ),1,0( L=n 化简得 ))(,()()( 1 nnnn xyxhfxyxy +≈+ 用 ny 近似代替 )( nxy ，所得结果作为 )( 1+nxy 的近似值，记为 1+ny ，则有 ),(1 nnnn yxhfyy +≈+ 按此式由初值 0y 逐次计算 1 2, ,Ly y 的方法即为 Euler 方法，上式被称为 Euler 公式。 2. 截断误差 假设 )( nn xyy = 是精确的，误差 111 )( +++ −= nnn yxyR 被称为局部截断误差，又简称 截断误差。 3. p 阶方法 如果某种数值方法的局部截断误差为 1( )+pO h ，则称该方法是 p 阶方法或具有 p 阶 精度。 4. 梯形公式 ．144 ． 实用数值分析解指导 对初值问题 ⎪⎩ ⎪⎨ ⎧ = = )( ),( 00 yxy yxf dx dy 的方程两端积分可得 dxxyxfxyxy n n x xnn ))(,()()( 11 ∫ +=−+ ),1,0( L=n 用梯形公式计算右端积分，即有 ))](,())(,([ 2 ))(,( 11 1 +++≈∫ + nnnnxx xyxfxyxfhdxxyxfnn 用 1, +nn yy 代替 )(),( 1+nn xyxy ，所得计算公式 )],(),([ 2 111 +++ ++= nnnnnn yxfyxfhyy 即为求解常微分方程初值问题的梯形公式。 5. 改进的 Euler 法 用 Euler 公式求 1+ny 的一个初步近似值 1+ny 作为预测值，再用梯形公式校正求近似 值 1+ny ，即 ⎪⎩ ⎪⎨ ⎧ ++= += +++ + )],(),([ 2 ),( 111 1 nnnnnn nnnn yxfyxfhyy yxhfyy 校正 预测 上式称为由 Euler 公式和梯形公式得到的预测－校正(Predictor-Corrector)系统，又叫改进 的 Euler 法。 6. R-K 方法的构造 一般地，设近似公式 ⎪⎪⎩ ⎪⎪⎨ ⎧ =++= = += ∑ ∑ − = = + ),,3,2( ),( ),( 1 1 1 1 1 piKbhyhaxfK yxfK Kchyy i j jijnini nn p i iinn 称为 p 阶显式 Runge-Kutta 方法，简称 R-K 方法。其中 iiji cba ,, 都是参数，确定它们的 原则是使近似公式在 ),( nn yx 处的 Taylor 展开式与 )(xy 在 nx 处的 Taylor 展开式的前面的 第六章 常微分方程数值解法 ．145． 项尽可能多地重合，使近似公式有尽可能高的精度。 7. 阿达姆斯（Adams）公式 对 R-K 方法取 r ＝3，并令 01321 ==== −βααα ，由方程组 ⎪⎪⎩ ⎪⎪⎨ ⎧ ==−+− = ∑ ∑ ∑ = −= − = 3 1 3 1 1 3 0 )4,3,2,1( 1)1()( 1 i i i k i k i i kki βα α 可解得 24 9, 24 37, 24 59, 24 55,1 32100 −==−=== ββββα 相应的线性多步公式为 )9375955( 24 3211 −−−+ −+−+= nnnnnn ffffhyy 因为 1 0β− = ，上式称为 Adams 显式公式，它是四阶公式，局部截断误差为 )(])(5)(1[ !5 6)5( 3 1 3 1 45 5 1 hOyii hR n i i iin +−−−−= ∑ ∑ = −= + βα )( 720 251 6)5(5 hOyh n += 如果令 03321 ==== βααα ，由方程组 )9375955( 24 3211 −−−+ −+−+= nnnnnn ffffhyy 可得解 24 1, 24 5, 24 19, 24 9,1 21010 =−==== − ββββα 相应的线性多步公式 )5199( 24 21111 −−+−+ +−++= nnnnnn ffffhyy 即称为四阶 Adams 隐式公式，其局部截断误差为 )( 720 19 6)5(5 1 hOyhR nn +−=+ ．146 ． 实用数值分析解题指导 二、典型例题选解 例 6.1 用 Euler 方法求初值问题 ⎪⎩ ⎪⎨ ⎧ = −= 0)0(y yx dx dy )10( ≤≤ x 的数值解，取步长 1.0=h 。 解 Euler 公式为 ),(1 nnnn yxhfyy +=+ 于是有 )(1.01 nnnn yxyy −+=+ 由 0)0( 0 == yy ，计算得 xi 0.1 0.2 0.3 0.4 0.5 y(xi) 0.000000 0.010000 0.029000 0.056100 0.090490 xi 0.6 0.7 0.8 0.9 1.0 y(xi) 0.131441 0.178297 0.230467 0.287420 0.348678 例 6.2 用 Euler 方法求初值问题 ⎪⎩ ⎪⎨ ⎧ = ++−= 1)0( 1 y xy dx dy )10( ≤≤ x 的数值解（取步长 1.0=h ）。 解 Euler 公式为 ),(1 nnnn yxhfyy +=+ 于是有 )1(1.01 ++−+=+ nnnn xyyy 由 1)0( 0 == yy ，计算得 第六章 常微分方程数值解法 ．147． xi 0.1 0.2 0.3 0.4 0.5 y(xi) 1.000000 1.010000 1.029000 1.056100 1.090490 xi 0.6 0.7 0.8 0.9 1.0 y(xi) 1.131441 1.178297 1.230467 1.287420 1.348678 例 6.3 用改进 Euler 法求初值问题 ⎪⎩ ⎪⎨ ⎧ = −= 0)0(y yx dx dy )10( ≤≤ x 的数值解，取 10,1.0 == nh 。 解 改进的为 Euler 公式为 ⎪⎩ ⎪⎨ ⎧ ++= += +++ + )]~,(),([2 ),(~ 111 1 nnnnnn nnnn yxfyxfhyy yxhfyy 于是有 ⎩⎨ ⎧ +−−−= −−= +++ + )~(05.0 )(1.0~ 111 1 nnnnnn nnnn xyyxyy yxyy 由 0)1( 0 == yy ，计算得 xi 0.1 0.2 0.3 0.4 0.5 y(xi) 0.005000 0.019025 0.041218 0.070802 0.107076 xi 0.6 0.7 0.8 0.9 1.0 y(xi) 0.149404 0.197211 0.249976 0.307228 0.368541 例 6.4 用欧拉预-校方法求解初值问题 ⎩⎨ ⎧ = =++′ 1)1( 0sin2 y xyyy 要求取步长 2.0=h ，计算 )2.1(y 及 )4.1(y 的近似值，小数点后至少保留 5 位。 解 欧拉预-校格式为 ．148 ． 实用数值分析解题指导 ⎪⎩ ⎪⎨ ⎧ ++= += +++ + )]~,(),([2 ),(~ 111 1 nnnnnn nnnn yxfyxfhyy yxhfyy 于是有 ⎩⎨ ⎧ +++−= +−= ++++ + )sin~~sin(1.0 )sin(2.0~ 1 2 11 2 1 2 1 nnnnnnnn nnnnn xyyxyyyy xyyyy 由 1)1( 0 == yy ，计算得 ⎩⎨ ⎧ =≈ = 715488.0)2.1( 63171.0~ 1 1 yy y ⎩⎨ ⎧ =≈ = 52611.0)4.1( 47696.0~ 2 2 yy y 例 6.5 常微分方程的初值问题 )10( 2.1)0( 2 ≤≤ ⎪⎩ ⎪⎨ ⎧ −= = x y xydx dy 用欧拉法分别取步长 1.0=h 和 05.0=h ，并和方程的准确解 235 6)( x xy +−= 做比较。 解 取 1.0=h ， 2.10 −=y ，欧拉计算公式： 2 1 1.0 nnnn yxyy +=+ 取 05.0=h ， 2.10 −=y ，欧拉计算公式： 2 1 ~05.0~~ nnnn yxyy +=+ 计算结果见下 nx ny )( nn xyy − nx ny~ )(~ nn xyy − 0.0 -1.2 0.0 0.1 -1.2 0.3 0.2 -1.1856 0.00724294 0.3 -1.15749 0.0143879 0.0 -1.2 0.05 -1.2 0.0 0.1 -1.1964 0.0018027 0.15 -1.18924 0.00359981 0.2 -1.17864 0.0053799 0.25 -1.16474 0.00713093 0.3 -1.14779 0.0088403 第六章 常微分方程数值解法 ．149． 0.4 -1.11729 0.0212262 0.5 -106736 0.0275306 0.6 -1.0104 0.0330812 0.7 -0.949143 0.0376992 0.8 -0.886082 0.04412751 0.9 -0.823271 0.0437814 1.0 -0.762271 0.0452661 0.35 -1.12803 0.0104949 0.4 -1.10576 0.0120815 0.45 -1.0813 0.0135872 0.5 -1.055 0.0149996 0.55 -1.02717 0.0163077 0.6 -0.998156 0.0175022 0.65 -0.968266 0.0185757 0.7 -0.937796 0.0195232 0.75 -0.907015 0.0203419 0.8 -0.876165 0.0210315 0.85 -0.845458 0.0215938 0.9 -0.815079 0.0220327 0.95 -0.785183 0.0223537 1.0 -0.755899 0.0225636 例 6.6 用梯形公式计算积分 ∫ −= x t dtey 0 2 在 1,75.0,5.0=x 时的近似值(至少保留 6 位小数)。 解 依本题特点，可采用步长 25.0=h 。梯形公式为 )],(),([2 111 +++ ++= nnnnnn yxfyxf hyy 代入 2),( xeyxf −= ，有 )(2 2 1 2 1 +−−+ ++= nn xxnn eehyy 当 0)0( 0 == yy 时，计算得 242427.0)25.0( 1 == yy 457204.0)50.0( 2 == yy 625777.0)75.0( 3 == yy ．150 ． 实用数值分析解题指导 742985.0)00.1( 4 == yy 。 例 6.7 用欧拉法的外推公式计算 )5.00( 5.0)0( 2 ≤≤ ⎪⎩ ⎪⎨ ⎧ = += x y yxdx dy 方程的准确解： 225.22)( xxexy x −−+−= 。 解 外推公式： )(~)2(2 2 hy hyz nnn −= 计算结果见下表 nx ny~ ny2 nz nn yxy 2)( − nn zxy −)( 0.1 0.55 0.551375 0.55275 0.0015523 0.000177295 0.2 0.606 0.609541 0.613082 0.00396596 0.00042502 0.3 0.6706 0.677244 0.683888 0.00740314 0.000759252 0.4 0.74666 0.757511 0.768363 0.0120504 0.00119898 0.5 0.837326 0.853681 0.870037 0.0181219 0.00176658 可见外推值 nz 比 )2 (2 hy n 计算效果更好。 附 Mathematica 程序： }]5,1,{2"",1,"",,""],2[2,""],[1int[Pr ]];1.0[]2[2[1];]1.0[2];[1]2[22[ 2][5.22:][ }]10,1,{2]];1[2,[2]1[2][2[ ;0;05.02;5.0]0[2 }]6,1,{1]];1[1,[1]1[1][1[ ;0;1.01;5.0]0[1 ];0,10,2[2 ]0,6,1[1 :],[ 2 2 kzzzkyky kfkyAbszzkAbsfzkykyzDo xxxExpxf khxxkyxfhkykyDo xhy khxxkyxfhkykyDo xhy yArrayY yArrayY yxyxf −=−=−= −−+−= +=−+−= === +=−+−= === = = += − −− 例 6.8 用二阶 R-K 方法求解下列初值问题，并画出近似解草图。 )05.0,05.22.1( 2 .3)2.1( cos =≤≤ ⎪⎩ ⎪⎨ ⎧ = = hx y yxdx dy 解 计算公式 第六章 常微分方程数值解法 ．151． ⎪⎪⎩ ⎪⎪⎨ ⎧ ++= = ++= → ⎪⎪⎩ ⎪⎪⎨ ⎧ ++= = ++= ++ 12 1 211 12 1 211 )cos( cos )(2 ),( ),( )(2 hkyhxk yxk kkhyy hkyhxfk yxfk kkhyy nn nn nn nn nn nn 计算结果见下表。 n nx ny n nx ny 1 1.25 3.23038 2 1.3 3.25662 3 1.35 3.2786 4 1.4 3.29623 5 1.45 3.30943 6 1.5 3.31813 7 1.55 3.3223 8 1.6 3.32192 9 1.65 3.31699 10 1.7 3.30752 11 1.75 3.29358 12 1.8 3.27521 13 1.85 3.2525 14 1.9 3.22555 15 1.95 3.19449 16 2 3.15945 附 Mathematica 程序： f[x_,y_]:=Cos[x]Sqrt[y]; xylist={{1.2,3.2}};h=0.05; Do[xn=xylist[[n]][[1]];yn=xylist[[n]][[2]]; k1=f[xn,yn]; k2=f[xn+h,yn+h k1]; d=(k1+k2)h/2;Print[n," ",xn+h," ",yn+d]; xylist=Append[xylist,{xn+h,yn+d}],{n,1,17}] ListPlot[xylist,PlotJoined->True] 运行结果见下表： 1 1.25 3.23038 2 1.3 3.25662 3 1.35 3.2786 4 1.4 3.29623 5 1.45 3.30943 6 1.5 3.31813 7 1.55 3.3223 8 1.6 3.32192 ．152 ． 实用数值分析解题指导 9 1.65 3.31699 10 1.7 3.30752 11 1.75 3.29358 12 1.8 3.27521 13 1.85 3.2525 14 1.9 3.22555 15 1.95 3.19449 16 2. 3 .15945 17 2.05 3.12059 近似解草图如图： 例 6.9 由微分方程 ⎪⎩ ⎪⎨ ⎧ = += 0)0( y baxdx dy 导出欧拉公式，改进欧拉公式的表达式，并与准确解 bxaxxy += 22 1)( 比较。 解 (1) 欧拉公式： hxabxxa bxnnah bxkhah bhnxxxah hbaxhbaxhbaxy hbaxhbaxy hbaxyyxhfyy nnn n n n k n nnn nnnnnn 11 2 1 1 2 1 0 210 100 11 1 22 2 )1( )1()( )()()( )()( )(),( +++ + + = −− + −+= ++= += +++++= +++++++= = ++++= ++=+= ∑ L L L 第六章 常微分方程数值解法 ．153． 或 hxabxxay nnnn 22 2 −+= 对确定的 nxnhx =, ，有 xhayxy n 2)( =− 欧拉公式解收敛于准确解。 (2) 改进欧拉公式 2 1 )2 1( )2 1( )1()(2 2)(2)(2 )(2 )]()[(2 ))],(,(),([2 1 2 1 1 22 0 1 1 2 0 0 10 111 1 1 11 ++ + + = = + +−− + + ++ += +++= ++++= ++++= = +++++= ++= ++++= +++= ∑ ∑ nn n n n k n k kk nnnnn nnn nnn nnnnnnnn bxax bxnahy bxnkahy bhnxxahy bhxxahxxahy bhxxahy baxbaxhy yxfyxfyxfhyy L 或 nnn bxx ay += 22 对确定的 nxnhx =, ，有 bxxaxyyn +== 22)( 即差分是微分方程的准确解。事实上，改进欧拉公式就是二阶 R-K 方法，因此对方程解 bxxaxy += 22)( 是精确的。 例 6.10 试推导求解初值问题 00 )(),( yxyxyfy ==′ 的如下数值格式 ．154 ． 实用数值分析解题指导 ),2,1,0( )],()[,(2),( 2 1 L=+′++=+ nyxfxyyxfhyxhfyy nnnnnnnnnn 并说明它是多少阶格式。 解 泰勒展开是建立数值格式并推导其截断误差的重要方法。由泰勒展开有 32 1 !3 )( !2 )()()()( hyhxyhxyxyxy nnnnn ξ′′′+′′+′+=+ 因为 ))()(())((),( xyxfyxyfyxyxyfyxyfy +′=′+′=′′=′ 于是上式成为 3 2 1 !3 )())(()())(((2))(()()( h yxyxfxxyxyxfhxyxhfxyxy nnnnnnnnnnn ξ′′′++′++=+ 为考虑局部截断误差，设 )( nn xyy = ，于是 )(!3 )()( 33111 hOh yyxyR nnnn =ξ′′′=−= +++ 所以所给数值格式是二阶格式，它由上述式子中舍掉局部截断误差 3!3 )( hy nξ′′′ 而得。 例 6.11 用如下 4 步四阶的阿达姆斯显格式 24/)9375955( 3211 −−−+ −+−+= nnnnnn ffffhyy 求解初值问题 0.50 1)0( 23 ≤≤ ⎪⎩ ⎪⎨ ⎧ = −= x y yx dx dy 取步长 1.0=h ，小数点后至少保留六位。 解 可先用四阶单步法如四阶龙格-库塔法求 321 ,, yyy ，由于 yxyxf 23),( −= ，于 是 ⎪⎪ ⎪⎪ ⎩ ⎪⎪ ⎪⎪ ⎨ ⎧ ++=++++= +−=++= +−=++= +−=++= −== + 01405.0818733.02719.0)22(6 1 0273.0163.0245.0),( 0135.0182.0273.0)2,2( 015.018.027.0)2,2( 2.03.0),( 43211 34 2 3 1 2 1 nnnn nnnn nnnn nnnn nnnn yxKKKKyy yxKyhxhfK yxKyhxhfK yxKyhxhfK yxyxhfK 第六章 常微分方程数值解法 ．155． 由 1)0( 0 == yy ，计算得 660429.0)3.0( 273067.0)2.0( 832783.0)1.0( 3 2 1 =≈ =≈ =≈ yy yy yy 再由 4 步阿达姆斯显格式得 24/)]1827()74111()118177()110165[(1.0 3322111 −−−−−−+ +−+−++−+−+= nnnnnnnnnn yxyxyxyxyy 从而 636466.0 24/)]1827()74111()118177()110165[(1.0)4.0( 0011223334 = +−+−++−+−+=≈ yxyxyxyxyyy 643976.0 24/)]1827()74111()118177()110165[(1.0)5.0( 1122334445 = +−+−++−+−+=≈ yxyxyxyxyyy 三、基于 Mathematica 的数值计算实例 例 1 求解微分方程 xexy 23'' += ，并作出相应的积分曲线。 解 Mathematica 程序： DSolve[y''[x] == 2x + Exp[x], y[x], x] 运行得方程的通解为： 21 3 2 2 ) 4 21 2 ()( CxCxexxxy x +++++−= g1 = Table[ Plot[(-x/2+(1+2x)/4)E^(2x) + x^3/2 + c1 + x*c2, {x, -5, 5}, DisplayFunction -> Identity], {c1, -10, 10, 10}, {c2, -5, 5, 10}]; Show[g1, DisplayFunction -> $DisplayFunction] 积分曲线为： ．156 ． 实用数值分析解题指导 例 2 用 Euler 法求微分方程初值问题： ⎪⎩ ⎪⎨ ⎧ = += 1)0( 52' y y xyy 的数值解（步长h 取 0.1），求 解范围为区间 ]1,0[ 。 解 用 Euler 方法得算式 ⎪⎩ ⎪⎨ ⎧ = ++=+ 1)0( ) 5 2(1 y y xyhyy n n nnn Mathematica 程序为： f[x_, y_] := 2y + 5x/y; x = 0; y = 1; h = 0.1; While[x + h < 1, y = y + h*f[x, y]; x = x + h; Print["x=",x, " ", "y=", y]] 运行结果见下表： x=0.1 y=1.2 x=0.2 y=1.48167 x=0.3 y=1.84549 x=0.4 y=2.29587 x=0.5 y=2.84216 x=0.6 y=3.49855 x=0.7 y=4.28401 x=0.8 y=5.22251 x=0.9 y=6.3436 x=1. y=7.68326 第六章 常微分方程数值解法 ．157． 例 3 用改进 Euler 法求初值问题 ⎪⎩ ⎪⎨ ⎧ = = 2.3)2.1( sin3 2 y xy dx dy )31( ≤≤ x 的数值解，取 05.0=h 。 解 Mathematica 程序为： Clear[f,x,y] f[x_,y_]:=Sin[x^2]Sqrt[3y]; xyL={{1.1,3.0}};h=0.05; Do[xn=xyL[[n]][[1]];yn=xyL[[n]][[2]]; K1=f[xn,yn]; K2=f[xn+h,yn+h*K1]; d=(K1+K2)*h/2; xyL=Append[xyL,{xn+h,yn+d}],{n,1,17}] xyL//TableForm ListPlot[xyL,PlotStyle->PointSize[0.03]] ListPlot[xyL,PlotJoined->True] 运行结果： 1.1 3. 1.15 3.14455 1.2 3.29689 1.25 3.4553 1.3 3.61755 1.35 3.78087 1.4 3.94198 1.45 4.0971 1.5 4.24209 1.55 4.37254 1.6 4.48396 1.65 4.57206 1.7 4.63298 1.75 4.66358 1.8 4.66177 1.85 4.62671 1.9 4.55908 1.95 4.4611 ．158 ． 实用数值分析解题指导 近似图形为： 例 4 求微分方程组 ⎪⎩ ⎪⎨⎧ += = ttty tttx cos2)(' sin)(' 3 2 当 ⎩⎨ ⎧ = = 1)0( 1)0( y x 时的特解，并画出解的图形。 解 求特解的 Mathematic 程序如下： Clear[x,t,y] DSolve[{x'[t]==t^2*Sin[t],y'[t]==2t^3+Cos[t], x[0]==1,y[0]==1},{x[t],y[t]},t] 运行结果得特解： {{x[t] -> -1 + 2 Cos[t] – t2Cos[t] + 2 t Sin[t], y[t] -> 1 + 2 4t + Sin[t]}} 特解图形的 Mathematic 程序如下： x[t_]:=1+2Cos[t]-t^2*Cos[t]+2t*Sin[t]; y[t_]:=1+t^4/2+Sin[t]; ParametricPlot[{x[t],y[t]},{t,0,4Pi}] 运行结果得解的图形为： 第六章 常微分方程数值解法 ．159． 例 5 在区间 ]10,0[ 上，作出微分方程 ⎩⎨ ⎧ == =++ 0)0(,0)0(' sin2x5y3y''2y' 2 yy 的解的图形，以及在 2.3,1.1 == xx 处的数值解，精确到 510− 。 解 Mathematic 程序如下 NDSolve[{2y''[x]+3y'[x]+5y[x]==Sin[2x^2],y[0]==0,y'[0]==0},y,{x,0,10}, AccuracyGoal->5] Plot[Evaluate[y[x]/.%],{x,0,10},PlotRange->All] y[1.1]/.%% y[3.2]/.%%% 运行结果： (1) 图形 (2) 数值解 y[1.1]=0.064775； y[3.2]=-0.016441 例 6 画出微分方程 1)0(',2)0(,04)(')5)('()('' ===+++ yyytytyty 在 ]10,0[∈t 上的解的曲线，并求方程在 6.8,5.0=t 处的解。 解 Mathematic 程序如下： solution=NDSolve[{y''[t]+(y'[t]+5)y'[t]+4y[t]==0, ．160 ． 实用数值分析解题指导 y'[0]==1,y[0]==2},y[t],{t,0,10}] Plot[y[t]/.solution,{t,0,10},PlotRange->All] %%/.t->0.5//N %%%/.t->8.6//N 运行结果： (1) 解的曲线为 (2) 解为：{{y[0.5] -> 1.64131}}；{{y[8.6] -> 0.000208144}} 四、基础知识练习 6.1 用欧拉公式解初值问题： 1.00.1 0.1)1.0( 2 ≤≤ ⎪⎩ ⎪⎨ ⎧ = −= x y y xy dx dy 取 1.0=h ，并用外推方法求出具有 )( 2hO 精度的数值解。 6.2 用梯形公式： )],(),([2 111 +++ ++= nnnnnn yxfyxf hyy 解初值问题： 1.01.0,0 2.1)0( =≤≤ ⎪⎩ ⎪⎨ ⎧ = += hx y yxdx dy 取 第六章 常微分方程数值解法 ．161． 方程的准确解为： xyxxy ++= 25.02.1)( 6.3 用二阶 R-K 方法解常微分方程组： ⎪⎪⎩ ⎪⎪⎨ ⎧ = = +=′ −=′ 0)0( 1)0( cos)()( sin)()( y x xtxty xtytx 6.4 用二阶 R-K 方法求初值问题 ⎩⎨ ⎧ = −+=′ 0)0( 2 y yxxy 10 << x 的数值解，取 1.0=h 。 6.5 取 2.0=h ，用精典四阶 R-K 方法求初值问题 （1）⎩⎨ ⎧ = +=′ 1)0(y yxy （ 10 << x ） （2） ⎪⎩ ⎪⎨ ⎧ =+ =′ 1)0( 1 3 y x yy （ 10 << x ） 五、练习参考 6.1 取 0.1)1.0(,1.0 == yh ，欧拉计算公式 )2(1.01 n n nnn y xyyy −+=+ 外推公式： )()2(2 2 hy hyz nnn −= nx ny nz 0.2 1.08 1.07515 0.3 1.05096 1.14119 0.4 1.21393 1.19874 0.5 1.26942 1.24795 0.6 1.31759 1.28858 0.7 1.35827 1.32001 0.8 1.39102 1.34112 0.9 1.4151 1.35066 1.0 1.42941 1.34617 ．162 ． 实用数值分析解题指导 6.2 nx ny nx ny 0.1 0.2 0.3 0.4 0.5 1.33158 1.48753 1.67043 1.88311 2.1287 0.6 0.7 0.8 0.9 1.0 2.41067 2.73284 3.09946 3.51519 3.98521 6.3 n nt nx ny 1 2 3 4 5 6 7 8 9 10 0.1 0.979073 0.198254 0.2 0.984164 0.39293 0.3 0.999005 0.585653 0.4 1.02383 0.776467 0.5 1.05887 0.965471 0.6 1.10448 1.15288 0.7 1.16111 1.33895 0.8 1.22932 1.52408 0.9 1.30978 1.70882 1.0 1.40328 1.89383 6.4 nx : 0.1 0.3 0.5 ny :0.0055 0.058181 0.14500 6.5 （1）1.2428，1.5836， 2.0442， 2.6510 ，3.4365 （2）1.7276 ，2.7430 ,，4.0942 ,，5.8292 ， 7.9960