某些常微分方程數(shù)值解法與程序實現(xiàn)

1、 畢業(yè)設計說明書題 目：某些常微分方程數(shù)值解法與程序實現(xiàn)學 院： 年級、專業(yè)： 學 生： 學 號： 指導教師： 完成日期： 目 錄摘 要4引 言61 準備知識61.1 微分方程的解61.2 常微分方程的初值問題及其解的存在唯一性71.3 常微分方程的收斂性與穩(wěn)定性71.4 MATLAB開發(fā)工具簡介81.5 用戶圖形界面GUI制作92 常微分方程數(shù)值解92.1 Euler法92.2 改進Euler法102.3 Runge-Kutta法113 算法及程序流程圖153.1 Euler方法153.1.1 Euler法的算法設計153.1.2 Euler法的程序流程圖163.2 改進Euler法183.

2、2.1 改進Euler法的算法設計183.2.2 改進Euler法的程序流程圖193.3 Runge-Kutta方法203.3.1 Runge-Kutta法的算法設計203.3.2 Runge-Kutta法的程序流程圖224 基于MATLAB的程序界面設計254.1 用戶圖形界面GUI說明及詳細操作過程254.2.1 Euler向前法264.2.2 Euler向后法274.2.3 改進Euler法284.2.4 二階Runge-Kutta方法304.2.5 三階Runge-Kutta方法314.2.5 四階Runge-Kutta方法324.2.6 綜合比較344.2用戶圖形界面源程序355結論

3、41總結與體會42謝 辭42參考文獻43摘 要隨著社會的不斷發(fā)展與進步，微分方程逐漸在越來越多的現(xiàn)實生活中得到了非常廣泛的應用，但是，我們根據(jù)問題建立的微分方程模型中，只有非常少的模型能夠得到對應的解析解，因此，通過對常微分方程數(shù)值解的剖析，我們能夠更加清楚的認識和了解常微分方程。數(shù)值分析是通過研究、開發(fā)以及分析各種數(shù)學問題來得到數(shù)值解的算法。本文主要是闡述幾種常見的常微分方程數(shù)值解的方法，其中包括Euler向前法、Euler向后法、改進歐拉法、Runge-kutta法。介紹了這幾種方法的基本原理及推導過程，并且用matlab實現(xiàn)了這些方法，包括它們算法的思想，程序流程以及具體的實現(xiàn)。最后，利

4、用MATLAB中的圖形用戶界面GUI的功能來實現(xiàn)，讓用戶操作起來更加的得心應手。關鍵詞：Euler 改進歐拉 Runge-kutta MATLAB 圖形用戶界面AbstractWith the continuous development and progress of society, more and more gradually in the differential equation in real life has been applied widely, but according to the differential equation model established in

5、the model only very few can get the corresponding analytical solution, therefore, for we can know more about the differential equation, we must study the differential the number of equations solution.The numerical analysis is through research, development and analysis of various mathematical problem

6、s to obtain the numerical solution of the algorithm. This paper is mainly for numerical solutions of ordinary differential equations of several common, including the Euler Euler forward method and backward method, the improved Euler method, Runge-kutta method is introduced. The basic principle and t

7、he derivation of these methods, and these the method is implemented by MATLAB, including their algorithm, program flow and implementation. Finally, the use of MATLAB graphical user interface GUI function to achieve, allowing users to operate more handy.Keywords: Euler improved Euler Runge-kutta MATL

8、AB graphical user interface引 言微分方程現(xiàn)在在許多方面已經得到了普遍應用，許多實際問題所建立的數(shù)學模型也是微分方程的定解問題。在描述系統(tǒng)的動態(tài)演變時，如物種的增長和蛻變、物體的運動、電路的震動瞬間、化學反應過程等都能表示為以時間t為變量的常微分方程或方程組。在實際問題中的研究中，很多時候我們很難求出微分方程的解析解或者是求解過程過于繁瑣，因此我們的研究是非常的有意義。常微分方程是微分方程中只有一個自變量的函數(shù)，本文主要是對微分方程數(shù)值解的方法進行分析，并且利用MATLAB軟件將這些理論知識進行實現(xiàn)，既能得到微分方程的數(shù)值解，又能將這些理論知識運用到實際生活中。本文首

9、先對了歐拉向前算法、歐拉向后格式以及改進歐拉格式進行了分析和總結，歐拉算法相對于其他幾種方法來說更易于理解，但是歐拉向前算法的弱點在于其精度不高，所以為了提高歐拉算法的精度，提出了歐拉向后算法和改進歐拉算法。其次，本文介紹了龍格-庫塔法，其中包括二階、三階、四階龍格-庫塔法，龍格-庫塔法計算過程穩(wěn)定，程序也易實現(xiàn)，總的來說精度比歐拉算法的高。最后，本文給出了實例來證明了這幾種算法的是可行的，并且得到幾種算法是正確的。1 準備知識1.1 微分方程的解 一個函數(shù)如果代入微分方程后能使方程為恒等式，則稱為微分方程的解。所以，對于某個區(qū)間上的所有來說，能使微分方程恒成立，那么為常微分方程的解。1.2

10、常微分方程的初值問題及其解的存在唯一性討論常微分的初值問題： （1.2.1） 或記為 只有一些特殊形式的，才能找到它的解析解；對于大多數(shù)常微分方程的初值問題，只能計算它的數(shù)值解。在計算中約定表示常微分方程準確解的值，表示的近似值。由微分方程理論，為保證微分方程解的存在唯一性，通常對(1.2.1)中對應滿足Lipschitz條件，即對任意、存在常數(shù)，有 （1.2.2）1.3 常微分方程的收斂性與穩(wěn)定性(一)收斂性若對任意固定的，當 (同時)時，有數(shù)值解，則稱該常微分方程（1.2.1）是收斂的。常微分方程初值問題的形式還可以表示為：若單步法有階精度，且增量函數(shù)關于滿足Lipschitz條件，即存在

11、 （1.3.1）則其數(shù)值解的整體截斷誤差。（二）穩(wěn)定性由于計算數(shù)值解的過程中舍入誤差總會存在，因此我們需要討論其穩(wěn)定性。但由于差分方法數(shù)值穩(wěn)定性問題很難作一般性討論，一般我們都只用方程 （1.3.2）作討論，當計算某一步有舍入誤差時，但從此之后在計算過程中誤差不會逐步擴大，我們就稱這種穩(wěn)定性是絕對穩(wěn)定性。1.4 MATLAB開發(fā)工具簡介MATLAB是美國MathWorks公司出品的商業(yè)數(shù)學軟件，它集合了數(shù)值分析、矩陣運算、信號處理和圖形顯示于一體，構成了一個方便、有著良好界面的用戶環(huán)境，讓更多不同專業(yè)的人可以方便使用。Matlab是一個非常有用的數(shù)學軟件，它擁有，它不僅能夠解決我們數(shù)學中遇到的

12、各種問題，而且可以根據(jù)所求結果繪畫出相應的函數(shù)圖形，方便我們直接觀察。無論你身處何方，從事著什么職業(yè)，只要你遇到了數(shù)學中的難題，你都可以利用Matlab，它可以為你解決你所遇到的問題，讓你不用通過復雜的計算和推導來得到結果，真正實現(xiàn)了讓用戶用起來方便。MATLAB的基本單位是矩陣。它的表達式與數(shù)學、工程計算中常用的形式十分相似，MATLAB之所以能夠得到廣大用戶的喜愛，主要是因為它對于用戶來說非常實用。由于MATLAB的許多特點，比如：MATLAB操作起來方便，使初學者易于學習和掌握；它用有非常多的工具，能夠大量滿足我們的需求；它能夠實現(xiàn)程序的可視化，更用戶能夠更直接的使用和觀察等等?？偟膩碚f

13、，對于大部分人來說MATLAB比起其他語言更適用，特別是對于那些對語言不是特別精通的人來講。1.5 用戶圖形界面GUI制作圖形用戶界面是用于用戶和軟件相互溝通的一個界面，換句話說就是我們輸入一些數(shù)據(jù)，然后調用程序，通過程序內部運算得出結果，并且把這個得到的結果顯示在界面上，從而讓用戶更加直觀、清晰的看到結果，從而方便用戶的使用。用戶則通過界面上的相關內容的提示，在相對應的文本框中輸入數(shù)據(jù)，然后通過點擊相對應的按鈕，讓計算機產生某種特定的動作或變化，比如實現(xiàn)計算、繪圖等。本文這個設計是通過GUI來制作一個一階常微分方程數(shù)值解法的界面。這個界面由按鈕、文本、圖形等組成。其主界面由六個按鈕組成，每個

14、按鈕對應著一個相應算法程序，每個算法程序可以實現(xiàn)一種求常微分方程的差分方法，每個按鈕要求在相應的文本中輸入的內容，得到對應的圖形是不同的。2 常微分方程數(shù)值解2.1 Euler法Euler法的思想就是：在每一個小區(qū)間上，用一條切線來代替原函數(shù)曲線。從整體上看就是用一條折線來代替解（或積分）曲線，并以此來求取一系列離散結點處函數(shù)近似值。對于常微分方程初值問題（1.2.1），在求解區(qū)間上作等距剖分，步長，記=，。用數(shù)值差商的方法，即用差商近似導數(shù)求解常微分方程。由數(shù)值微分向前差商公式 得 式（1.2.1）實際上給出 故得 再由推導出向前的Euler公式： （2.1.1）同理，利用數(shù)值微分向后差商公

15、式，可以推導出向后的Euler公式： （2.1.2）式（2.1.2）是關于的非線性方程，它是隱式Euler公式。對于歐拉公式，對在作Taylor展開得到： （2.1.3）Euler公式（2.1.1）是由（2.1.3）展開式中截斷所得到，因此Euler公式的截斷誤差為，故Euler公式是一階方法。2.2 改進Euler法為了提高精度，從另一角度來對初值問題進行研究。對微分方程： 從到對兩邊進行積分，得到：即： （2.2.1）從（2.2.1）中，可以看出，只要近似地算出積分：，就可得到微分方程的另一種求解格式。為了減少誤差，由梯形公式： 可以得出： 進一步有： 再用，分別代替，就得到了新的求解格式

16、： （2.2.2）公式（2.2.2）是一個隱式方程，并且計算量大。于是我們往往先用歐拉向前公式（2.1.1）求出一個初步的近似值，記作，稱之為預報值，再將代替(2.2.2)右端的進行計算，得到校正值，于是就有： （2.2.3）方法（2.2.3）稱為改進的歐拉方法。改進歐拉法的截斷誤差為，歐拉方法每一步只要對調用一次，但是因為增加了校正過程，計算量比歐拉方法增加了一倍，付出這種代價是為了提高精度。2.3 Runge-Kutta法通過對Euler法和改進Euler法的對比，我們可發(fā)現(xiàn)Euler公式（2.1.1）僅取一個點的斜率作為平均斜率的近似值，計算精度較低。而改進Euler公式（2.2.1）則

17、是利用了兩個點和的斜率和的平均值作為平均斜率的近似值，即，其中是通過已知的值，然后再利用公式（2.1.1）所求得，此時計算精度提高到了二階。通過上面的分析，我們可以找出改進Euler法比Euler法計算更精確的原因。那就是在確定平均斜率時，多取了一個點的斜率值。因此，受到啟發(fā)，如果想辦法在區(qū)間上多找出幾個點的斜率值，然后把這幾個點的斜率值做加權平均，所得到的值作為的近似值，因此就有可能構造出更高精確度的數(shù)值計算公式，這就是Runge-Kutta法的基本思想。（一）二階Runge-Kutta法在區(qū)間上取兩點分別為、，以及該兩點處的斜率值和的加權平均(或稱線性組合)來求取平均斜率的近似值，即 其中

18、點的斜率值 點的斜率值 令：, (2.3.1)其中和為待定常數(shù)。對在處進行二階Taylor展開，則有 (2.3.2)對在處進行一階Taylor展開，得 (2.3.3)將式（2.3.3）代入式（2.3.1）中，得 （2.3.4）比較式（2.3.2）和式（2.3.4），得 (2.3.5)由于方程組（2.3.5）中有三個未知量，但只有兩個方程，方程的個數(shù)比未知量的個數(shù)少，因此它存在著無窮多組解。如果去。則有，用和分別代替和，此時計算公式（2.3.1）變?yōu)?這一多方程就是改進Euler公式。凡是滿足條件式（2.3.5）的計算格式都統(tǒng)稱為二階Runge-Kutta格式。改進Euler公式就是眾多的二階R

19、unge-Kutta公式中的一種特殊格式。若取，則，則得到另一種形式的二階Runge-Kutta法的計算公式： (2.3.6)其中點為區(qū)間的中點，所以公式（2.3.6）又稱為中點公式。（二）三階、四階Runge-Kutta法為了進一步提高數(shù)值計算的精度，在區(qū)間上除了和兩點之外，還需要再加一點，以及在這三點處的斜率值和的加權平均值來求平均斜率的近似值，此時得到一個計算公式為 (2.3.7)其中和為待定常數(shù)。為了求出點的斜率，必須要先找出點處在位置對應的函數(shù)值，我們在區(qū)間上可以利用二階Runge-Kutta公式，即得到 其中和為待定常數(shù)。因此就有 于是得到一個計算公式為 (2.3.8)把在處進行二

20、階Taylor展開，并代入（2.3.8）中的第一個方程式，在對在處進行三階Taylor展開，并且比較兩式的系數(shù)，就可以得到參數(shù)滿足的方程組 (2.3.9)方程組（2.3.9）中總共有5個方程以及7個未知量，方程個數(shù)少于未知量的個數(shù)，因此解不唯一。常見的三階Runge-Kutta公式為 (2.3.10)如果還需要再進一步提高計算精度到四階，用類似上述的構造方法就可以推導出四階Runge-Kutta法。但因為推導過程復雜，在這里就不進行推導了，四階Runge-Kutta法也不止一個，下面給出常用的四階方法中最常用的一個計算公式： (2.3.11)3 算法及程序流程圖3.1 Euler方法1 變量說

21、明表1 Euler方法的符號說明變量相應的說明的初始條件的初始條件迭代精度區(qū)間的左端點區(qū)間的右端點區(qū)間的步長迭代次數(shù)第步的計算結果第步的計算結果3.1.1 Euler法的算法設計歐拉向前算法3.1（1）輸入：，函數(shù)，步長，初值（2）賦值：（3）計算 輸出（4）如果，則，轉3，否則停止運行歐拉向后算法3.2（1）輸入，函數(shù)，步長，初值，精度（2）賦值：（3）計算 其中用迭代法求出，輸出（4）如果且，則，轉3，否則停止運行3.1.2 Euler法的程序流程圖圖3-1是Euler向前方法的程序流程，圖3-2是Euler向后方法的程序流程。 圖3-1 圖3-23.2 改進Euler法1 符號說明表2

22、改進Euler方法的符號說明變量相應的說明的初始條件的初始條件區(qū)間的左端點區(qū)間的右端點區(qū)間的步長迭代次數(shù)第步的計算結果第步的計算結果3.2.1 改進Euler法的算法設計為了方便編程設計，我們通常把（2.2.3）改為 算法3.3（1）輸入：區(qū)間，函數(shù)，步長，初值（2）賦值：（3）計算 輸出（4）如果，則，轉3；否則停止運行3.2.2 改進Euler法的程序流程圖圖3-3是改進Euler方法的程序流程。 圖3-33.3 Runge-Kutta方法1 符號說明表3 Runge-Kutta方法的符號說明變量相應的說明的初始條件的初始條件區(qū)間的左端點區(qū)間的右端點區(qū)間的步長迭代次數(shù)第步的計算結果第步的計

23、算結果3.3.1 Runge-Kutta法的算法設計二階Runge-Kutta算法3.4（1）輸入：，函數(shù)，步長，初值（2）賦值：（3）計算 輸出（4）如果，則，轉3；否則停止.三階Runge-Kutta算法3.5（1）輸入：，函數(shù)，步長，初值（2）賦值：（3）計算 輸出（4）如果，則，轉3；否則停止四階階Runge-Kutta算法3.6（1）輸入：，函數(shù)，步長，初值（2）賦值：（3）計算 輸出（4）如果，則，轉3；否則停止3.3.2 Runge-Kutta法的程序流程圖圖3-4是二階Runge-Kutta方法的程序流程，圖3-5是三階Runge-Kutta方法的程序流程，圖3-6是四階Run

24、ge-Kutta方法的程序流程。 圖3-4 圖3-5 圖3-64 基于MATLAB的程序界面設計4.1 用戶圖形界面GUI說明及詳細操作過程隨著現(xiàn)在軟件的更新，大部分程序都是用圖形界面來實現(xiàn)。圖形用戶界面利用窗口、按鈕、文本框、文字等等對象來構成一個界面，大大方便了用戶的使用。根據(jù)題目的需求，在用戶圖形界面在設計了幾個功能按鈕，每一個按鈕對應一種微分方程數(shù)值解法，并且在界面顯示了不同解法與精確值的比較圖。詳細界面設計如下：主用戶圖形界面由七個按鈕，一個文本組成，每一個按鈕又對應一個界面，其界面由兩個按鈕、五個需要輸入數(shù)據(jù)的文本框和一個控制圖形的控件組成。在執(zhí)行每種不同的數(shù)值解法時，都需要輸入指

25、定類型的數(shù)據(jù)，包括常微分方程、的取值范圍、初始值以及步長，當輸入數(shù)據(jù)之后，點擊確認，可以求出每步迭代之后的值和對應的值，根據(jù)對應的數(shù)據(jù)生成對應的圖形，通過圖形我們可以觀察到數(shù)值解與精確值的差別。開始運行程序時，我們直接進入的是如下用戶圖形開始操作界面：圖4-1在用戶圖形界面中，我們可以選擇不同的常微分方程的數(shù)值解法，根據(jù)我們選擇不同的解法，界面將跳轉到相應的界面。4.2.1 Euler向前法當點擊Euler向前這個按鈕時，程序界面將跳轉到如圖4-2這個界面，在Euler向前這個界面中，我們看到我們需要輸入一些相關的數(shù)據(jù)，如：需要求的函數(shù)，的取值范圍，的初始值，以及所取的步長。圖4-2然后，我們

26、求的數(shù)值解，的取值范圍為，步長=0.01，將程序得到的數(shù)值解與精確解對比得出圖形，結果為圖4-3，當步長取值較大時，誤差較大。由對比圖可以看出精確解與數(shù)值解重疊效果比較好，說明其誤差小，結果是正確的。圖4-34.2.2 Euler向后法當點擊Euler向后這個按鈕時，程序界面將跳轉到如圖4-4這個界面，在Euler向后這個界面中，我們看到我們需要輸入一些相關的數(shù)據(jù)，如：需要求的函數(shù)，的取值范圍，的初始值，以及所取的步長和代數(shù)精度。圖4-4然后，我們求的數(shù)值解，的取值范圍為，步長=0.01，代數(shù)精度為0.001，將程序得到的數(shù)值解與精確解對比得出圖形，結果為圖4-5，當步長取值較大時，誤差較大。

27、由對比圖可以看出精確解與數(shù)值解重疊效果比較好，說明其誤差小，結果是正確的。圖4-54.2.3 改進Euler法當點改進Euler這個按鈕時，程序界面將跳轉到如圖4-6這個界面，在改進Euler這個界面中，我們看到我們需要輸入一些相關的數(shù)據(jù)，如：需要求的函數(shù)，的取值范圍，的初始值，以及所取的步長。圖4-6然后，我們求的數(shù)值解，的取值范圍為，步長=0.1，將程序得到的數(shù)值解與精確解對比得出圖形，結果為圖4-7。由對比圖可以看出精確解與數(shù)值解重疊效果比較好，說明其誤差小，結果是正確的。圖4-74.2.4 二階Runge-Kutta方法當點擊二階Runge-Kutta這個按鈕時，程序界面將跳轉到如圖4

28、-8這個界面，在二階Runge-Kutta這個界面中，我們看到我們需要輸入一些相關的數(shù)據(jù)，如：需要求的函數(shù)，的取值范圍，的初始值，以及所取的步長。圖4-8然后，我們求的數(shù)值解，的取值范圍為，步長=0.1，將程序得到的數(shù)值解與精確解對比得出圖形，結果為圖4-9。由對比圖可以看出精確解與數(shù)值解重疊效果比較好，說明其誤差小，結果是正確的。圖4-94.2.5 三階Runge-Kutta方法當點擊三階Runge-Kutta這個按鈕時，程序界面將跳轉到如圖4-10這個界面，在三階Runge-Kutta這個界面中，我們看到我們需要輸入一些相關的數(shù)據(jù)，如：需要求的函數(shù)，的取值范圍，的初始值，以及所取的步長。圖

29、4-10然后，我們求的數(shù)值解，的取值范圍為，步長=0.1，將程序得到的數(shù)值解與精確解對比得出圖形，結果為圖4-11。由對比圖可以看出精確解與數(shù)值解重疊效果比較好，說明其誤差小，結果是正確的。圖4-114.2.5 四階Runge-Kutta方法當點擊四階Runge-Kutta這個按鈕時，程序界面將跳轉到如圖4-12這個界面，在四階Runge-Kutta這個界面中，我們看到我們需要輸入一些相關的數(shù)據(jù)，如：需要求的函數(shù)，的取值范圍，的初始值，以及所取的步長。圖4-12然后，我們求的數(shù)值解，的取值范圍為，步長=0.1，將程序得到的數(shù)值解與精確解對比得出圖形，結果為圖4-13。由對比圖可以看出精確解與數(shù)

30、值解重疊效果比較好，說明其誤差小，結果是正確的。圖4-134.2.6 綜合比較當點擊綜合比較這個按鈕時，程序界面將跳轉到如圖4-14這個界面，在綜合比較這個界面中，我們看到我們需要輸入一些相關的數(shù)據(jù)，如：需要求的函數(shù)，的取值范圍，的初始值，以及所取的步長。圖4-14然后，我們求的數(shù)值解，的取值范圍為，步長=0.1，將程序得到的幾種方法的數(shù)值解與精確解的圖形，結果為圖4-15。圖4-154.2用戶圖形界面源程序這幾種常微分數(shù)值解在用戶圖形界面中實現(xiàn)的主要源程序代碼如下：（一）Euler向前法：function pushbutton1_Callback(hObject, eventdata, ha

31、ndles)fun=inline(get(handles.edit1,String);xmin=str2num(get(handles.edit2,String);xmax=str2num(get(handles.edit3,String);y0=str2num(get(handles.edit4,String);h=str2num(get(handles.edit5,String);x=xmin:h:xmax;N=(xmax-xmin)/h;y(1)=y0;for n=1:N y(n+1)=y(n)+h*feval(fun,x(n),y(n);endset(handles.text13,St

32、ring,num2str(x);set(handles.text14,String,num2str(y);x=x;y=y;y1=dsolve(Dy=-y+x+1,y(0)=1,x);ezplot(y1,xmin,xmax);hold on;plot(x,y,-r);legend(解析解,數(shù)值解);hold off;% hObject handle to pushbutton1 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles structure with handles a

33、nd user data (see GUIDATA)（二）Euler向后法：function pushbutton2_Callback(hObject, eventdata, handles)fun=inline(get(handles.edit1,String);xmin=str2num(get(handles.edit2,String);xmax=str2num(get(handles.edit3,String);y0=str2num(get(handles.edit4,String);h=str2num(get(handles.edit5,String);a=str2num(get(ha

34、ndles.edit6,String);x=xmin:h:xmax;N=(xmax-xmin)/h;y(1)=y0;k=1;for n=1:N z=y(n)+h*feval(fun,x(n),y(n); while kN z1=y(n)+h*feval(fun,x(n+1),z); if(abs(z1-z)a break; end z=z1; end y(n+1)=z1;endset(handles.text9,String,num2str(x);set(handles.text10,String,num2str(y);x=x;y=y;y1=dsolve(Dy=-y+x+1,y(0)=1,x)

35、;ezplot(y1,xmin,xmax);hold on;plot(x,y,-r);legend(解析解,數(shù)值解);hold off;% hObject handle to pushbutton2 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles structure with handles and user data (see GUIDATA（三）改進Euler法：function pushbutton2_Callback(hObject, eventdata, han

36、dles)fun=inline(get(handles.edit1,String);xmin=str2num(get(handles.edit2,String);xmax=str2num(get(handles.edit3,String);y0=str2num(get(handles.edit4,String);h=str2num(get(handles.edit5,String);N=(xmax-xmin)/h;y(1)=y0;x(1)=xmin;x(N)=0;y(N)=0;for i=1:N x(i+1)=xmin+i*h; y1=y(i)+h*feval(fun,x(i),y(i); y

37、2=y(i)+h*feval(fun,x(i+1),y1); y(i+1)=(y1+y2)/2;endset(handles.text8,String,num2str(x);set(handles.text9,String,num2str(y);y1=dsolve(Dy=-y+x+1,y(0)=1,x);ezplot(y1,xmin,xmax);hold on;plot(x,y,-r);legend(解析解,數(shù)值解);hold off;% hObject handle to pushbutton2 (see GCBO)% eventdata reserved - to be defined i

38、n a future version of MATLAB% handles structure with handles and user data (see GUIDATA)（四）二階Runge-Kutta方法：function pushbutton2_Callback(hObject, eventdata, handles)fun=inline(get(handles.edit1,String);xmin=str2num(get(handles.edit2,String);xmax=str2num(get(handles.edit3,String);y0=str2num(get(handl

39、es.edit4,String);h=str2num(get(handles.edit5,String);x=xmin:h:xmax;N=(xmax-xmin)/h;y(1)=y0;for n=1:N k1=feval(fun,x(n),y(n); k2=feval(fun,x(n)+h/2,y(n)+h/2*k1); y(n+1)=y(n)+h*k2;endset(handles.text8,String,num2str(x);set(handles.text9,String,num2str(y);x=x;y=y;y1=dsolve(Dy=-y+x+1,y(0)=1,x);ezplot(y1

40、,xmin,xmax);hold on;plot(x,y,-r);legend(解析解,數(shù)值解);hold off; % hObject handle to pushbutton2 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles structure with handles and user data (see GUIDATA)（五）三階Runge-Kutta方法：function pushbutton2_Callback(hObject, eventdata, hand

41、les)fun=inline(get(handles.edit1,String);xmin=str2num(get(handles.edit2,String);xmax=str2num(get(handles.edit3,String);y0=str2num(get(handles.edit4,String);h=str2num(get(handles.edit5,String);x=xmin:h:xmax;N=(xmax-xmin)/h;y(1)=y0;for n=1:N k1=feval(fun,x(n),y(n); k2=feval(fun,x(n)+h/2,y(n)+h/2*k1);

42、k3=feval(fun,x(n)+h,y(n)-h*k1+2*h*k2); y(n+1)=y(n)+h/6*(k1+4*k2+k3);endset(handles.text10,String,num2str(x);set(handles.text11,String,num2str(y);y1=dsolve(Dy=-y+x+1,y(0)=1,x);ezplot(y1,xmin,xmax);hold on;x=x;y=y;plot(x,y,-r);legend(解析解,數(shù)值解);hold off;% hObject handle to pushbutton2 (see GCBO)% eventd

43、ata reserved - to be defined in a future version of MATLAB% handles structure with handles and user data (see GUIDATA)（六）四階Runge-Kutta方法：function pushbutton2_Callback(hObject, eventdata, handles)fun=inline(get(handles.edit1,String);xmin=str2num(get(handles.edit2,String);xmax=str2num(get(handles.edit

44、3,String);y0=str2num(get(handles.edit4,String);h=str2num(get(handles.edit5,String);x=xmin:h:xmax;N=(xmax-xmin)/h;y(1)=y0;for n=1:N 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)+h,y(n)+h*k3); y(n+1)=y(n)+h/6*(k1+2*k2+2*k3+k4);endset(

45、handles.text8,String,num2str(x);set(handles.text9,String,num2str(y);x=x;y=y;y1=dsolve(Dy=-y+x+1,y(0)=1,x);ezplot(y1,xmin,xmax);hold on;plot(x,y,-r);legend(解析解,數(shù)值解);hold off;% hObject handle to pushbutton2 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles structure with handles and user data (see GUIDATA)5結論本次的設計主要對一階常微分方程數(shù)值解的基本思想及數(shù)學算法原理進行了描述，分別對Euler法、改進Euler法、Runge-kutta法的基本數(shù)學方法和計算公式進行了分析。本次設計還通過MATLAB對Euler法、改進Euler法、Runge-kutta法進行了程序的實現(xiàn)，包括算法思想、程序框圖及具體的程序。更是利用MATLAB的用戶圖形界面GUI實現(xiàn)了這幾種方法的具體實現(xiàn)，讓用戶擁有了一個方便、