數(shù)值計算方法實驗指導講解

《數(shù)值計算方法》

實驗指導（Matlab版）肇慶學院數(shù)學與統(tǒng)計學學院

計算方法課程組《數(shù)值計算方法》實驗1報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗1算法設計原則驗證(之相近數(shù)相減、大數(shù)吃小數(shù)和簡化計算步驟)實驗題目取z=1。16,計算-云和1/(t7軍1+云)，驗證兩個相近的數(shù)相減會造成有效數(shù)字的損失.按不同順序求一個較大的數(shù)(123)與1000個較小的數(shù)(3X10-15)的和，驗證大數(shù)吃小數(shù)的現(xiàn)象.分別用直接法和秦九韶算法計算多項式P(x)=axn+axn-1++ax+a在x=1.00037處的值.驗證簡化計算步驟能減少運算時間.n+1，一,對于第(3)題中的多項式P(x),直接逐項計算需要n+(n-1)+???+2+1=——次乘法和n次加法，使用秦九韶算法P(x)=(((ax+a)x+a)x—+a)x+a則只需要n次乘法和n次加法.實驗目的驗證數(shù)值算法需遵循的若干規(guī)則.基礎理論設計數(shù)值算法時，應避免兩個相近的數(shù)相減、防止大數(shù)吃小數(shù)、簡化計算步驟減少運算次數(shù)以減少運算時間并降低舍入誤差的積累.兩相近的數(shù)相減會損失有效數(shù)字的個數(shù)，用一個大數(shù)依次加小數(shù)，小數(shù)會被大數(shù)吃掉，乘法運算次數(shù)太多會增加運算時間.實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程直接計算并比較；法1：大數(shù)逐個加1000個小數(shù)，法2：先把1000個小數(shù)相加再與大數(shù)加；將由高次項到低次項的系數(shù)保存到數(shù)組A[n]中，其中n為多項式次數(shù).結(jié)果與分析(1)計算的k：z+1—z=，1/(yZ+1+tz).分析：(2)123逐次加1000個3x10-6的和是，先將1000個3X10-6相加，再用這個和與123相加得.分析：(3)計算次的多項式：直接計算的結(jié)果是，用時;用秦九韶算法計算的結(jié)果是，用時分析：8.附錄：程序清單兩個相近的數(shù)相減.%**個**********************************************************%*程序名：ex1_1.m*%*程序功能：驗證兩個相近的數(shù)相減會損失有效數(shù)字個數(shù)%*************************************************************z=1e16;%x=氣Z+1—tZ；%y=1/"E+.&)x,y(2)大數(shù)吃小數(shù)%*************************************************************%*程序名：ex1_2.m*%*程序功能：驗證大數(shù)吃小數(shù)的現(xiàn)象.%*************************************************************clc;%清屏clearall;%釋放所有內(nèi)存變量formatlong;%按雙精度顯示浮點數(shù)z=123;%大數(shù)t=3e-15;%小數(shù)x=z;%大數(shù)依次加小數(shù)%重復1000次給x中加上t

y=0；%先累加小數(shù)%重復1000y=0；%先累加小數(shù)%重復1000次給y中加上ty=z+y;%再加到大數(shù)(3)秦九韶算法%**個**********************************************************%*程序名：ex1_3.m*%*程序功能：驗證秦九韶算法可節(jié)省運行時間.%*************************************************************clc;%清屏clearall;%釋放所有內(nèi)存變量formatlong;%按雙精度顯示浮點數(shù)A=[8,4,-1,-3,6,5,3,2,1,3,2,-1,4,3,1,-2,4,6,8,9,50,-80,12,35,7,-6,42,5,6,23,74,65,55,80,78,77,98,56];A(10001)=0;%擴展到10001項，后面的都是分量0%A為多項式系數(shù)，從高次項到低次項x=1.00037;n=9000;%n為多項式次數(shù)%直接計算begintime=clock;P=0:鼬i=「：-l』T=1:奴k=l:it=t*京.endp=p-A(n-i-L)%;end.endtime=clock;time1=etime(endtime,begintime);disp('直接計算');disp(['p(',num2str(x),')=',num2str(p)]);%開始執(zhí)行的時間%求x的i次冪%累加多項式的i次項%結(jié)束執(zhí)行的時間%運行時間disp(['運行時間：',num2str(time1),'秒秒]);

%%開始執(zhí)行的時間%求x的i次冪%累加多項式的i次項%結(jié)束執(zhí)行的時間%運行時間%秦九韶算法計算%開始執(zhí)行的時間%累加秦九韶算法中的一項begintime=clock;fori=0;p-=p-A{i-1):endendtime=clock;disp(['運行時間：',num2str(time2),秒]);《數(shù)值計算方法》實驗1報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗1算法設計原則驗證（之數(shù)值穩(wěn)定性）實驗題目計算定積分I=j1心心一i&,n=0,1,10，分別用教材例1-7推導出的算法A和B，n0其中：算法A：算法A：<nI0=1一nIn-1r0.63216II一(1-I)算法B：'n-1nn110^0驗證算法不穩(wěn)定時誤差會擴大.實驗目的驗證數(shù)值算法需遵循的若干規(guī)則.基礎理論設計數(shù)值算法時，應采用數(shù)值穩(wěn)定性好的算法.數(shù)值穩(wěn)定的算法，誤差不會放大，甚至會縮??；而數(shù)值不穩(wěn)定的算法會放大誤差.實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程分別用數(shù)組IA[]和IB[]保存兩種算法計算的結(jié)果.8.附錄：程序清單%**個**********************************************************%*程序名：ex1_4.m%*程序功能：驗證數(shù)值穩(wěn)定性算法可控制誤差.%*************************************************************clc;%清屏clearall;%釋放所有內(nèi)存變量formatlong;%按雙精度顯示浮點數(shù)I=[0.63212055882856,0.36787944117144,0.26424111765712,0.20727664702865,0.17089341188538,0.14553294057308,0.12680235656154,0.11238350406938,.0.10093196744492,0.09161229300662,0.08387707010843];%保留14位小數(shù)的精確值，?是Matlab中的續(xù)行符%算法AIA(1)=0.6321;%Matlab下標從1開始，所以要用IA(n+1)表示原問題中的I(n)forn=1:10IA(n-L)=1-n*lA(n):end.%算法B0:forn=10:-1:1/n:enddisp('n算法A算法B精確值')；forn=1:11fprintf('%2d%14.6f%14.6f%14.6f\n',n-1,IA(n),IB(n),I(n));end%n顯示為2位整數(shù)，其它顯示為14位其中小數(shù)點后顯示6位的小數(shù)《數(shù)值計算方法》實驗1報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗1算法設計原則（除數(shù)絕對值不能太小）實驗題目將線性方程組增廣矩陣利用初等行變換可化為(a11a12a22r—%r2an1T01a12a22r—―r1。22Tf(a11a12a22r—%r2an1T01a12a22r—―r1。22Tfa11〔0由此可解得T叩an，%分別解增廣矩陣為實驗目的驗證數(shù)值算法需遵循的若干規(guī)則.基礎理論設計數(shù)值算法時，應避免除數(shù)絕對值遠小于被除數(shù)絕對值的除法，否則絕對誤差會被放大，使結(jié)果失真.實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程用二維數(shù)組A和B存放方程組的增廣矩陣，利用題目所給初等行變換求解方程組.結(jié)果與分析第1種順序的方程組的解為x=，y=；第2種順序的方程組的解為x=，y=.分析：8.附錄：程序清單%**個**********************************************************%*程序名：ex1_5.m*%*程序功能：驗證除數(shù)的絕對值太小可能會放大誤差.%*************************************************************

clc;A=[1e-16,1,1;2,1,2];B=[2,1,2;1e-16,1,1];%方程組A1);m=-A(l:2)/A(2;2);A(1:：)=A(1:A(l;3)=A(L3)/A(l=1);%增廣矩陣%m=-a_{21}/a_{11}是第2行加第1行的倍數(shù)%消去a_{21}%m=-a_{12}/a_{22}是第1行加第2行的倍數(shù)%消去a_{12},系數(shù)矩陣成對角線%未知數(shù)x1的值%未知數(shù)x2的值disp(['方程組A的解:x1=',num2str(A(1,3)),',x2=',num2str(A(2,3))]);disp('');%方程組B%m=-b_{21}/b_{11}是第2行加第1clc;A=[1e-16,1,1;2,1,2];B=[2,1,2;1e-16,1,1];%方程組A1);m=-A(l:2)/A(2;2);A(1:：)=A(1:A(l;3)=A(L3)/A(l=1);%增廣矩陣《數(shù)值計算方法》實驗2報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗2非線性方程的迭代解法(之簡單迭代法)實驗題目用簡單迭代法求方程尤3+4x2-10=0在區(qū)間［1,2］內(nèi)的一個實根，取絕對誤差限為10-4.實驗目的掌握非線性方程的簡單迭代法.基礎理論簡單迭代法：將方程f(x)=0改寫成等價形式x=9(x)，從初值x0開始，使用迭代公式x*=9(xk)可以得到一個數(shù)列，若該數(shù)列收斂，則其極限即為原方程的解.取數(shù)列中適當?shù)捻椏勺鳛榻平?實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程結(jié)果與分析附錄：程序清單《數(shù)值計算方法》實驗2報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗2非線性方程的迭代解法(之Newton迭代法)實驗題目用Newton迭代法求方程尤3+4尤2-10=0在區(qū)間［1,2］內(nèi)的一個實根，取絕對誤差限為10-4.實驗目的掌握求解非線性方程的Newton迭代法.基礎理論Newton迭代法：解方程/(x)=0的Newton迭代公式為x=x-k.Ef'(x)k實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程結(jié)果與分析附錄：程序清單《數(shù)值計算方法》實驗2報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗2非線性方程的迭代解法(之對分區(qū)間法)實驗題目用對分區(qū)間法求方程x3-x-1=0在區(qū)間［1,1.5］內(nèi)的一個實根，取絕對誤差限為10-4.實驗目的掌握求解非線性方程的對分區(qū)間法.基礎理論對分區(qū)間法：?。踫b］的中點p,若f(p)=0或b-a<們則p為方程f(x)=0的近似解；若f(a)f(p)<0，則說明根在區(qū)間?。踑,p］中；否則，根在區(qū)間?。踦,b］中.將新的有根區(qū)間記為［a”b1］,對該區(qū)間不斷重復上述步驟，即可得到方程的近似根.實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程用宏定義函數(shù)f(x)；為了循環(huán)方便，得到的新的有根區(qū)間始終用［a，b］表示；由于新的有根區(qū)間可能仍以a為左端點，這樣會反復使用函數(shù)值f(a)，為減少運算次數(shù)，將這個函數(shù)值保存在一個變量fa中；同樣在判斷新的有根區(qū)間時用到函數(shù)值f(p)，若新的有根區(qū)間以p為左端點，則下一次用到的f(a)實際上就是現(xiàn)在的f(p)為減少運算次數(shù)，將這個函數(shù)值保存在一個變量fp中.算法的偽代碼描述：Input:區(qū)間端點a，b；精度要求(即誤差限)&函數(shù)f(x)；最大對分次數(shù)NOutput:近似解或失敗信息行號偽代碼注釋1n—1；對分次數(shù)計數(shù)器2fa—f(a)；左端點的函數(shù)值3whilen<Ndo4p—(a+b)/2；區(qū)間中點5fp—Hp)；中點的函數(shù)值6iffp=0or(b-a)/2<￡then函數(shù)值為0或半?yún)^(qū)間長不超￡7returnp；輸出近似解并退出程序8endif9n—n+1；計數(shù)器加一10iffa-fp>0then若中點與左端點函數(shù)值同號11a—p；新區(qū)間取右半?yún)^(qū)間12fa—fp；13else否則14b—p；新區(qū)間取左右半?yún)^(qū)間

15endif16enddo17return錯誤信息輸出錯誤信息并結(jié)束程序結(jié)果與分析附錄：程序清單說明：源程序中帶有數(shù)字的空行，對應著算法描述中的行號rvf““““““““““““““““““““““““““““““““““““““““““““““““““““““““““%**個*******************************************************%*程序名：Bisection.m*%*程序功能：使用二分法求解非線性方程.rvf““““““““““““““““““““““““““““““““““““““““““““““““““““““““““%**********************************************************f=inline('xAf=inline('xA3-x-1');%定義函數(shù)f(x)a=input('有根區(qū)間左端點：a=');b=input('右端點:b=');epsilon=input('誤差限:epsilona=');N=input('最大對分次數(shù):N=');%對分次數(shù)計數(shù)器n置%對分次數(shù)計數(shù)器n置1%左端點的函數(shù)值給變量fafprintf('\nkpfprintf('b(k)f(p)a(k)f(a(k))');b-a\n');%顯示表頭fprintf('%2d%36.6f%12.6f%12.6f%12.6f\n',0,a,fa,b,b-a);%占2位其中0位小數(shù)顯示步數(shù)0,共12位其中小數(shù)6位顯示各值3%whilen<N4%取區(qū)間中點p5%求p點函數(shù)值給變量fpfprintf('%2d%12.6f%12.6f,n,p,fp);%輸出迭代過程中的中點信息p和f(p)6fprintf('\n\n近似解為:％f\n',p);return;8%如果f(p)-0或b-a的一半小于誤差限￡%則輸出近似根p(7)%并結(jié)束程序(7)%計數(shù)器加1%若f(a)與f(p)同號11%則取右半?yún)^(qū)間為新的求根區(qū)間，即a取作p12%保存新區(qū)間左端點的函數(shù)值13%否則14%左半?yún)^(qū)間為新的求根區(qū)間，即b取作p15fprintf('%12.6f%12.6f%12.6f%12.6f\n',a,fa,b,b-a);%顯示新區(qū)間端點及左端函數(shù)值、區(qū)間長度16fprintf('\n\n經(jīng)過％d次迭代后未達到精度要求.\n',N);%輸出錯誤信息(行17)《數(shù)值計算方法》實驗2報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗2非線性方程的迭代解法(之Aitken-Steffensen加速法)實驗題目用Aitken-Steffensen加速法求方程尤3+4尤2—10=0在區(qū)間［1,2］內(nèi)的一個實根，取絕對誤差限為10-4.實驗目的熟悉求解非線性方程的Aitken-Steffensen加速法.基礎理論將方程f(x)=0改寫成等價形式x=9(x)，得到從初值七開始的迭代公式氣+］=甲(xk)后，基于迭代公式xk+1=9(xk)的Aitken-Steffensen加速法是通過“迭代-再迭代-加速”完成迭代的，具體過程為J=9(x),乙=93),x=x—(七七*.

kkkkk+1kz—2y+x實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程為了驗證Aitken-Steffensen加速法可以把一些不收斂的迭代加速成迭代收斂，我們使用將方程組變形為x=1<10—x3，取迭代函數(shù)9(x)=；(10—x3，并利用宏定義出迭代函數(shù).由于不用保存迭代過程，所以用x0表示初值同時也存放前一步迭代的值，y和z是迭代過程中產(chǎn)生的恥和zk，x存放新迭代的結(jié)果.算法的偽代碼描述：Input：初值x0；精度要求(即誤差限)&迭代函數(shù)p(x)；最大迭代次數(shù)NOutput:近似解或失敗信息行號偽代碼注釋1n—1；迭代次數(shù)計數(shù)器2whilen<Ndo3y—p(x0)；迭代4z—p(y)；再迭代5x—x0—(y-x0)2/(z-2y+x0)加速6if1x一x0|<￡then如果達到精度要求7returnx；則輸出近似值并退出程序8endif

9n—n+1；計數(shù)器加一10x0—x；新近似值給x0做下次的初值11enddo12return錯誤信息輸出錯誤信息并結(jié)束程序7.結(jié)果與分析8.附錄：程序清單rvf“““““““““““““““““““““““““““““““““““““““““““““““““““““““““““““%**個**********************************************************%*程序名：Aitken_Steffensen.m%*程序功能：用Aitken-Steffensen加速法求方程.rvf“““““““““““““““““““““““““““““““““““““““““““““““““““““““““““““%*************************************************************clc;clearall;phi=inline('0.5*sqrt(10-xA3)');%迭代函數(shù)x0=input(初值:x0=');epsilon=input('誤差限：epsilon=');N=input('最大迭代次數(shù):N=');disp('n迭代中間值y(n-1)再迭代結(jié)構(gòu)z(n-1)加速后的近似值x(n)');fprintf('%2d%54.6f\n',0,x0);%占2位整數(shù)顯示步數(shù)0,為了對齊，占54位小數(shù)6位顯示x0%n是計數(shù)器%n是計數(shù)器%whilen<=N%迭代%再迭代%加速z是中間變量,x是下一步的近似值y=3;z=3;x=3;%x0初值及前一步的近似值,y和fprintf('%2d%18.6f%18.6f%18.6f\n',n,y,z,x);%顯示中間值和迭代近似值6%如果與上一步近似解差的絕對值不超過誤差限fprintf('\n\n近似解x~x(%d)~%f\n',n,x);%則輸出近似根(7)，可簡略為:fprintf('\n\n近似解x=%f,x);return;%并結(jié)束程序(7)8%相當于endif9%計數(shù)器加110%新近似值x作為下一次迭代的初值11fprintf('\n迭代％d次還不滿足誤差要求.\n\n',N);%輸出錯誤信息(12)《數(shù)值計算方法》實驗2報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗2非線性方程的迭代解法(之Newton下山法)實驗題目用Newton下山法求方程尤3+4尤2-10=0在區(qū)間［1,2］內(nèi)的一個實根，取絕對誤差限為10-4.實驗目的熟悉非線性方程的Newton下山法.基礎理論Newton下山法：Newton下山法公式為^=^-人H「，使I/(x)1<1/(x)1,k+1kk廣(x)k+1kk其中0<Xk<1.實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程定義函數(shù)f(x)和df(x)，其中df(x)是f(x)的導函數(shù).每步迭代時先取下山因子為1,嘗試迭代，判斷嘗試結(jié)果是否滿足下山因子，若滿足則作為這步的迭代結(jié)果；否則將下山因子減半，然后再嘗試.為防止當前的xk是極小值點，附近不會有滿足下述條件的其它點，使嘗試陷入死循環(huán)，同時計算機中能表示出的浮點數(shù)也有下界，因此我們設置了最大嘗試次數(shù).當超過最大嘗試次數(shù)時，不再進行下山嘗試.由于反復嘗試迭代且要判斷下山條件，所以f(x0)和f(x0)會反復使用，為避免重復計算浪費運行時間，將這兩個值分別保存在變量fx0和dfx0.而嘗試產(chǎn)生的節(jié)點，判斷下山條件時要用到它的函數(shù)值，若嘗試成功，這個點會作為下一步的初值再使用，所以把該點的函數(shù)值也保存在變量fx中.算法的偽代碼描述：Input:初值x0；精度要求(即誤差限)&函數(shù)及其導函數(shù)f(x)和f(x)；最大迭代次數(shù)N；K下山嘗試最大次數(shù)Output:近似解或失敗信息行號偽代碼注釋1n—1；迭代次數(shù)計數(shù)器2F0—f(x0)；3whilen<Ndo4F0’—g)；5ifF0’=0then若該點導數(shù)為06returnFalse；則無法進行迭代，結(jié)束程序7endif

8lambda—1；下山因子入從1試起9k—1；嘗試次數(shù)計數(shù)器10whilek<Kdo11x—x0-lambda*F0/F0’Newton下山公式12Fx—f（x）；13ifIFxl<IF0Ithen判斷下山條件14break；退出嘗試循環(huán)15endif16lambda—lambda/2；下山因子減半17k—k+1；嘗試次數(shù)計數(shù)器加118endwhile19ifk>Kthen如果因超過嘗試次數(shù)退出循環(huán)20returnFalse；則提示錯誤信息并結(jié)束程序21endif否則時嘗試成功退出上邊循環(huán)22ifIx一x0|<￡then如果達到精度要求23returnx；則輸出近似值并退出程序24endif25x0—X；新近似值給x0做下次的初值26F0—Fx；所求函數(shù)值下次也用到27n—n+1；計數(shù)器加一28enddo29return錯誤信息輸出錯誤信息并結(jié)束程序7.結(jié)果與分析8.附錄：程序清單rvf“““““““““““““““““““““““““““““““““““““““““““““““““““““““““““““%**個**********************************************************%*程序名：NewtonDownhill.m*%*程序功能：用Newton下山法求解非線性方程.m“““““““““““““““““““““““““““““““““““““““““““““““““““““““““““““%*************************************************************clc;clearall;%函數(shù)f(x)%函數(shù)f(x)的導函數(shù)f=inline（'xA3-x-1'）;df=inline（'3*xA2-1'）;x0=input（初值:x0='）;epsilon=input（'誤差限：epsilon='）;N=input（'最大迭代次數(shù)%函數(shù)f(x)%函數(shù)f(x)的導函數(shù)迭代次數(shù)計數(shù)器存x0點函數(shù)值fprintf('\n\nnx(n)f(x(n))\n');fprintf('\n\nnx(n)f(x(n))\n');顯示表頭fprintf('%2d%14.6f%14.6f\n',0,x0,fx0);%2位整數(shù)顯示0,共14位小數(shù)6位顯示x0和fx03%whilen<Ndisp('');%換行顯示下山嘗試過程的表頭disp('下山因子嘗試x(n)對應f(x(n))滿足下山條件');disp('');4%存x0點導數(shù)值，每次下山嘗試不用重新計算ifdfx0==0%導數(shù)為0不能迭代disp(‘無法進行Newton迭代’)；return;endlambda=1.0;%下山因子從1開始嘗試k=1;%k下山嘗試次數(shù)計數(shù)器whilek<=K%下山最多嘗試K次k=xO-lambda*fftO'dhtO:%下山公式fx=f(x);fprintf('%22.6f%14.6f%14.6f,lambda,x,fx);%函數(shù)值%顯示嘗試結(jié)果if(abs(fx)<abs(fx0))fprintf('滿足4');%判斷是否滿足下山條件break;elsefprintf('不滿足4');end%是，則退出下山嘗試的循環(huán)lambda=lambda/2;%不是，則下山因子減半k=k+1;endifk>K%計數(shù)器加1fprintf('\n下山條件無法滿足，迭代失敗.\n\n');return;endfprintf('%2d%14.6f%14.6f\n',n,x,fx);%2位整數(shù)顯示步數(shù)n,共14位小數(shù)6位顯示下步迭代結(jié)果22%達到精度要求否fprintf('\n\n方程的近似解為：x^f\n',x);%(23)return;%達到，則顯示結(jié)果并結(jié)束程序(23)22end%(24)%用x0,fx0存放前一步的近似值和它的函數(shù)值，進行循環(huán)迭代25262728fprintf('\n迭代％d次還不滿足誤差要求.\n\n',N);《數(shù)值計算方法》實驗2報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗2非線性方程的迭代解法（之弦截法）實驗題目用弦截法求方程X3+4X2-10=0在區(qū)間［1,2］內(nèi)的一個實根，取絕對誤差限為10-4.實驗目的熟悉非線性方程的弦截法.基礎理論將Newton迭代法中的導數(shù)用差商代替，得到弦截法（或叫正割法）公式X=X一-k__-^f(X).*kf(X)-f(X)kkk-1實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程不保存迭代過程，所以始終以X0和X1分別存放Xk-1和Xk，而x存放新產(chǎn)生的迭代值Xk+1，這樣，下一次迭代時需要把上一步的X1（即Xk）賦值于X0（做新的Xk-1）.這些點的函數(shù)值會重復用到，在迭代公式中也要用到，上一步的X］作為下一步的X0也會再一次用它的函數(shù)值，為減少重新計算該點函數(shù)值的運行時間，將X］點的函數(shù)值保存在變量fx1中.算法的偽代碼描述：Input:初值x0，x1；精度要求（即誤差限）￡；函數(shù)f（x）；最大迭代次數(shù)NOutput:近似解或失敗信息行號偽代碼注釋1fx0—f(x0)；初值點的函數(shù)值2fx1—f(x,)；3n—2；迭代次數(shù)計數(shù)器（從2開始）4whilen<Ndo5x—x1-(x1-x0)*f(x1)/(f(x,)-f(x0))；迭代6if1x一x0|<￡then如果達到精度要求7returnx；則輸出近似值并退出程序8endif9x0—x1；（按順序）最新的2個近似值，10x1—x；分別做下次的初值x0，x111fx0—fx1；下次初值點的函數(shù)值12fx1—f(x1)；13n—n+1；14enddo15return錯誤信息輸出錯誤信息并結(jié)束程序7.結(jié)果與分析8.附錄：程序清單%**個**********************************************************%*程序名：SecantMethod.m*%*程序功能：用弦截法求解非線性方程.%*************************************************************clc;clearall;f=inline('2*xA3-5*x-1');%函數(shù)f(x)x0=input('第一初值:x0=');x1=input('第二初值:x1=');epsilon=input('誤差限：epsilon=');N=input('最大迭代次數(shù):N=');fprintf('\nnx(n)\n');fprintf('%2d%14.6f\n',0,x0);fprintf('%2d%14.6f\n',1,x1);%顯示表頭%占2位顯示步數(shù)0,共14位其中小數(shù)6位顯示x0%占2位顯示步數(shù)1,共14位其中小數(shù)6位顯示x11%存x0點函數(shù)值2%存x1點函數(shù)值3%迭代計數(shù)器4%whilen<N:―"三土頊%弦截法公式fprintf（'%2d%14.6f\n',n,x）;%顯示迭代過程6%達到精度要求否fprintf（'\n\n方程的近似解為:x~%f\n\n',x）;return;%達到，則顯示結(jié)果并結(jié)束程序89%原x1做x0為前兩步的近似值W%現(xiàn)x做x1為一兩步的近似值11%x0點函數(shù)值12%計算x1點函數(shù)值，為下一次循環(huán)13%計數(shù)器加114fprintf（'\n迭代％d次還不滿足誤差要求.\n\n',N）;

《數(shù)值計算方法》實驗3報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗3解線性方程組的直接法（之Gauss消去法）實驗題目用Gauss消去法求解線性方程組'0.0012.000-1.000'0.0012.000-1.0003.712廠2.0001.072實驗目的掌握解線性方程組的Gauss消去法.基礎理論3.000卜1]4.623x25.643人x3J[L000]=2.0000.000/Gauss消去法是通過對增廣矩陣的初等行變換，將方程組變成上三角方程組，然后通過回代，從后到前依次求出各未知數(shù).Gauss消去法的第k步（1MkMn-1）消元:若九豐0,則依次將增廣矩陣第k行的-七/九倍加到第1行（k+1<i<n）,將第k列對角線下的元素都化成0.實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程結(jié)果與分析附錄：程序清單

《數(shù)值計算方法》實驗3報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗3解線性方程組的直接法（之Gauss列主元消去法）實驗題目用Gauss列主元消去法求解線性方程組'0.001-1.000廠2.0002.0003.7121.072'0.001-1.000廠2.0002.0003.7121.0723.000卜1]4.623x25.643人x3J[L000]=2.0000.000/掌握解線性方程組的Gauss列主元消去法.基礎理論Gauss列主元消去法也是通過對增廣矩陣的初等行變換，將方程組變成上三角方程組，然后通過回代，從后到前依次求出各未知數(shù).Gauss列主元消去法的第k步（1MkMn-1）消元:先在九,a^,,ak中找絕對值最大的，將它所在的行與第k行交換，然后將第k行的-七/%倍加到第i行（k+1<i<n），將第k列對角線下的元素都化成0.實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程結(jié)果與分析附錄：程序清單

《數(shù)值計算方法》實驗3報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗3解線性方程組的直接法（之Doolittle分解）實驗題目對矩陣A進行Doolittle分解，其中(621-1、2410A=114-1[-10-13；3.實驗目的掌握矩陣的Doolittle分解.4.基礎理論矩陣的Doolittle分解是指將矩陣A=（七）呻可以分解為一個單位下三角矩陣和一個上三角矩陣的乘積.若設(100…0、'uuu…u、1112131nl10…00uu…u2122232nL=ll1…0,U=00u…u3132333n?????.....?----.?...???.?12ln3…1^V000…u/nn則可依如下順序公式計算1「九=akj-如l*r=1j=k,k+1,…，n1l=(aikik—Jlu)/u,ittkkkt=1i=k+1,k+2,…，n其中k=1,2,…，n.5.實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程（1）按計算公式依次計算一行u同時計算一列l(wèi)；（2）因為計算完u..（或l..）后，a..就不再使用，為節(jié)省存儲空間，將計算的u..（和l..）ijljljlj.j仍存放在矩陣A中的相應位置；（3）使用L矩陣和U矩陣時需要根據(jù)元素所在位置取固定值或A中相應位置的值.L對角線上的元素為1,上三角部分為0,下三角部分為A中對應的元素；U的下三角部分為0,上三角部分為A中對應的元素.

算法的偽代碼描述：Input:階數(shù)n；矩陣AOutput：矩陣L和U（合并存儲在數(shù)組A中）行號偽代碼注釋1forkj1ton2forjjkton3u=a-*11u,kjkjkrrjr=1計算U的第k行（計算結(jié)果存放在aj4endfor5forijkton6l=(a-^11u)/uikikittkkkt=1計算L的第k列（計算結(jié)果存放在aj7endfor8endfor9returnL,U輸出L和U結(jié)果與分析附錄：程序清單rvf““““““““““““““““““““““““““““““““““““““““““““““““““““%**個*************************************************%程序名:Doolittle.m*%程序功能：矩陣LU分解中的Doolittle分解.rvf““““““““““““““““““““““““““““““““““““““““““““““““““““%****************************************************clc;%矩陣階數(shù)clearall;n=4;A=[621-1;2410;114-1;-10-13]disp('A=');disp(A);%LU分解(Doolittle分解)fork=1:n%計算矩陣U的元素u_{kj}%（可參照下面l_{ik}的公式填寫%矩陣階數(shù)知i=k-1:nfcri=1:k-1Ag=A(ik)-A&n項均％計算矩陣L的元素l_(ik}endend%L在A下三角,U在上三角(對角線為1)enddisp('分解結(jié)果：');disp('L=');fori=1:nforj=1:nifi>j%在下三角部分，則取A對于的元素顯示fprintf('%8.4f,A(i,j));elseifi==j%在對角線上，則顯示1fprintf('%8d',1);else%在上三角部分，則顯示0fprintf('%8d',0);endendfprintf('\n');%換行enddisp('U=');fori=1:nforj=1:nifi<=j%在上三角部分或?qū)蔷€上，則取A對于的元素顯示fprintf('%8.4f,A(i,j));else%在下三角部分，則顯示0fprintf('%8d',0);endendfprintf('\n');%換行end《數(shù)值計算方法》實驗3報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗3解線性方程組的直接法（之LU分解法）實驗題目用LU分解（Doolittle分解）法求解線性方程組x+2x+6x=1<2x+5x+15x=36x'+15；+46【=10

、123實驗目的熟悉解線性方程組LU分解法.基礎理論、若將矩陣A進行了Doolittle分解，A=LU,則解方程組杰=b可以分解求解兩個三角方程組Ly=b和Ux=M它們都可直接代入求解，其中Ly=b的代入公式為y=b一云ly,k=1,2,…,nj=1而ux=y的代入公式為x=（y-^Lux）/u,k=n,n-1,?…，1.kkkjjkkj=k+1實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程（1）Doolittle分解過程依次計算一行u同時計算一列l(wèi)完成，并將計算的勺（和l〃）仍存放在矩陣A中的相應位置；""（2）求解方程組的代入公式中用到的勺和1〃都直接在A的相應位置取值即可.算法的偽代碼描述：""Input：階數(shù)n；矩陣A；常數(shù)項向量bOutput:方程組Ax=b的解x行號偽代碼注釋1LU分解參加上一實驗程序2fork,1ton3y=b一支11y,k=1,2,…，nj=1代入法求解Ly=b4endfor

5fork,1ton6x=(y—才ux)/u,k=n,n-1,,1kkkjjkkj=k+1回代求解Ux=y7endfor8returnX],…,xn7.結(jié)果與分析8.附錄：程序清單rvf$$$$$$$力e“$$$$$$$$$$$$$$““$$$$$$$$e%****************************************************%程序名：LinearSystemByLU.m*%程序功能：利用LU分解(Doolittle分解)解方程組.rvf$$$$$$-力、e“$$$$$$$$$$$$$$““$$$$$$$$e%****************************************************clc;%矩陣階數(shù)clearall;n=3;A=[126;2515;61546];b=[1;3;10];%LU分解(Doolittle分解)fork=1:n%計算矩陣U的元素u_{kj}%（可參照下面l_{ik}的公式填寫%矩陣階數(shù)融i=k-1:nfcrj=1:k-1%計算矩陣L的元素l_{ik}AUk)=A(就-A(i:ij*A(j,k);endAMk)=AMk)'AfKk):%計算矩陣L的元素l_{ik}end%L在A下三角，U在上三角（對角線為1）endfork=1:ny(k)=b(k);%用代入法求解下三角方程組fork=1:ny(k)=b(k);%用代入法求解下三角方程組Ly=b%七=bk-11/.,k=1,2,…,nj=13enddisp('方程組Ly=b的解：y=');disp(y');fork=n:-1:1%回代求解上三角方程組Ux=y一.寸%Xk=y一乙x(k)=y(k);ux,k=n,n-1,?…,1j=k+1enddisp('原方程組的解：x=');一.寸%Xk=y一乙《數(shù)值計算方法》實驗3報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗3解線性方程組的直接法(之Cholesky分解)實驗題目對矩陣A進行Cholesky分解，其中f621-1]A=241°.114-1"-10-13J實驗目的理解矩陣的Cholesky分解.基礎理論矩陣的Cholesky分解是指將矩陣A=(%)心可以分解為一個下三角矩陣L和L轉(zhuǎn)置的乘積，即A=LLt，其中L各元素可依如下順序公式計算I七="_蕓%＜》r=1l=(a-習ll)/1,i=k+1,k+2,,nikikittkkkIt=1其中k=1,2,...,n.實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：VC++實驗過程按計算公式依次先計算一列對角線上的元素1kk，再計算這列其他元素履且對稱位置的元素也取同一個值；因為計算完1..后，a?就不再使用，為節(jié)省存儲空間，將計算的1..仍存放在矩陣A中的相應位置；““"使用L矩陣時需要根據(jù)元素所在位置取固定值或A中相應位置的值.L上三角部分為0，對角線和下三角部分為A中對應的元素.算法的偽代碼描述：Input:階數(shù)n；矩陣AOutput：矩陣L(合并存儲在數(shù)組A中)行號偽代碼注釋1fork—1ton

2l=\a-t^112kk\|kkkr'r=13fori,kton41=(a-t111)/1ikikittkkkt=1計算結(jié)果存放在a..5endfor6endfor7returnL輸出L結(jié)果與分析附錄：程序清單%**個**********************************************************%*程序名：Cholesky.m*%*程序功能：對稱正定矩陣的Cholesky分解.%*************************************************************n=4;%矩陣階數(shù)A=[6,2,1,-1;2,4,1,0;1,1,4,-1;-1,0,-1,3];disp('A=');fori=1:nforj=1:nfprintf('%10.4f,A(i,j));%共占14位endfprintf('\n');%一行結(jié)束換行endfork=1:n%計算對角線上的Zfork=1:n%計算對角線上的Z_{kk}%計算其他的匚{ik}%和匚{ki}for戶1:虹]A(k,k)=A-(l￡nk)-A(Kjr2：endA'tki=，qrt[A[哄由血i=k+l:nfcrj=l:k-lA(i,k)=A(i,k)-AM+A03k)：endAftk)=Aftt)/A(ktXA(Ki)=A(i,k)：end

end%L在A下三角，Lend%L在A下三角，LAT在上三角forj=1:nifi>=j%在下三角部分或?qū)蔷€上，則取A對于的元素顯示fprintf('%10.4f,A(i,j));else%在上三角部分，則顯示0fprintf('%10d',0);endendfprintf('\n');%換行end《數(shù)值計算方法》實驗3報告班級：20xx級XXXXx班學號：20xx2409xxxx姓名：XXX成績：實驗名稱實驗3解線性方程組的直接法(之改進的Cholesky分解)實驗題目對矩陣A進行改進的Cholesky分解，其中(621-1)2410A=114-1-1V10-13°7實驗目的理解矩陣改進的Cholesky分解.基礎理論矩陣的改進的Cholesky分解是指將矩陣A=(%)心可以分解為一個單位下三角矩陣L和對角矩陣D及L轉(zhuǎn)置的乘積，即A=LDLt，其中L和D各元素可依如下順序公式計算d=a一如d12kkkrkrr=11=(a-切dl1)/d,i=k+1,k+2,,nikiktitktkIt=1其中k=1，2，…，n.5.實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：VC++6.實驗過程按計算公式依次先計算D的一個元素dk，再計算L中這列的元素匕/且對稱位置的元素也取同一個值；因為計算完dk和/〃后，akk或％就不再使用，為節(jié)省存儲空間，將計算的％或1..仍存放在矩陣A中的相應位置；使用L矩陣時需要根據(jù)元素所在位置取固定值或A中相應位置的值.L對角線和上三角部分為0，下三角部分為A中對應的元素；D對角線為A中對應的元素，其余都是0.算法的偽代碼描述：Input:階數(shù)n；矩陣AOutput：矩陣L(合并存儲在數(shù)組A中)行號偽代碼注釋1forkj1ton

2d=a一切d12r=13fori,kton41=(a-云dl1)/dikiktitktkt=1計算結(jié)果存放在a--IJ5endfor6endfor7returnL,D輸出L和D7.結(jié)果與分析8.附錄：程序清單%**個**********************************************************%*程序名：ImprovedCholesky.m*%*程序功能：對稱正定矩陣的改進的Cholesky分解.%*************************************************************n=4;%矩陣階數(shù)A=[6,2,1,-1;2,4,1,0;1,1,4,-1;-1,0,-1,3];disp('A=');fori=1:nforj=1:nfprintf('%10.4f,A(i,j));%共占14位endfprintf('\n');%一行結(jié)束換行end%Cholesky分解fork=1:n%計算D對角線上的w_{kk}%計算L的元素Z_{ik}%和L轉(zhuǎn)置的元素Z_{ki}%計算D對角線上的w_{kk}%計算L的元素Z_{ik}%和L轉(zhuǎn)置的元素Z_{ki}endfori=k+l:nA(U)=A{U)-Ag*削ij)*endA(Lk)=A(iX)..0?k〕；A(ti)=A(U)；end

end%L在A下三角,Dend%L在A下三角,D在對角線forj=1:nifi>j%在下三角部分，則取A對于的元素顯示fprintf('%10.4f,A(i,j));elseifi==j%在對角線上，則顯示1fprintf('%10d',1);else%在上三角部分，則顯示0fprintf('%10d',0);endendfprintf('\n');%換行enddisp('D=');fori=1:nforj=1:nifi==j%在對角線上，則取A對于的元素顯示fprintf('%10.4f,A(i,j));else%其余顯示0fprintf('%10d',0);endendfprintf('\n');%換行end《數(shù)值計算方法》實驗3報告班級：20xx級XXXXx班姓名：XXX學號：20xx2409xxxx成績：對于系數(shù)矩陣為三對角矩陣的方程組其Crout分解可分解為(b1a2C1b2a對于系數(shù)矩陣為三對角矩陣的方程組其Crout分解可分解為(b1a2C1b2an-1n-1ancn-1bnJ(s1a2asn-1n-1asnnt11t2實驗3解線性方程組的直接法（之追趕法）實驗題目用追趕法求解線性方程組-13-10實驗目的熟悉解線性方程組的追趕法.基礎理論這樣，解方程組可以由如下2步完成:“追”：s=s=b,t=c/s,s=b—at,i-1i-1i-1iiii-1"/s1，J=(f—ay)/s,i=2,3,???,n,iiii-1i其中：（匕，…，頃為方程組的常數(shù)項，八沒用；“趕”：Xn=yn,x=y-1x,i=n-1,n-2,—,1.iiii+1實驗環(huán)境操作系統(tǒng)：Windowsxp；程序設計語言：Matlab實驗過程在“追”的過程中，向量s和y都有n個元素，t只有n-1個元素，又七和"的計算公式與其它s和y不同，所以先單獨計算s和y，然后在一個n-1次循環(huán)中，求其它s和y以ii11ii

由于在“追”的過程中，b,c和f在分別計算完對應的s,t.和y后就不再使用，所iiiiIi以借用數(shù)組b，c和f存儲向量s，t和y；同樣在“趕”的過程中，七在計算完對應的X后就不再使用，所以再一次借用數(shù)組f存儲向量x追趕法算法的偽代碼描述：Input:階數(shù)n；三對角矩陣的三條對角線向量a，b，c，常數(shù)項向量fOutput:方程組的解x行號偽代碼注釋1s1。b1；追2y】。f1/s1；3forkj2ton4tk-1jck-1/sk-1；5skjbk-akL；6ykj(fk-akyk-1)/sk；7endfor8