matlab計(jì)算方法實(shí)驗(yàn)指導(dǎo)誤差分析

1、實(shí)驗(yàn)一誤差分析實(shí)驗(yàn)1（病態(tài)問題）實(shí)驗(yàn)?zāi)康模核惴ㄓ小皟?yōu)”與“劣”之分，問題也有“好”與“壞”之別。對數(shù)值方法的研究而言，所謂壞問題就是問題本身對擾動敏感者，反之屬于好問題。通過本實(shí)驗(yàn)可獲得一個(gè)初步體會。數(shù)值分析的大部分研究課題中，如線性代數(shù)方程組、矩陣特征值問題、非線性方程及方程組等都存在病態(tài)的問題。病態(tài)問題要通過研究和構(gòu)造特殊的算法來解決，當(dāng)然一般要付出一些代價(jià)（如耗用更多的機(jī)器時(shí)間、占用更多的存儲空間等）。問題提出：考慮一個(gè)高次的代數(shù)多項(xiàng)式p（x）=（x，l）（x，2）（x，20）=（x，k）（1.1）k=1顯然該多項(xiàng)式的全部根為1,2,20共計(jì)20個(gè)，且每個(gè)根都是單重的?，F(xiàn)考慮該多項(xiàng)式的一

2、個(gè)擾動p（x）+x19=0（1.2）其中是一個(gè)非常小的數(shù)。這相當(dāng)于是對（1.1）中x19的系數(shù)作一個(gè)小的擾動。我們希望比較（1.1）和（1.2）根的差別，從而分析方程（1.1）的解對擾動的敏感性。實(shí)驗(yàn)內(nèi)容:為了實(shí)現(xiàn)方便，我們先介紹兩個(gè)MATLAB函數(shù):“roots”和“poly”。u=roots（a）其中若變量a存儲n+1維的向量，則該函數(shù)的輸出u為一個(gè)n維的向量。設(shè)a的元素依次為a,a，,a，則輸出u的各分量是多項(xiàng)式方程12n+1axn+axn-1+ax+a=012nn+1的全部根；而函數(shù)b=poly（v的輸出b是一個(gè)n+1維向量，它是以n維向量v的各分量為根的多項(xiàng)式的系數(shù)?？梢姟皉oots

3、”和“poly”是兩個(gè)互逆的運(yùn)算函數(shù)。ess=0.000000001ve=zeros（1,21）;ve（2）=ess;roots（poly（1:20）+ve）上述簡單的MATLAB程序便得到（1.2）的全部根，程序中的“ess”即是（1.2）中的。實(shí)驗(yàn)要求:選擇充分小的ess，反復(fù)進(jìn)行上述實(shí)驗(yàn)，記錄結(jié)果的變化并分析它們。如果擾動項(xiàng)的系數(shù)很小，我們自然感覺(1.1)和(1.2)的解應(yīng)當(dāng)相差很小。計(jì)算中你有什么出乎意料的發(fā)現(xiàn)？表明有些解關(guān)于如此的擾動敏感性如何？將方程(1.2)中的擾動項(xiàng)改成x18或其它形式，實(shí)驗(yàn)中又有怎樣的現(xiàn)象出現(xiàn)？(選作部分)請從理論上分析產(chǎn)生這一問題的根源。注意我們可以將方程

4、(1.2)寫成展開的形式，p(x,)二X20-Ox19,0(1.3)同時(shí)將方程的解X看成是系數(shù)的函數(shù)，考察方程的某個(gè)解關(guān)于的擾動是否敏感，與研究它關(guān)于的導(dǎo)數(shù)的大小有何關(guān)系？為什么？你發(fā)現(xiàn)了什么現(xiàn)象，哪些根關(guān)于的變化更敏感？思考題一：(上述實(shí)驗(yàn)的改進(jìn))在上述實(shí)驗(yàn)中我們會發(fā)現(xiàn)用roots函數(shù)求解多項(xiàng)式方程的精度不高，為此你可以考慮用符號函數(shù)solve來提高解的精確度，這需要用到將多項(xiàng)式轉(zhuǎn)換為符號多項(xiàng)式的函數(shù)poly2sym，函數(shù)的具體使用方法可參考MATLAB的幫助。思考題二：(二進(jìn)制產(chǎn)生的誤差)用MATLAB計(jì)算0.1-100。結(jié)果居然有誤差！因?yàn)閺氖M(jìn)制數(shù)角度分析,i=1這一計(jì)算應(yīng)該是準(zhǔn)確的。

5、實(shí)驗(yàn)反映了計(jì)算機(jī)內(nèi)部的二進(jìn)制本質(zhì)。思考題三：(一個(gè)簡單公式中產(chǎn)生巨大舍入誤差的例子)可以用下列式子計(jì)算自然對數(shù)的底數(shù)1e=e1=lim(1+)nngn這個(gè)極限表明隨著n的增加，計(jì)算e值的精度是不確定的。現(xiàn)編程計(jì)算1f(n)=(1+-)n與exp(1)值的差。n大到什么程度的時(shí)候誤差最大？你能解釋其中n的原因嗎？相關(guān)MATLAB函數(shù)提示：poya求給定的根向量a生成其對應(yīng)的多項(xiàng)式系數(shù)(降序)向量roots(p)求解以向量p為系數(shù)的多項(xiàng)式(降序)的所有根poly2sym(p)將多項(xiàng)式向量p表示成為符號多項(xiàng)式(降序)sym(arg)將數(shù)字、字符串或表達(dá)式arg轉(zhuǎn)換為符號對象symsarg1arg2a

6、rgk將字符arg1,arg2,argk定義為基本符號對象solve(eq1)求符號多項(xiàng)式方程eq1的符號解實(shí)驗(yàn)二非線性方程求根實(shí)驗(yàn)2(迭代法、初始值與收斂性)實(shí)驗(yàn)?zāi)康模撼醪秸J(rèn)識非線性問題的迭代法與線性問題迭代法的差別，探討迭代法及初始值與迭代收斂性的關(guān)系。問題提出：迭代法是求解非線性方程的基本思想方法，與線性方程的情況一樣，其構(gòu)造方法可以有多種多樣，但關(guān)鍵是怎樣才能使迭代收斂且有較快的收斂速度。實(shí)驗(yàn)內(nèi)容：考慮一個(gè)簡單的代數(shù)方程X2X10針對上述方程，可以構(gòu)造多種迭代法，如XX21(7.1)n+1nX1+(7.2)n+1XnXn+1=X+1n(7.3)在實(shí)軸上取初始值X。，請分別用迭代(7.1

7、)-(7.3)作實(shí)驗(yàn)，記錄各算法的迭代過程。實(shí)驗(yàn)要求：取定某個(gè)初始值，分別計(jì)算(7.1)-(7.3)迭代結(jié)果，它們的收斂性如何？重復(fù)選取不同的初始值，反復(fù)實(shí)驗(yàn)。請自選設(shè)計(jì)一種比較形象的記錄方式(如利用MATLAB的圖形功能)，分析三種迭代法的收斂性與初值選取的關(guān)系。(2)對三個(gè)迭代法中的某個(gè)，取不同的初始值進(jìn)行迭代，結(jié)果如何？試分析迭代法對不同的初值是否有差異？(3)線性方程組迭代法的收斂性是不依賴初始值選取的。比較線性與非線性問題迭代的差異，有何結(jié)論和問題。思考題一：用Newton法求方程X3一X一10在區(qū)間-3,3上誤差不大于10-5的根，分別取初值x1.5,0,1進(jìn)行計(jì)算，比較它0們的迭

8、代次數(shù)。相關(guān)MATLAB函數(shù)提示：x=fzero(fun,xO)返回一元函數(shù)fun的一個(gè)零點(diǎn)，其中fun為函數(shù)句柄，x0為標(biāo)量時(shí)，返回在x0附近的零點(diǎn)；x0為向量a,b時(shí)，返回函數(shù)在a,b中的零點(diǎn)x,f,h=fsolve(fun,x0)返回一元或多元函數(shù)x0附近fun的一個(gè)零點(diǎn)，其中fun為函數(shù)句柄，x0為迭代初值；f返回fun在x的函數(shù)值，應(yīng)該接近0；h返回值如果大于0,說明計(jì)算結(jié)果可靠，否則不可靠實(shí)驗(yàn)三解線性方程組的迭代法實(shí)驗(yàn)3.1（病態(tài)的線性方程組的求解）問題提出：理論的分析表明，求解病態(tài)的線性方程組是困難的。實(shí)際情況是否如此，會出現(xiàn)怎樣的現(xiàn)象呢？實(shí)驗(yàn)內(nèi)容：考慮方程組Hx=b的求解，其中

9、系數(shù)矩陣H為Hilbert矩陣，H(h),i,jnxnh,i,j1,2,ni,jI,j1這是一個(gè)著名的病態(tài)問題。通過首先給定解（例如取為各個(gè)分量均為1）再計(jì)算出右端b的辦法給出確定的問題。實(shí)驗(yàn)要求：（1）選擇問題的維數(shù)為6,分別用Gauss消去法、J迭代法、GS迭代法和SOR迭代法求解方程組，其各自的結(jié)果如何？將計(jì)算結(jié)果與問題的解比較，結(jié)論如何？（2）逐步增大問題的維數(shù)，仍然用上述的方法來解它們，計(jì)算的結(jié)果如何？計(jì)算的結(jié)果說明了什么？思考題一：討論病態(tài)問題求解的算法Jacobi迭代法與Gauss-seidel迭代法的比較：用Jacobi迭代法與Gauss-seidel迭代法解下列方程組：10(

10、1)-1111X15110 x2-7-3丿XV3丿v-171(2)0.5、0.50.510.505、0.5x1x20.5I-25丿、12取X（0）（0,0,0）T,10-6，你能得出什么結(jié)論？思考題二：SOR（超松馳）迭代法松馳因子對收斂性及速度的影響試用SOR（超松馳）迭代法求解下列方程組：43340-101x12411-1x2304丿xV3丿V-24取X（0）（1,1,1）t,10-6，選擇松馳因子0.8,091,1.1,1.2等,試看對算法收斂性的影響，并找出所選用的松馳因子的最佳者。要求編制矩陣迭代求解的函數(shù)文件。實(shí)驗(yàn)3.2方程組性態(tài)討論（1）求Axb的解向量，其中：129，1，A=3

11、00020001000b二20004/1063/1062/106丿3/106丿求系數(shù)矩陣A的條件數(shù)；將a改為3/106,b改為4/106,求解向量壬;333令P=diag(1,10-3,106),求解PAx=Pb,并求系數(shù)矩陣PA的條件數(shù);對PA中的a和Pb中的b給以10-6的擾動，求解向量X333結(jié)合上述結(jié)果，討論兩個(gè)方程組的性態(tài)實(shí)驗(yàn)3.2大型稀疏方程組的數(shù)值解法設(shè)n階方陣：1-41-2-31-41-21-2-1-41-21-43-1-231-23丿B為A的各元素之和，顯然Ax=b的解為x=(1,1,丄)丁，用下面三種方法對于階數(shù)n=100,200,.,500，誤差限為=10-2,10-3,

12、丄0-6的各種組合求解，分析收斂速度選取列主元Gauss消元法；Jacobi迭代法Gauss-seidel迭代法算法的改進(jìn)：由于是稀疏矩陣，請考慮改進(jìn)矩陣的存儲方式，減少存儲空間的使用來提高計(jì)算效率。也就是用稀疏存儲方式。實(shí)驗(yàn)四解線性方程組的直接方法實(shí)驗(yàn)4(主元的選取與算法的穩(wěn)定性)問題提出：Gauss消去法是我們在線性代數(shù)中已經(jīng)熟悉的。但由于計(jì)算機(jī)的數(shù)值運(yùn)算是在一個(gè)有限的浮點(diǎn)數(shù)集合上進(jìn)行的，如何才能確保Gauss消去法作為數(shù)值算法的穩(wěn)定性呢？Gauss消去法從理論算法到數(shù)值算法，其關(guān)鍵是主元的選擇。主元的選擇從數(shù)學(xué)理論上看起來平凡，它卻是數(shù)值分析中十分典型的問題。實(shí)驗(yàn)內(nèi)容：考慮線性方程組Ax

13、二b,ARn,n,bRn又能手動選取主元的求解線性方程組的Gauss編制一個(gè)能自動選取主元消去過程。實(shí)驗(yàn)要求：61-786115,b=861158614(1)取矩陣A=，則方程有解x*二(1,1,，,1)T。取n=10計(jì)算矩陣的條件數(shù)。讓程序自動選取主元，結(jié)果如何？現(xiàn)選擇程序中手動選取主元的功能。每步消去過程總選取按模最小或按模盡可能小的元素作為主元，觀察并記錄計(jì)算結(jié)果。若每步消去過程總選取按模最大的元素作為主元，結(jié)果又如何？分析實(shí)驗(yàn)的結(jié)果。取矩陣階數(shù)n=20或者更大，重復(fù)上述實(shí)驗(yàn)過程，觀察記錄并分析不同的問題及消去過程中選擇不同的主元時(shí)計(jì)算結(jié)果的差異，說明主元素的選取在消去過程中的作用。將上述矩陣A中的主元改為0.00006再重新作一次數(shù)值實(shí)驗(yàn)看看。選取其他你感興趣的問題或者隨機(jī)生成矩陣，計(jì)算其條件數(shù)。重復(fù)上述實(shí)驗(yàn)，觀察記錄并分析實(shí)驗(yàn)結(jié)果。相關(guān)MATLAB函數(shù)提示：zeros(m,n)生成m行，n列的零矩陣ones(m,n)生成m行，n列的元素全為1的矩陣eye(n)生成n階單位矩陣rand(m,n)生成m行,n列(0,1)上均勻分布的隨機(jī)矩陣diag(x)返回由向量x的元素構(gòu)