数值分析第一次大作业 .docx

《数值分析第一次大作业 .docx》由会员分享，可在线阅读，更多相关《数值分析第一次大作业 .docx（7页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、数值分析第一次大作业王昕飏 一、 算法的设计方案A中多处为0，所以可将A501501以压缩存储为C5501来节省空间。由题中120，则501为a，1=b，否则相反。由于s是A中按模最小的特征值，所以可由反幂法得到lamdas。采用doolittle分解法将A分解成L、U。A的条件数cond(A)2=|A的按模最大特征值/A的按模最小特征值|=|a/lamdas|。在1中将A进行LU分解时，得到的矩阵U的对角线乘积就是A的行列式。求与每个k最接近的特征值，通过原点平移的方法，B=A-k*I。则此时与其最接近的特征值ik变成B按模最小的特征值，于是通过反幂法求解再加上k于是就得到ik。二、 源代码

2、#include#include#include#includeusing namespace std;int Max(int value1,int value2);int Min(int value1,int value2);void Transform(double A5501);double mifa(double A5501);void daizhuangdoolite(double A5501,double x501,double b501);double fanmifa(double A5501);double Det(double A5501);/*定义2个判断大小的函数便于以后

3、调用*/int Max(int value1,int value2) return(value1value2)?value1:value2);int Min(int value1,int value2) return (value1value2)?value1:value2);void Transform(double A5501,double b,double c) int i=0,j=0; Aij=0,Aij+1=0; for(j=2;j=500;j+) Aij=c; i+; j=0; Aij=0; for(j=1;j=500;j+) Aij=b; i+; for(j=0;j=500;j+

4、) Aij=(1.64-0.024*(j+1)*sin(0.2*(j+1)-0.64*exp(0.1/(j+1); i+; for(j=0;j=499;j+) Aij=b; Aij=0; i+; for(j=0;j=498;j+) Aij=c; Aij=0,Aij+1=0;double mifa(double A5501) int s=2,r=2,m=0,i,j; double b2,b1=0,sum,u501,y501; for (i=0;i=500;i+) ui = 1.0; do sum=0; if(m!=0)b1=b2; m+; for(i=0;i=500;i+) sum+=ui*ui

5、; for(i=0;i=500;i+) yi=ui/sqrt(sum); for(i=0;i=500;i+) ui=0; for(j=Max(i-r,0);j=Min(i+s,500);j+) ui=ui+Ai-j+sj*yj; b2=0; for(i=0;i=exp(-12); return b2;/带状DOOLITE分解并且求解出方程组的解void daizhuangdoolite(double A5501,double x501,double b501) int i,j,k,t,s=2,r=2; double B5501,c501; for(i=0;i=4;i+) for(j=0;j=5

6、00;j+) Bij=Aij; for(i=0;i=500;i+) ci=bi; for(k=0;k=500;k+) for(j=k;j=Min(k+s,500);j+) for(t=Max(0,Max(k-r,j-s);t=k-1;t+) Bk-j+sj=Bk-j+sj-Bk-t+st*Bt-j+sj; for(i=k+1;i=Min(k+r,500);i+) for(t=Max(0,Max(i-r,k-s);t=k-1;t+) Bi-k+sk=Bi-k+sk-Bi-t+st*Bt-k+sk; Bi-k+sk=Bi-k+sk/Bsk; for(i=1;i=500;i+) for(t=Max(

7、0,i-r);t=0;i-) xi=ci; for(t=i+1;t=Min(i+s,500);t+) xi=xi-Bi-t+st*xt; xi=xi/Bsi; /用于求解模最大的特征值反幂法double fanmifa(double A5501) int s=2,r=2,m=0,i; double b2,b1=0,sum=0,u501,y501; for (i=0;i=500;i+) ui = 1.0; do if(m!=0)b1=b2; m+; sum=0; for(i=0;i=500;i+) sum+=ui*ui; for(i=0;i=500;i+) yi=ui/sqrt(sum); da

8、izhuangdoolite(A,u,y); b2=0; for(i=0;i=fabs(b1)*exp(-12); return 1/b2;/行列式的LU分解U的主线乘积即位矩阵的DETdouble Det(double A5501) int i,j,k,t,s=2,r=2; for(k=0;k=500;k+) for(j=k;j=Min(k+s,500);j+) for(t=Max(0,Max(k-r,j-s);t=k-1;t+) Ak-j+sj=Ak-j+sj-Ak-t+st*At-j+sj; for(i=k+1;i=Min(k+r,500);i+) for(t=Max(0,Max(i-r

9、,k-s);t=k-1;t+) Ai-k+sk=Ai-k+sk-Ai-t+st*At-k+sk; Ai-k+sk=Ai-k+sk/Ask; double det=1; for(i=0;i0) /判断求得最大以及最小的特征值.如果K0则它为最大特征 /并以它为偏移量再用一次幂法求得新矩阵最大特征值即为最大 /与最小的特征值的差 lamda501=k; for(i=0;i=500;i+) A2i=A2i-k; lamda1=mifa(A)+lamda501; for(i=0;i=500;i+) A2i=A2i+k; else /如果K=0则它为最小特征值值并以它为偏移量再用一次 /求得新矩阵最大特

10、征值即为最大与最小的特征值的差 lamda1=k; for(i=0;i=500;i+) A2i=A2i-k; lamda501=mifa(A)+lamda1; for(i=0;i=500;i+) A2i=A2i+k; lamdas=fanmifa(A); printf(1=%.12en,lamda1); printf(501=%.12en,lamda501); printf(s=%.12enn,lamdas); printf(t要求接近的值ttt实际求得的特征值n); /*for(i=1;i=39;i+) /反幂法求得与给定值接近的特征值 p=lamda1+(i+1)*(lamda501-lamda1)/40; for(j=0;j=500;j+) A2j=A2j-p; q=fanmifa(A)+p; for(j=0;j|2|n|为其对应特征值。设0为初始迭代向量=a1*X1+a2*X2+an*Xn。k=1k*a1*X1+a2*(2/1)k*X2+an*(n/1)k*Xn。如果初始迭代向量0选取使得其在特征值1上的分量为零，但是由于计算中存在误差，迭代次数充分大后理论上使得k在特征值1上的分量不为零于是以k为新的初始迭代向量仍可以得到结果。但是实际中由于精度的要求和算法允许的计算量无法做到迭代次数充分大，于是在实际选取中初始迭代向量最好使得其在按模最大特征值的分量上不为零。