數(shù)值分析實驗報告

1、數(shù)值分析實驗報告姓 名鄒昊學 號20081138班 級軟一0801指 導 教 師馮男實驗名稱數(shù)值分析開 設(shè) 學 期2009 2010第二學期評 定 成 績評定人簽字評 定 日 期東北大學軟件學院2010年5月實驗一一、課題名稱解線性方程組的迭代法二、目的和意義 了解各迭代法的基本原理和特點，判斷雅克比迭代、高斯-塞德爾迭代對任意初始向量的收斂性，完成雅克比迭代、高斯-塞德爾迭代算法的程序?qū)崿F(xiàn)三、計算公式l 雅可比 xi(k+1)=1/aii(bi-aijxj(k)l 高斯-塞德爾xi(k+1)=1/aii(bi-aijxj(k+1)-aijxj(k) l 超松弛迭代 xi(k+1)=（1-w）

2、xi(k)+w*(bi-aijxj(k+1)-aijxj(k) /aii四、結(jié)構(gòu)程序設(shè)計1、主程序流程圖 圖1-1 主程序流程圖2、jacobi迭代算法流程：圖1-2 jacobi 迭代法流程圖3、gauss-seidel迭代算法流程同jacobi算法:4、sor迭代算法：圖1-4 sor迭代法流程圖5、代碼:#includeusing namespace std; #define n 40const int n=10;int jacobi(float *p,float b,float x,float x,int n);int gs(float *p,float b,float x,float

3、 x,int n);int sor(float *p,float b,float x,float x,int n);void print(float *a,int r);void main()float a1010=4,2,-3,-1,2,1,0,0,0,0, 8,6,-5,-3,6,5,0,1,0,0,4,2,-2,-1,3,2,-1,1,9,4,0,-2,1,5,-1,3,-1,1,9,4,-4,2,6,-1,6,7,-3,3,2,3,8,6,-8,5,7,17,2,6,-3,5,0,2,-1,3,-4,2,5,3,0,1,16,10,-11,-9,17,34,2,-1,2,2,4,6,2

4、,-7,13,9,2,0,12,4,0,0,-1,8,-3,-24,-8,6,3,-1; float a10=7,5,-13,2,6,-12,14,-4,5,-5; float x110=0,0,0,0,0,0,0,0,0,0; float x110; float b88=4,2,-4,0,2,4,0,0, 2,2,-1,-2,1,3,2,0, -4,-1,14,1,-8,-3,5,6, 0,-2,1,6,-1,-4,-3,3, 2,1,-8,-1,22,4,-10,-3, 4,3,-3,-4,4,11,1,-4, 0,2,5,-3,-10,1,14,2,0,0,6,3,-3,-4,2,19;

5、 float b8=0,-6,6,23,11,-22,-15,45; float x28=0,0,0,0,0,0,0,0; float x28; float c1010=4,-1,0,0,0,0,0,0,0,0, -1,4,-1,0,0,0,0,0,0,0, 0,-1,4,-1,0,0,0,0,0,0, 0,0,-1,4,-1,0,0,0,0,0, 0,0,0,-1,4,-1,0,0,0,0, 0,0,0,0,-1,4,-1,0,0,0, 0,0,0,0,0,-1,4,-1,0,0, 0,0,0,0,0,0,-1,4,-1,0, 0,0,0,0,0,0,0,-1,4,-1, 0,0,0,0,0

6、,0,0,0,-1,4; float c10=7,5,-13,2,6,-12,14,-4,5,-5; float x310; float x310=0,0,0,0,0,0,0,0,0,0; float *p3; p0=&a00; p1=&b00; p2=&c00; coutjacobi迭代法解第一個方程:endl; jacobi(p0,a,x1,x1,10); coutjacobi迭代法解第二個方程:endl; jacobi(p1,b,x2,x2,8); coutjacobi迭代法解第三個方程:endl; jacobi(p2,c,x3,x3,10); coutgauss-seidel迭代法解第

7、一個方程:endl; gs(p0,a,x1,x1,10); coutgauss-seidel迭代法解第二個方程:endl; gs(p1,b,x2,x2,8); coutgauss-seidel迭代法解第三個方程:endl; gs(p2,c,x3,x3,10); coutsor迭代法解第一個方程:endl; sor(p0,a,x1,x1,10); coutsor迭代法解第二個方程:endl; sor(p1,b,x2,x2,8); coutsor迭代法解第三個方程:endl; sor(p2,c,x3,x3,10); int jacobi(float *p,float b,float x,float

8、 x,int n) int k,i,j; float m,r,r,e; coute; for(k=0;kn;k+) r=0; for(i=0;in;i+) m=0; for(j=0;jn;j+) m=m+(*(p+i*n+j)*xj; xi=xi+(bi-m)/(*(p+i*n+i); r=xi-xi; if(rr) r=r; if(r0) print(x,n); cout迭代次數(shù)為:k+1endl; return k; for(j=0;j10;j+) xj=xj; print(x,n);cout迭代次數(shù)為:kendl; cout方程解發(fā)散,無法用jacobi方法解此方程!endl; retu

9、rn 0;int gs(float *p,float b,float x,float x,int n) int i,j,k;float t,r,r,e; coute; for(k=0;kn;k+) for(i=0;in;i+) t=0; for(j=0;jn;j+) if(ji) t+=(*(p+i*n+j)*xj; xi=(bi-t)/(*(p+i*n+i); for(i=0;i10;i+) r=xi-xi; if(rr) r=r; if(r0) print(x,n); cout迭代次數(shù)為:k+1endl; return k; for(j=0;j8;j+) xj=xj; print(x,n)

10、;cout迭代次數(shù)為:kendl; cout方程解發(fā)散,無法用gauss-seidel方法解此方程!endl; return 0;int sor(float *p,float b,float x,float x,int n)int i,j,k;float t,r,r,e,w; cout請輸入松弛因子w(0ww;coute; for(i=0;in;i+) xi=xi; for(k=0;kn;k+) r=0; for(i=0;in;i+) t=0; for(j=0;jn;j+) t+=(*(p+i*n+j)*xj; r=w*(bi-t)/(*(p+i*n+i); xi+=r; if(rr) r=r

11、; if(r0) print(x,n); cout迭代次數(shù)為:k+1endl; return k; print(x,n);cout迭代次數(shù)為:kendl;cout方程解發(fā)散,無法用xor方法解此方程!endl; return 0;void print(float *a,int n) int j; float *t=a; coutx=( ; for(j=0;jn-1;j+) cout*(t+j),; cout*(t+j)endl;五、結(jié)果討論和分析1.程序截圖：jacobi: gauss-seidel: sor: 2.算法總結(jié)：與直接法相比，迭代法適用于稀疏矩陣的線性方程組。在試驗中，不同的系數(shù)

12、矩陣對上述三種迭代方法有很大影響，會導致結(jié)果發(fā)散無法得到正常結(jié)果。三種算法的收斂，sor方法最大，其次gauss-seidel方法，jacobi方法最小，松弛因子的不同也和收斂速度密切相關(guān)。實驗二一、課題名稱曲線擬合的最小二乘法二、目的和意義掌握曲線擬合的最小二乘法；了解最小二乘法亦可以用于解超定線性方程組；探索擬合函數(shù)的選擇與擬合精度間的關(guān)系。三、計算公式e 22=2i=(xi)-f(xi)2四、結(jié)構(gòu)程序設(shè)計#include stdio.h#include math.h#define num 10float neiji(float bnum,float cnum) int p;float n

13、j=0;for(p=1;pnum;p+)nj+=cp*bp;return nj;float snum,xnum,fainumnum,afanum;float beidanum,anum,xfainum,ydnum,max,pcpfh;void main() int i,j,k,n,index,flag;char conti;conti= ;printf(請輸入已知點的個數(shù)n=n);scanf(%d,&n);printf(請輸入x和y:);for(i=1;i=n;i+) printf(x%d=,i);scanf(%f,&xi);printf(y%d=,i);scanf(%f,&yi);while

14、(conti= ) printf(請輸入擬和次數(shù)=);scanf(%d,&index);pcpfh=0;afa1=0;a0=0;for(i=1;i=n;i+) afa1+=xi;a0+=ydi;fai0i=1;afa1=afa1/n;a0=a0/n;for(i=1;i=n;i+)fai1i=xi-afa1;a1=neiji(fai1,yd)/neiji(fai1,fai1);for(k=1;kindex;k+) for(i=1;i=n;i+)xfaii=xi*faiki;afak+1=neiji(faik,xfai)/neiji(faik,faik);beidak=neiji(faik,fai

15、k)/neiji(faik-1,faik-1);for(j=1;j=n;j+)faik+1j=(xj-afak+1)*faikj-beidak*faik-1j;ak+1=neiji(faik+1,yd)/neiji(faik+1,faik+1);printf(%d次擬和結(jié)果為n,index);for(i=0;i=index;i+)printf(a%d=%fn,i,ai);for(i=1;i=index;i+)printf(afa%d=%fn,i,afai);for(i=1;iindex;i+)printf(beida%d=%fn,i,beidai);for(i=1;i=n;i+) for(k=

16、0;k=index;k+)si+=ak*faiki;ydi=fabs(ydi-si);pcpfh+=ydi*ydi;si=0;max=0;for(i=1;imax)max=ydi;flag=i;printf(當x=%f時,偏差最大=%f,偏差平方和為%fn,xflag,max,pcpfh);printf(繼續(xù)擬和請按space,按其他鍵退出);conti=getchar();conti=getchar();五、結(jié)果討論和分析1.程序截圖：請輸入已知點的個數(shù)n=10請輸入x和y：x1=0y1=0請輸入擬合次數(shù)=55次擬合結(jié)果為當x=0.000000時，偏差最大=6706185.000000，偏差

17、平方和為449729146126336.000000.實驗三一、課題名稱數(shù)值積分二、目的和意義1、 深刻認識數(shù)值積分法的意義；2、 明確數(shù)值積分精度與步長的關(guān)系；3、 根據(jù)定積分的計算方法，可以考慮二重積分的計算問題。三、計算公式選用復(fù)合梯形公式，復(fù)合simpson公式，romberg算法，計算（1） i = （2） i = （3） i = （4） i = 四、結(jié)構(gòu)程序設(shè)計1、主程序流程圖圖2-1 主程序流程圖2.復(fù)合梯形公式的計算方法圖2-2 復(fù)合梯形公式程序流程圖在此算法中，主要思想是將區(qū)間端點的兩個值相加，然后將區(qū)間間的函數(shù)值乘以2進行相加，然后總體和乘以h/2.3.復(fù)合simpson公

18、式計算方法圖2-3 .復(fù)合simpson算法流程圖在此方法中，主要思想和上述梯形復(fù)合公式的算法相似，只是有些項加的數(shù)不同，就不具體說明了。4.romberg算法計算公式圖2-4 romberg算法在此算法中，利用一個二維數(shù)組分別存儲xk和xk+1，只要記錄下所謂t數(shù)表的上一行，就能得出下一行即xk+1的值，進而可根據(jù)t0j和t1j計算出t1j+1，根據(jù)遞推可以得出t1k-2。5.程序代碼#include #include math.husing namespace std;#define n 100float simpson(float a,float b);float romberg(flo

19、at a,float b);float trapezium(float a,float b);void trans(float (*p)(float),float a,float b,int n);void transromberg(float (*p)(float),float a,float b,int m);float f1(float x);float f2(float x);float f3(float x);float f4(float x);float fn;float (*fp)(float);int main() char m; float a,b; float i; cou

20、t輸入1選擇函數(shù)1!endl輸入2選擇函數(shù)2!endl輸入3選擇函數(shù)3!endl輸入4選擇函數(shù)4!endl其他退出!m; switch(m) case 1:fp=f1;a=0;b=0.785;cout計算函數(shù)1的積分:endl;break; case 2:fp=f2;a=0;b=1;cout計算函數(shù)2的積分:endl;break; case 3:fp=f3;a=0;b=1;cout計算函數(shù)3的積分:endl;break; case 4:fp=f4;a=0;b=1;cout計算函數(shù)4的積分:endl;break; default:exit(1); i=simpson(a,b); couti=ie

21、ndl; i=trapezium (a,b); couti=iendl; i=romberg(a,b); couti=iendl; return 0;float simpson(float a,float b) int j,n; float x,h; float f0,f1,f2,sn; cout復(fù)化simpson求積公式,請輸入?yún)^(qū)間劃分n.endln; trans(fp,a,b,2*n); h=(b-a)/(2*n); f0=f0+f2*n; f1=0; f2=0; for(j=1;j(2*n);j+) x=a+j*h;if(j%2=0) f2+=fj;elsef1+=fj; sn=(f0+

22、4*f1+2*f2)*h/3; return sn;float trapezium(float a,float b) int j,n; float x,h; float f0,f1,sn; cout復(fù)化梯形求積公式,請輸入?yún)^(qū)間劃分n.endln; trans(fp,a,b,n); h=(b-a)/n; f0=f0+fn; f1=0; for(j=1;jn;j+) x=a+j*h;f1+=fj; sn=(f0+2*f1)*h/2; return sn;void trans(float (*p)(float),float a,float b,int n) int i;float h;h=(b-a)/n; for(i=0;i=n;i+) fi=p(a+i*h);void transromberg(float (*p)(float),float a,float b,int m) int i,n,k;float t,h;n=1;h=b-a;f0=p(a)+p(b); for(k=1;k=m;k+) h/