淺談對矩陣理論、數(shù)值分析及Matlab應(yīng)用的理解 4

1、3數(shù)值分析數(shù)值分析是研究用計算機求解各種數(shù)學問題的數(shù)值方法及其理論的一門學科。數(shù)值分析也稱為數(shù)值計算方法、科學計算。隨著計算科學與技術(shù)的進步和發(fā)展，科學與工程計算的應(yīng)用范圍已擴大到許多的學科領(lǐng)域，新的、有效的數(shù)值方法不斷出現(xiàn)，形成了許多新型交叉學科。如：計算金融學、計算物理學、計算力學、計算化學等。繼理論方法和實驗方法之后，數(shù)值計算方法已成為科學研究的第三種基本手段，科學計算已經(jīng)成為了人們進行科學活動必不可少的科學方法和工具。所以，數(shù)值分析既是一個基礎(chǔ)性的，同時也是一個應(yīng)用性的數(shù)學學科，與其他學科的聯(lián)系十分緊密。 在計算機問世前, 由于受到計算工具的限制, 許多很好的數(shù)值方法難以應(yīng)用。計算機出

2、現(xiàn)后, 其在通用性、標準化、大容量、高精度等方面都取得了驚人的成就, 特別是近20 年高速巨型計算機的成功研制與發(fā)展, 為數(shù)值分析方法提供了極好的機會, 大大提高了作為計算方法理論基礎(chǔ)的數(shù)值分析方法課程的重要性與實用性。利用數(shù)值計算方法和計算機來解決科學與工程中的問題，通常稱為科學與工程計算，或簡稱科學計算?？茖W計算是人們在社會實踐中，特別是在科學和工程實踐中普遍面臨的問題：原問題 抽象為數(shù)學問題（建立原問題的數(shù)學模型） 求解該數(shù)學問題（必然涉及數(shù)值計算） 回答原問題。由于科學計算本身的迅速發(fā)展及不斷取得成效，使得科學計算與傳統(tǒng)的理論研究和實驗研究并列稱為當今研究的三大方法。而科學計算與具體科

3、學的交叉發(fā)展，又形成了諸如計算力學、計算物理、計算化學、計算生物學等新的邊緣科學。這些學科與“數(shù)值分析”息息相關(guān)。誤差來源誤差是實驗科學術(shù)語。指測量結(jié)果偏離真值的程度。對任何一個物理量進行的測量都不可能得出一個絕對準確的數(shù)值，即用測量技術(shù)所能達到的最完善的方法，測出的數(shù)值也和真實值存在差異，這種測量值和真實值的差異稱為誤差。分為絕對誤差和相對誤差。也可以根據(jù)誤差的來源分為系統(tǒng)誤差（又稱偏性）和隨機誤差（又稱機會誤差）計量或測定中的誤差是指測定結(jié)果與真實結(jié)果之間的差值。 數(shù)值計算普遍存在于科學研究和工程應(yīng)用中，由于誤差的存在，一般難以獲得精確的計算結(jié)果，產(chǎn)生誤差的原因主要有以下幾個方面：1.模型

4、誤差：數(shù)學模型對實際問題的僅是刻畫；基于對實際問題近似描述的數(shù)學模型進行數(shù)值計算，例如利用函數(shù)的n階Taylor展式計算函數(shù)值；2.觀測誤差：數(shù)學模型或計算公式中通常包含若干參數(shù)，這些參數(shù)往往是通過觀測或?qū)嶒灥玫降?，這樣得到的參數(shù)與其真值之間有一定的差異即所謂的觀測誤差，例如描述彈簧受迫振動的二階線性常系數(shù)微分方程中的質(zhì)量、阻尼系數(shù)和彈性系數(shù)等。更一般地：對物體的長寬高、電壓、溫度、速度的量測等。3.截斷誤差：許多數(shù)學運算是通過極限過程定義的，如微分、積分以及無窮級數(shù)求和等，由于計算機只能完成有限的算術(shù)運算和邏輯運算，所以在利用計算機進行計算是需要把無限的計算過程用有限的計算過程代替，由此產(chǎn)生

5、的誤差成為截斷誤差；4.舍入誤差：實際計算時只能按有限位進行，特別是里用計算機計算，由于計算機的有限位的限制，對參與運算的數(shù)據(jù)以及運算結(jié)果往往要進行舍入，例如利用公式計算圓的面積時， 需用有限的小數(shù)代替，由此產(chǎn)生的誤差成為舍入誤差。數(shù)值分析與計算主要研究截斷誤差和舍入誤差，研究誤差產(chǎn)生的原因、分析算法的誤差以及控制計算過程中誤差的擴散，由此把握計算結(jié)果的精度。a) 盡量減少計算步驟，簡化計算；減少運算次數(shù)，既可以提高解題速度，又有可能使計算中的舍入誤差積累減少。例：計算多項式的值。若采用逐項計算然后相加的算法所需的乘法次數(shù)為加法次數(shù)為n,若采用秦九韶算法則只需n次乘法和n次加法即可算出。b)

6、盡量避免兩個相近的數(shù)減；在數(shù)值計算中，兩個相近的數(shù)作減法時有效數(shù)字會損失。例：原始數(shù)據(jù)有9位有效數(shù)字，運算結(jié)果只有1位有效數(shù)字。c)盡量避免較大的數(shù)和較小的數(shù)相加；大數(shù)“吃掉”小數(shù)是指計算過程中，較小的數(shù)加不到較大的數(shù)中。這種現(xiàn)象有時會產(chǎn)生嚴重的后果。例：較小的數(shù)被較大的數(shù)“吃”掉了。d) 盡量避免絕對值較小的數(shù)做除數(shù)；在用計算機作運算時，絕對值很小的數(shù)做除數(shù)會溢出停機。而且當絕對值很小的除數(shù)稍微有一點誤差時，對計算結(jié)果影響很大。例：設(shè) ，計算結(jié)果起了很大的變化，因此，在計算中必須避免絕對值很小的數(shù)作除數(shù)。e) 選擇數(shù)值穩(wěn)定的算法定義：如果一個算法的計算結(jié)果對原始數(shù)據(jù)的誤差（擾動）以及計算過程

7、中的舍入誤差不敏感，則稱該算法是數(shù)值穩(wěn)定的，否則稱為數(shù)值不穩(wěn)定的。基本數(shù)值計算方法及其（結(jié)合Matlab相關(guān)函數(shù)）應(yīng)用示例線性代數(shù)方程組的解法 線性代數(shù)方程組的解法在數(shù)值分析中占有極其重要的位置，一方面，在工程技術(shù)領(lǐng)域，經(jīng)常以線性代數(shù)方程組作為其基本模型，例如，電學網(wǎng)絡(luò)問題，熱傳導問題、質(zhì)譜儀數(shù)據(jù)分析等；另一方面，在許多有效的數(shù)值方法中，求解線性代數(shù)方程組是其中關(guān)鍵的一步，比如樣條差值法、矩陣特征值問題、微分方程數(shù)值解法等，都離不開線性方程組的求解。 求解線性代數(shù)方程組的數(shù)值方法可分為兩大類：直接法和迭代法。直接法是指若沒有舍入誤差的影響，經(jīng)過有限次算數(shù)運算可求得方程組的精確解的方法。由于實際

8、計算時舍入誤差不可避免，故直接法得到的解并不精確。迭代法是利用迭代公式產(chǎn)生的向量序列去逼近精確解，當?shù)諗繒r，可得到滿足精確要求的近似解。一般的，對于低階的稠密線性方程組以及大型帶形方程組的解，采用直接法比較有效，而對于大型稀疏（非常形）方程組則用迭代法求解比較有利。.1高斯消元法高斯消去法的求解思路：把一般的線性方程組化成（上或下）梯形的形式。例：使用Matlab軟件，采用高斯消元法解方程Ax=b程序見附錄1。.2高斯列主元消元法 高斯消元法是按方程組和未知元的自然順序（即方程組對應(yīng)增廣矩陣的行和列的自然順序）逐步消元的。高斯消元法中若主元的絕對值很小，將使相應(yīng)的乘數(shù)的絕對值很大，從而導致

9、舍入誤差的嚴重擴散（放大），嚴重影響所得解的精確度。改進措施是，調(diào)整方程組或未知元的次序（即對增廣矩陣的行或列進行交換），使得系數(shù)矩陣中絕對值較大的元素作主元，這就是所謂的“選主元”技巧。列主元高斯消元法只是調(diào)整方程組（增廣矩陣的行）的次序，不改變未知元（增廣矩陣）的次序，與順序高斯消元法相比，只是在每一步消元前增加一步按列選主元的工作。例：使用Matlab軟件，采用列主元高斯消元法解方程Ax=b程序見附錄1。.3直接三角分解法 高斯消去法的消元過程，從代數(shù)運算的角度看就是用一個下三角矩陣左乘方程組的系數(shù)矩陣，且乘積的結(jié)果為上三角矩陣，即例：使用Matlab軟件，采用分解法解方程Ax=b程序見

10、附錄1。.4迭代法迭代法是解線性方程組的另一類方法，特別適用與解大型稀疏線性方程組（階數(shù)高、零元素多的方程組）。作為一種求解數(shù)值問題的通用方法，迭代法的基本思路是統(tǒng)一的，即利用設(shè)計好的迭代公式所產(chǎn)生的迭代序列去逐步逼近精確解。迭代法求解線性方程組的基本思想是a.不追求“一下子”得到方程組的解，而是在逐步逼近方程組的精確界的迭代過程中獲得滿足精度要求的近似值解，這一點與直接法不同；b.通過對問題的轉(zhuǎn)化，避免（困難的）矩陣求逆運算。由迭代格式確定如下的迭代算法：對于給定的線性方程組，可以寫成不同的（無窮多）迭代格式，有意義的（可用的）迭代格式應(yīng)具有收斂性生成的序列收斂于方程組的解；好的迭代法應(yīng)具有

11、較高的收斂速度。例：使用Matlab軟件，采用Jacobi迭代法驗證生成線性方程組是否收斂，若收斂求出譜半徑，并改變譜半徑的大小觀察其與收斂速度的關(guān)系。保證Jacobi迭代法，Gauss-Seidel迭代法都收斂的情況下，選擇最佳松弛因子計算線性方程組的解。經(jīng)過實驗發(fā)現(xiàn)運行時譜半徑反別為0.5984、0.2256、0.1215時，收斂速度分別為37、26、20由上面三次運算結(jié)果可以看出譜半徑越小收斂速度越快。SOR引入松弛因子后收斂速度明顯加快。程序見附錄2。3.4.2非線性代數(shù)方程組的解法一般非線性代數(shù)方程可有如下形式：設(shè) 在區(qū)間上連續(xù)。確定方程在該區(qū)間內(nèi)是否有實根，如果有實根，求某個或全部

12、的實根（近似值）或方程的近似解。.1 二分法設(shè)在區(qū)間上連續(xù)、單調(diào)并且 。如方程首先要確定適當?shù)陌膮^(qū)間，這可以依據(jù)閉區(qū)間上連續(xù)函數(shù)的介值定理來確定。注：二分法適用于一個方程的場合，收斂速度是線性的。迭代次數(shù)的估計：.2 黃金分割法在區(qū)間內(nèi)取對稱的兩點： 使得求根（方法）程序如下： 按這種方法選取點和 ，每次去掉的區(qū)間長度至少是原區(qū)間長度的0.618倍，上述適用于一個方程的場合，收斂速度是線性的。迭代次數(shù)的估計：例：使用Matlab軟件，采用黃金分割法解初始區(qū)間-3, 3， 收斂精度為的解。程序見附錄3。.3 牛頓迭代法牛頓迭代法構(gòu)造說明：設(shè)方程右端的函數(shù) 在所確定的包含唯一根的區(qū)間內(nèi)具有一

13、階連續(xù)導數(shù) ，那么由一階微分的性質(zhì)，有 由上式可得從而得到牛頓迭代格式 以及迭代法: 。Newton迭代法的幾何含義曲線在點 處的切線方程為 。 該切線與 軸的交點橫坐標就是Newton迭代公式確定的：例：使用Matlab軟件，利用函數(shù)Newton 迭代法求非線性方程組其中。的所有實根，誤差不超過，并在途中畫出所求實根對應(yīng)的各點。所得圖形如下，程序見附錄4。差值方法建立函數(shù)表達式、計算函數(shù)值等問題，都是數(shù)學中最基礎(chǔ)同時也很有實際需要的問題。差值問題就是這種問題中的一類近似處理方法。這些問題的產(chǎn)生，可能是已知函數(shù)太復雜，或者是已知函數(shù)只是一張表格，而這兩者都不便使用。這時，我們就希望通過已知函數(shù)

14、的一些離散的值，按照一定的要求，構(gòu)造出此函數(shù)的另外一個近似解表達式，一種常用的辦法是從某個性質(zhì)優(yōu)良、便于結(jié)算的函數(shù)類中選出一個函數(shù)，如代數(shù)多項式，尋找的方法就是插值法。定理：給定個曲線上點 ，如果 ，互不相同，那么，在所有次數(shù)不超過次的多項式函數(shù)中，存在唯一的多項式函數(shù)，滿足條件。.1拉格朗日插值法拉格朗日插值法是多項式插值中最基礎(chǔ)的一種。構(gòu)造插值多項式的基函數(shù) ， 拉格朗日插值多項式：拉格朗日插值法的不足：在實際問題中，觀測的數(shù)據(jù)可能會不斷增加，如果用拉格朗日插值公式構(gòu)造插值多項式，那么，每當增加數(shù)據(jù)就要重新計算多項式的系數(shù)，由此增加許多不必要的計算工作量。例：使用Matlab軟件利用Lag

15、range差值方法求溫度隨水深（近似）變化函數(shù)表達式。在夏季，較大湖泊的水體按深度被躍變層分為上部的變溫層和下部的均溫層。水體的分層化對環(huán)境工程中污染問題的研究具有重要的意義，例如，有機物的分解將導致被躍變層隔離的底部水體中氧急劇減少。按溫度隨水深的變化曲線，躍變層位于水深處：。 現(xiàn)有美國普拉特湖（Platte Lake） 的一組數(shù)據(jù)：深度(): 0溫度(): 所得圖形如下圖所示，程序見附錄5。所得方程為：.2三次樣條插值法三次樣條插值（簡稱Spline插值）是通過一系列形值點的一條光滑曲線，數(shù)學上通過求解三彎矩方程組得出曲線函數(shù)組的過程。早期工程師制圖時，把富有彈性的細長木條（所謂樣條）用壓

16、鐵固定在樣點上，在其他地方讓它自由彎曲，然后沿木條畫下曲線。成為樣條曲線。從數(shù)學上抽象就得到三次樣條函數(shù)這一概念。實際計算是還需要引入邊界條件才能完成計算。邊界通常有自然邊界（邊界點的導數(shù)為0），夾持邊界（邊界點導數(shù)給定），非扭結(jié)邊界（使兩端點的三階導與這兩端點的鄰近點的三階導相等）。a、確定三次樣條函數(shù)的條件根據(jù)三次樣條插值的要求，樣條插值函數(shù)在每個小區(qū)間 ， ，上都是三次多項式，每個三次多項式有4個系數(shù)，總共需要確定個系數(shù)，因此，需要個獨立條件才能保證唯一地確定滿足要求的三次樣條插值函數(shù)。已知：條件b給出了條件：， 條件c要求 在區(qū)間上具有二階連續(xù)導數(shù)，所以，插值函數(shù)在中間插值節(jié)點 ， ，

17、處，必須滿足條件：， ； ， ； ， ； 這樣已經(jīng)有個條件：還需要2個條件，才能保證對插值函數(shù)的唯一性要求；為此，通常在插值區(qū)間的端點處附加2個條件：第一類邊界條件： 固定端點的斜率：固定邊界條件： ，. ， 第二類邊界條件： 給定端點的的二階導數(shù)：自由邊界條件： 。 第三類邊界條件：周期性條件： ， 。 b、 三次樣條插值函數(shù)的構(gòu)造1. 設(shè)待求的三次樣條插值函數(shù) 在各插值節(jié)點處的一階、二階導數(shù)，由于插值函數(shù)在各小區(qū)間，上都是三次多項式，所以二階導數(shù) 是一次多項式，記 ， 利用Lagrange插值公式， 在各個區(qū)間上可表示為：， 此式的兩端積分兩次，得到樣條函數(shù)每個區(qū)間上的表達式： ， 在每個

18、小區(qū)間的左端點 處 ， 和右端點 處 ， ，利用插值條件：得到一階導數(shù)的表達式： ， (2由此確定了個系數(shù)： 。再利用一階導數(shù)在中間個點處的連續(xù)性：，；得到個等式（方程）： 經(jīng)進一步整理：， 最終歸結(jié)為求解 階的三對角線性方程組 的問題，其中： 對于自然邊界條件，有例：使用Matlab軟件利用三次樣條差值方法求溫度隨水深（近似）變化函數(shù)表達式。（題目如上題）所得圖形如下圖所示，程序見附錄6。所得方程為：。分析：根據(jù)常識我們可以知道，隨著向湖面或湖底延伸，水溫應(yīng)該呈現(xiàn)自然的延伸狀態(tài)，lagrange插值法沒有展現(xiàn)這種趨勢，而且在各個點之間變化比較大，而三次樣條插值法加上自己設(shè)置的邊界條件，所給出

19、的曲線較好地展示了這些特性，因此，相比較而言，三次樣條插值法更適用于此場合，更具有效性。.3曲線擬合的最小二乘法給出一組離散點，確定一個函數(shù)逼近原函數(shù)，插值是這樣的一種手段。在實際中，數(shù)據(jù)不可避免的會有誤差，插值函數(shù)會將這些誤差也包括在內(nèi)。因此，我們需要一種新的逼近原函數(shù)的手段：不要求過所有的點（可以消除誤差影響）；盡可能表現(xiàn)數(shù)據(jù)的趨勢，靠近這些點。設(shè)近似函數(shù)為： 函數(shù)值 與觀測值 之差稱為殘差 可以用殘差來衡量近似函數(shù)的好壞。擬合問題的幾何背景是尋求一條近似通過給定離散點的曲線，故稱曲線擬合問題。最小二乘法一般形式 為線性無關(guān)的基函數(shù)，求駐點，令 即令 則法方程組可寫成以下形式 函數(shù)空間的基

20、然后列出法方程 函數(shù)空間的基然后列出法方程 例：使用Matlab軟件利用最小二乘法擬合法按要求做下題。煉鋼是通過氧化降低含碳量的過程。鋼液含碳量的多少直接影響冶煉時間的長短，因此，必需了解它們之間的關(guān)系。下表給出的是一組熔畢碳 （爐料熔化完畢時）與精煉時間（從爐料熔化完畢到出爐）的數(shù)據(jù)。含碳量:104121134147150177180191190204時間: 分鐘 100125135155170185200205210235選擇不同次數(shù)的多項式，試用該組數(shù)據(jù)，畫出數(shù)據(jù)散點圖和擬合曲線圖形，比較擬合的效果，選擇你認為適當?shù)拇螖?shù)，并說明理由.所得圖形如下圖所示，程序見附錄6當n=1時 當n=2時

21、 當n=3時 當n=4時 結(jié)果分析：因為所給數(shù)據(jù)近似線性分布當n=1或2時數(shù)據(jù)點均勻分布在擬合線的兩側(cè)，即n=1或2就可滿足結(jié)果需求。.4最佳平方逼近多項式 所謂函數(shù)最佳逼近就是從指定的一類簡單的函數(shù)中尋找一個和給定的函數(shù)“最貼近”的函數(shù)。從幾何（空間）的角度看，函數(shù)最佳逼近就是從指定的函數(shù)集合（函數(shù)空間）中尋找一個函數(shù)（點），該函數(shù)和給定的函數(shù)（定點）最貼近（該點到定點的距離最?。?。為此，要給出函數(shù)之間距離的概念和定義，這時要把函數(shù)看做點（函數(shù)空間中的點）。區(qū)間 上連續(xù)函數(shù)組全體 構(gòu)成一個函數(shù)空間，可以看做該空間中的一個點；所有次數(shù)不超過次的多項式集合 也構(gòu)成函數(shù)空間。顯然有 例：使用Mat

22、lab軟件利用最佳平方逼近多項式求函數(shù)在區(qū)間上距離不超過 的最佳平方逼近多項式。在同一窗口畫出該函數(shù)和最佳逼近多項式曲線，寫出該多項式的表達式。所得圖形如下圖所示，程序見附錄7N =8整理得：數(shù)值積分和數(shù)值微分公式.1數(shù)值積分對于給定的函數(shù)，如果1. 的函數(shù)關(guān)系式比較復雜；2. 未知，而僅僅知道該函數(shù)在個點，處的函數(shù)值 .則希望能用相對簡單的計算方法，求得在上的定積分的近似值。a、 兩點公式（梯形公式）b、 三點公式（Simpson公式）c、復化的梯形求積公式d、復化的Simpson求積公式復化的求積公式是受“利用低次多項式分段插值”能夠提高逼近精度且計算簡單的事實啟發(fā)，將上述的“在整個積分區(qū)

23、間上用一個插值多項式代替被積函數(shù)”的做法變?yōu)椤坝梅侄尾逯刀囗検酱姹环e函數(shù)”。f、Monte Carlo 求積法 函數(shù)在區(qū)間上的平均值為當那么有由此，利用在上服從均勻分布的隨機數(shù)生成函數(shù)unifrnd(a,b,m,n)生成，代入。例：使用Matlab軟件利用復化的梯形求積公式、復化的Simpson求積公式、Monte Carlo 求積法求解下題。下面是一處地質(zhì)巖層斷面上部邊緣的深度測量數(shù)據(jù)。水平距離():0深度(): 1.1 試利用復化的梯形求積法求該組數(shù)據(jù)所在曲線與基準線（ 軸）在范圍 內(nèi)所圍成圖形面積.畫出數(shù)據(jù)散點圖和圖形的示意圖.所得圖形如下圖所示，程序見附錄81.2 試利用復化的Sim

24、pson求積法求該組數(shù)據(jù)所在曲線與基準線（ 軸）在范圍 內(nèi)所圍成圖形面積.畫出數(shù)據(jù)散點圖和圖形的示意圖.所得圖形如下圖所示，程序見附錄8利用Monte Carlo法求面積（注：各做 次試驗，取平均值作為面積值）。程序見附錄8結(jié)果分析：法一、選擇相鄰兩點為一組按照梯形的面積公式即可求出相鄰兩點的面積，然后相加即可。法二、選擇相鄰三點為一組，中間的點位插入點，然后按照Simpson求積法求出相鄰三點的面積，然后相加即可。.2數(shù)值微分對于給定的函數(shù)，如果1. 的函數(shù)關(guān)系式比較復雜；2. 未知，而僅僅知道該函數(shù)在個點，處的函數(shù)值 .則希望能用相對簡單的計算方法，求得導數(shù)的近似值?；谏鲜隹紤]，選擇的方

25、法之一是利用函數(shù)的插值多項式的導數(shù)作為函數(shù) 導數(shù)的近似值，例如Lagrange插值多項式，由因而有這里需要說明一點的是，盡管 和的函數(shù)值可能相差不多，但是它們的導數(shù)有可能相差很大。常微分方程的數(shù)值解法 微分方程（組）是科學研究和工程應(yīng)用中最常用的數(shù)學模型之一。如揭示質(zhì)點運動規(guī)律的Newton第二定律：和刻畫回路電流或電壓變化規(guī)律的基爾霍夫回路定律等，但是，只有一些簡單的和特殊的常微分方程及常微分方程組，可以求得用公式給出的所謂解析解或公式解，如一階線性微分方程的初值問題：的解為：但是，絕大多數(shù)在實際中遇到的常微分方程和常微分方程組得不到“解析解”，因此，基于如下的事實：I. 絕大多數(shù)的常微分方

26、程和常微分方程組得不到解析解；II. 實際應(yīng)用中往往只需要知道常微分方程（組）的解在（人們所關(guān)心的,感興趣的）某些點處的函數(shù)值（可以是滿足一定精度要求的近似值）；如果只需要常微分方程（組）的解在某些點處的函數(shù)值，則沒有必要非得通過求得解的表達式，然后再計算出函數(shù)值不可，事實上，我們可以采用下面將介紹的常微分方程（組）的初值問題的數(shù)值解法，就可以達到這一目的。一般的一階常微分方程初值問題是指如下的一階常微分方程的定解問題： 微分方程（組）的初值問題通常是對一動態(tài)過程演化規(guī)律的描述，求解常微分方程（組）的初值問題就是要了解和掌握動態(tài)過程演化規(guī)律。解的存在且唯一性條件（定理）：如果函數(shù)在區(qū)域上連續(xù)，

27、且存在（Lipschitz）常數(shù),對自變量 滿足不等式（Lipschitz條件）：那么初值問題（Cauchy 問題）存在唯一的連續(xù)可微解 。.1 Euler 方法假設(shè)要求在點（時刻） ， ， 處初值問題的解 的近似值。首先對式的兩端積分，得對其右邊，如果被積函數(shù)用積分下限處的函數(shù)值代替被積函數(shù)作積分（從幾何上的角度看，是用矩形面積代替曲邊梯形面積），則有進而得到下式給出的遞推算法Euler 方法改進的Euler 方法對于的右邊，如果被積函數(shù)用積分限和處的函數(shù)值的算術(shù)平均值代替（幾何上，是用梯形面積代替曲邊梯形面積），則有進而得到下式給出的遞推算法：通常算法改進的上式方法比Euler 方法的精度

28、高，但是，按算法Euler方法求 時要解（非線性）方程，這是上式方法不如Euler上式算法比Euler方法精度高的優(yōu)點；II 盡可能地利用Euler方法計算簡單的長處；人們采取了如下的稱之為改進的Euler 方法的折衷方案：預測 修正 .2 RungeKutta法RungeKutta法是按選取區(qū)間 上函數(shù) 變化率的個數(shù)的多少和截斷誤差的階數(shù) 來區(qū)分的一系列方法，如I. 二階的RungeKutta法（改進的Euler方法） II. 三階的RungeKutta法 III. 四階的RungeKutta法III.1 古典形式 III.2 Gill 公式(具有減小舍入例：使用Matlab軟件利用Eule

29、r 法，4階R-K公式解下題。給定如下常微分方程Cauchy 問題：1.分別取參數(shù)，初值和步長：，1.1 取 ，在區(qū)間 上分別利用Euler 法，4階R-K公式和matlab的ode45求問題（1）的數(shù)值解；的曲線；2.1 分別取 和 ，重復1.1的實驗和1.2的做法，給出你根據(jù)觀察，對微分方程（1）的解 變化趨勢的判斷。注：上述的實驗2.1和2.2 初值可按如下方法隨機確定：實驗結(jié)果及分析如下，程序見附錄9運行結(jié)果：只分別去了他們最后一個結(jié)果 結(jié)果分析：由（）分別取 和，重復1.1的實驗和1.2的做法得到的圖可以看出當區(qū)間變大時，得到的解趨于常數(shù)1.有圖像也可以看出由4階R-K公式得到的曲線

30、比Euler 法得到的圖形曲線等貼近Matlab得到的圖象。.3 Matlab法 ode是Matlab專門用于解微分方程的功能函數(shù)。該求解器有變步長（variable-step）和定步長（fixed-step）兩種類型。不同類型有著不同的求解器，其中ode45求解器屬于變步長的一種，采用Runge-Kutta算法；和他采用相同算法的變步長求解器還有ode23。ode45表示采用四階，五階Runge-Kutta單步算法,截斷誤差為(x)3。解決的Nonstiff(非剛性)常微分方程。解決數(shù)值解問題的首選方法是ode45，若長時間沒結(jié)果，應(yīng)該就是剛性的，可換用ode23試試。Matlab中求微分方

31、程數(shù)值解的函數(shù)有七個：ode45，ode23，ode113，ode15s，ode23s，ode23t，ode23tb 。ode45 此指令以 Runge-Kutta (4,5) 演算法為基礎(chǔ)，運算上是屬單步驟計算運y(tn)，過僅需時間點前之解y(tn-1)。 ode45函數(shù)指令可作為第一階段評估的數(shù)據(jù)，然后再思考是否采取下一步驟。ode23與ode45函數(shù)同以 Runge-Kutta 方程式為基礎(chǔ)，并以Bogacki與Shampine配對。在溫和剛性的情況，容許度較粗略時比ode45 更有效率。此指令也是單步驟解題。ode15s以數(shù)值微分方程(NDFs).為基礎(chǔ)之可變順序指令。可選擇后向微分

32、方程BDFs, (或稱為Gear 法)。ode15s是一個多步驟指令。若執(zhí)行前認為問題特性屬于剛性，亦或使用過ode45指令但無法達成目的時，可以試用f ode15s指令。ode23s以Rosenbrock修正方程式二階為基礎(chǔ)。由于它是單步驟解題指令，故在粗略容許范圍下可能比ode15s 更有效率。它可能解一些ode15s效率低或無法解的剛性問題。ode23t采用自由間機的梯形法。在問題屬于中等剛性或不需含括數(shù)值阻尼的特性時可以使用此指令解題。ode23tb執(zhí)行TR-BDF2初期利用梯形法含Runge-Kutta 方程式而第二階段采用后向微分法二階。如同ode23s一樣，在粗略容許度下此法比o

33、de15s有效率。語法T,Y = ode45(odefun,tspan,y0) T,Y = ode45(odefun,tspan,y0,options)T,Y,TE,YE,IE = ode45(odefun,tspan,y0,options) sol = ode45(odefun,t0tf,y0.)T,Y = ode45(odefun,tspan,y0)odefun 是函數(shù)句柄,可以是函數(shù)文件名,匿名函數(shù)句柄或內(nèi)聯(lián)函數(shù)名tspan 是區(qū)間 t0 tf 或者一系列散點t0,t1,.,tfy0 是初始值向量 Simulink與信號處理T 返回列向量的時間點 Y 返回對應(yīng)T的求解列向量T,Y = o

34、de45(odefun,tspan,y0,options)options 是求解參數(shù)設(shè)置,可以用odeset在計算前設(shè)定誤差,輸出參數(shù),事件等T,Y,TE,YE,IE =ode45(odefun,tspan,y0,options)在設(shè)置了事件參數(shù)后的對應(yīng)輸出TE 事件發(fā)生時間 YE 事件解決時間IE The index i of the event functionthat vanishes. sol =ode45(odefun,t0 tf,y0.) sol 結(jié)構(gòu)體輸出結(jié)果例：使用Matlab軟件利用matlab法解下列題目二階常微分方程是一類典型的、有十分廣泛應(yīng)用背景的數(shù)學模型。令例1、考慮如下二階常微分方程Cauchy 問題：取求解（仿真）的終止時刻 中的某個值。令是如下特征方程的根: 。 參數(shù)分別取如下五組：， （4-1）， （4-2）， （4-3）， （4-4） ， （4-5） 試分別用Matlab的ode45求上述問題的數(shù)值解。并按每組參數(shù)（4-1）-（4-5），分別在在同一個圖形窗口畫出所求數(shù)值解在區(qū)間 的圖像。1.2 針對（4-1）-（4-5）中的每組參數(shù)，觀察所求數(shù)值解曲線隨時間的演化趨勢，