常微分方程數(shù)值求解問題的實(shí)習(xí)報告

1、 計(jì)算機(jī)與信息學(xué)院（數(shù)學(xué)類課程）課程實(shí)習(xí)報告課程名稱：常微分方程課程實(shí)習(xí)實(shí)習(xí)題目：常微分方程數(shù)值求解問題的實(shí)習(xí)姓 名：系：信息與計(jì)算科學(xué)專 業(yè)：信息與計(jì)算科學(xué)年 級：學(xué) 號：指導(dǎo)教師：職 稱：講師2010 年 12 月 1 日目 錄1.實(shí)習(xí)的目的和任務(wù)12.實(shí)習(xí)要求13.實(shí)習(xí)地點(diǎn)14.主要儀器設(shè)備15.實(shí)習(xí)內(nèi)容15.1 用不同格式對同一個初值問題的數(shù)值求解及其分析.15.1.1求精確解15.1.2用歐拉法求解45.1.3用改進(jìn)歐拉法求解65.1.4用4級4階龍格庫塔法求解8 5.1.5 問題討論與分析11 5.2 一個算法不同不長求解同一個初值問題及其分析. 146.結(jié)束語15參考文獻(xiàn)15常微

2、分方程課程實(shí)習(xí)1. 實(shí)習(xí)的目的和任務(wù)目的：通過課程實(shí)習(xí)能夠應(yīng)用matlab軟來計(jì)算微分方程（組）的數(shù)值解；了解常微分方程數(shù)值解。任務(wù)：通過具體的問題，利用matlab軟件來計(jì)算問題的結(jié)果，分析問題的結(jié)論。2. 實(shí)習(xí)要求能夠從案例的自然語言描述中，抽象出其中的數(shù)學(xué)模型；能夠熟練應(yīng)用所學(xué)的數(shù)值解計(jì)算方法；能夠熟練使用matlab軟件；對常微分方程數(shù)值解有所認(rèn)識，包括對不同算法有所認(rèn)識和對步長有所認(rèn)識。3. 實(shí)習(xí)地點(diǎn)數(shù)學(xué)實(shí)驗(yàn)室、學(xué)生宿舍4. 主要儀器設(shè)備計(jì)算機(jī)、microsoft windows xpmatlab 6.55. 實(shí)習(xí)內(nèi)容5.1 用歐拉方法，改進(jìn)歐拉方法，4階龍格庫塔方法分別求下面微分方

3、程的初值dy/dx= y(0)=1 x0，65.1 .1求精確解首先可以求得其精確解為：y=sinx+1/(x+1)程序代碼：x=0:0.2:6; y= sinx+1./(x+1) plot(x,y,'b*-');data=x',y'>> aay = columns 1 through 9 1.0000 1.0320 1.1037 1.1896 1.2729 1.3415 1.3866 1.4021 1.3842 columns 10 through 18 1.3310 1.2426 1.1210 0.9696 0.7933 0.5981 0.391

4、1 0.1797 -0.0283 columns 19 through 27 -0.2251 -0.4035 -0.5568 -0.6793 -0.7664 -0.8151 -0.8238 -0.7923 -0.7222 columns 28 through 31 -0.6165 -0.4798 -0.3175 -0.1366data = 0 1.0000 0.2000 1.0320 0.4000 1.1037 0.6000 1.1896 0.8000 1.2729 1.0000 1.3415 1.2000 1.3866 1.4000 1.4021 1.6000 1.3842 1.8000 1

5、.3310 2.0000 1.2426 2.2000 1.1210 2.4000 0.9696 2.6000 0.7933 2.8000 0.5981 3.0000 0.3911 3.2000 0.1797 3.4000 -0.0283 3.6000 -0.2251 3.8000 -0.4035 4.0000 -0.5568 4.2000 -0.6793 4.4000 -0.7664 4.6000 -0.8151 4.8000 -0.8238 5.0000 -0.7923 5.2000 -0.7222 5.4000 -0.6165 5.6000 -0.4798 5.8000 -0.31756.

6、0000 -0.1366 5.1.2 用歐拉法求解程序如下：建立函數(shù)文件cwfa1.mfunction x,y=cwfa1(fun,x_span,y0,h)x=x_span(1):h:x_span(2);y(1)=y0;for n=1:length(x)-1 y(n+1)=y(n)+h*feval(fun,x(n),y(n);endx=x'y=y'在matlab輸入以下程序：fun=inline('-y2+2*y*sin(x)+cos(x)-sin(x)2');x,y=cwfa1(fun,0,6,1,0.2); x,yplot(x,y,'r*-'

7、)結(jié)果及其圖象ans = 0 1.0000 0.2000 1.0000 0.4000 1.0676 0.6000 1.1598 0.8000 1.2540 1.0000 1.3358 1.2000 1.3950 1.4000 1.4246 1.6000 1.4200 1.8000 1.3788 2.0000 1.3006 2.2000 1.1867 2.4000 1.0404 2.6000 0.8663 2.8000 0.6703 3.0000 0.4594 3.2000 0.2411 3.4000 0.0235 3.6000 -0.1854 3.8000 -0.3780 4.0000 -0.

8、5471 4.2000 -0.6866 4.4000 -0.7915 4.6000 -0.8581 4.8000 -0.8842 5.0000 -0.8692 5.2000 -0.8141 5.4000 -0.7214 5.6000 -0.5950 5.8000 -0.4401 6.0000 -0.2631 5.1.3用改進(jìn)歐拉法求解：程序如下：建立函數(shù)文件cwfa2.mfunction x,y=cwfa2(fun,x_span,y0,h)x=x_span(1):h:x_span(2);y(1)=y0;for n=1:length(x)-1 k1=feval(fun,x(n),y(n); y(

9、n+1)=y(n)+h*k1; k2=feval(fun,x(n+1),y(n+1); y(n+1)=y(n)+h*(k1+k2)/2;endx=x'y=y'在matlab輸入以下程序： fun=inline('-y2+2*y*sin(x)+cos(x)-sin(x)2'); x,y=cwfa2(fun,0,6,1,0.2); x,y plot(x,y,'r*-')結(jié)果及其圖象：ans = 0 1.0000 0.2000 1.0338 0.4000 1.1050 0.6000 1.1898 0.8000 1.2716 1.0000 1.3388

10、1.2000 1.3827 1.4000 1.3972 1.6000 1.3785 1.8000 1.3248 2.0000 1.2362 2.2000 1.1146 2.4000 0.9635 2.6000 0.7877 2.8000 0.5933 3.0000 0.3871 3.2000 0.1767 3.4000 -0.0302 3.6000 -0.2260 3.8000 -0.4034 4.0000 -0.5557 4.2000 -0.6772 4.4000 -0.7636 4.6000 -0.8117 4.8000 -0.8199 5.0000 -0.7882 5.2000 -0.7

11、180 5.4000 -0.6125 5.6000 -0.4760 5.8000 -0.31426.0000 -0.1338 5.1.4 用4階龍格庫塔求解程序如下：建立函數(shù)文件cwfa3.mfunction x,y=cwfa3(fun,x_span,y0,h)x=x_span(1):h:x_span(2);y(1)=y0;for n=1:length(x)-1 k1=feval(fun,x(n),y(n); k2=feval(fun,x(n)+h/2,y(n)+h/2*k1); k3=feval(fun,x(n)+h/2,y(n)+h/2*k2); k4=feval(fun,x(n+1),y

12、(n)+h*k3); y(n+1)=y(n)+h*(k1+2*k2+2*k3+k4)/6;endx=x'y=y'在matlab輸入以下程序：fun=inline('-y2+2*y*sin(x)+cos(x)-sin(x)2'); x,y=cwfa3(fun,0,6,1,0.2); x,y plot(x,y, 'r*-')結(jié)果及其圖象： 0 1.0000 0.2000 1.0320 0.4000 1.1037 0.6000 1.1896 0.8000 1.2729 1.0000 1.3415 1.2000 1.3866 1.4000 1.4021

13、1.6000 1.3842 1.8000 1.3310 2.0000 1.2426 2.2000 1.1210 2.4000 0.9696 2.6000 0.7933 2.8000 0.5981 3.0000 0.3911 3.2000 0.1797 3.4000 -0.0283 3.6000 -0.2251 3.8000 -0.4035 4.0000 -0.5568 4.2000 -0.6793 4.4000 -0.7664 4.6000 -0.8151 4.8000 -0.8238 5.0000 -0.7923 5.2000 -0.7222 5.4000 -0.6165 5.6000 -0

14、.4798 5.8000 -0.31766.0000 -0.1366 5.1.5 問題討論與分析由以上數(shù)值分析結(jié)果繪制表格：精確解歐拉方法改進(jìn)的歐拉方法四階龍格-庫塔方法xiyiyi誤差yi誤差yi誤差011010100.21.03210.0321.03380.00181.03200.41.10371.06760.03611.1050.00131.103700.61.18961.15980.02981.18980.00021.189600.81.27291.2540.01891.27160.00131.2729011.34151.33580.00571.33880.00271.341501.2

15、1.38661.3950.00841.38270.00391.386601.41.40211.42460.02251.39720.00491.402101.61.38421.420.03581.37850.00571.384201.81.3311.37880.04781.32480.00621.331021.24261.30060.0581.23620.00641.242602.21.1211.18670.06571.11460.00641.12102.40.969581.04040.070820.963490.006090.969571e-052.60.793280.86630.073020

16、.787710.005570.793271e-052.80.598150.670310.072160.593280.004870.598141e-0530.391120.459370.068250.387120.0040.391111e-053.20.179720.241120.06140.17670.003020.179711e-053.4-0.028270.0235210.051789-0.030240.001976-0.028288e-063.6-0.22513-0.185410.03972-0.226030.0009-0.225141e-053.8-0.40352-0.377990.0

17、2553-0.403380.00014-0.403531e-054-0.5568-0.547120.00968-0.555670.00113-0.556811e-054.2-0.67927-0.686640.00737-0.677240.00203-0.679281e-054.4-0.76642-0.791530.02511-0.763620.0028-0.766431e-054.6-0.81512-0.858120.043-0.811710.00341-0.815131e-054.8-0.82375-0.884230.06048-0.81990.00385-0.823772e-055-0.7

18、9226-0.869240.07698-0.788150.00411-0.792282e-055.2-0.72216-0.814110.09195-0.717990.00417-0.722193e-055.4-0.61651-0.721370.10486-0.612450.00406-0.616543e-055.6-0.47975-0.594960.11521-0.475970.00378-0.479772e-055.8-0.31754-0.440110.12257-0.314210.00333-0.317573e-056-0.13656-0.263130.12657-0.13380.0027

19、6-0.136582e-05x=0,0.2,0.4,0.6,0.8,1,1.2,1.4,1.6,1.8,2,2.2,2.4,2.6,2.8,3,3.2,3.4,3.6,3.8,4,4.2,4.4,4.6,4.8,5,5.2,5.4,5.6,5.8,6;y1=1,1.032,1.1037,1.1896,1.2729,1.3415,1.3866,1.4021,1.3842,1.331,1.2426,1.121,0.96958,0.79328,0.59815,0.39112,0.17972,-0.028268,-0.22513,-0.40352,-0.5568,-0.67927,-0.76642,-

20、0.81512,-0.82375,-0.79226,-0.72216,-0.61651,-0.47975,-0.31754,-0.13656;y2=1,1,1.0676,1.1598,1.254,1.3358,1.395,1.4246,1.42,1.3788,1.3006,1.1867,1.0404,0.8663,0.67031,0.45937,0.24112,0.023521,-0.18541,-0.37799,-0.54712,-0.68664,-0.79153,-0.85812,-0.88423,-0.86924,-0.81411,-0.72137,-0.59496,-0.44011,-

21、0.26313;y3=1,1.0338,1.105,1.1898,1.2716,1.3388,1.3827,1.3972,1.3785,1.3248,1.2362,1.1146,0.96349,0.78771,0.59328,0.38712,0.1767,-0.030244,-0.22603,-0.40338,-0.55567,-0.67724,-0.76362,-0.81171,-0.8199,-0.78815,-0.71799,-0.61245,-0.47597,-0.31421,-0.1338;y4=1,1.032,1.1037,1.1896,1.2729,1.3415,1.3866,1

22、.4021,1.3842,1.331,1.2426,1.121,0.96957,0.79327,0.59814,0.39111,0.17971,-0.028276,-0.22514,-0.40353,-0.55681,-0.67928,-0.76643,-0.81513,-0.82377,-0.79228,-0.72219,-0.61654,-0.47977,-0.31757,-0.13658;plot(x,y1,'r+-')hold on,plot(x,y2,'b-') title('精確解與歐拉方法比較')plot(x,y1,'r+ -') hold on,plot(x,y3,'b-') title('精確解與改進(jìn)歐拉方法比較')plot(x,y1,'r+-')hold on,plot(x,y4,'b-')title('精確解與4階龍格庫塔法比較')由上表和圖可以看出歐拉法誤差最大，而改進(jìn)歐拉和龍格庫塔方法誤差相對較小，并且龍格庫塔方法誤差最小且大部