数值分析报告1.doc

数 值 分 析 计 算 报 告第 一 次 作 业系 （院）： 宇航学院 专 业： 航空宇航推进理论与工程 学生姓名： 向纪鑫 学 号： SY1415313 2014年11月北京航空航天大学大学数值分析大作业一计算方案的设计1.初始化并压缩存储矩阵A定义矩阵，由矩阵性质可知。根据计算实习说明书要求，所有矩阵的零元素都不存储，所以将矩阵压缩存储为矩阵。在矩阵c中检索矩阵的带内元素的方法为：。2.求和的值首先，用幂法求出矩阵按模最大的特征值。然后，根据原点平移法，同样根据幂法求出矩阵按模最大的特征值。当时，；当时，；当时，或，若，则，若，则。最后，综上所述，只要求出和，通过比较就可以求出和。3. 求,的值用反幂法就可以直接求出。结合上一步求出的，从而很容易求出。其中，用反幂法求的过程中用到了Doolittle分解，即把矩阵写成,是上三角矩阵，所以的值就是所有矩阵主对角线元素的乘积。4. 求的与最接近的特征值用带平移的反幂法求出的值然后加上就是，即：求出矩阵按模最小的特征值，。2 计算源程序File1 初始化的函数（shushihua）#include#includevoid chushihua(double A5501) int i;for(i=0;i501;i+) A2i=(1.64-0.024*(i+1)*sin(0.2*(i+1)-0.64*exp(0.1/(i+1); A0i=A4i=-0.064; A1i=A3i=0.16; A00=A01=A10=A4500=A4499=A3500=0;File 2 幂法的函数 （mifa）#include#includedouble mifa(double A5501) int i,j,t,k=0; double g,B=1,c=2,u501,y501; for(i=0;i501;i+) ui=1; do g=0; c=B; for(i=0;i501;i+) g=ui*ui+g; g=sqrt(g); for(i=0;i501;i+) yi=ui/g; for(i=0;i501;i+) ui=0; for(j=0;j=0&t=4) ui=Atj*yj+ui; B=0; for(i=0;i1e-12); return(B);File 3 Doolittle分解 （doolittle）#include#includevoid doolittle(double A5501)int k,i,j,t,ik;double S501,Max,B=1;for(k=1;k=501;k+) /*用选主元Doolittle分解法解方程组（LU分解）*/ ik=k; Sk-1=0; for(t=1;t=k-1;t+) if(fabs(k-t)3) Sk-1=Ak-t+2t-1*At-k+2k-1+Sk-1; Sk-1=A2k-1-Sk-1; Max=Sk-1;for(i=k;i=501;i+) Si-1=0; for(t=1;t=k-1;t+) if(fabs(i-t)3&fabs(t-k)3) Si-1=Ai-t+2t-1*At-k+2k-1+Si-1; if(fabs(i-k)3) Si-1=Ai-k+2k-1-Si-1; else Si-1=0-Si-1; A2k-1=Max;for(j=k+1;j=501;j+) /*求Ukj*/B=0;for(t=1;t=k-1;t+)if(fabs(k-t)3&fabs(t-j)3)B=Ak-t+2t-1*At-j+2j-1+B; if(fabs(k-j)3) Ak-j+2j-1=Ak-j+2j-1-B;for(i=k+1;i=501;i+) /*求Lik*/if(fabs(i-k)3)Ai-k+2k-1=Si-1/A2k-1;File 4 反幂法的函数 （fanmifa）#include#includedouble fanmifa(double A5501)extern void doolittle(double A5501);int k,i,j,t,ik,K=0;double S501,x501,u501,y501,H,D,Max,g,B=1,c=2,det;doolittle(A);for(i=0;i501;i+)ui=1;do /*求按模最小特征值*/ g=0; c=B;for(i=0;i501;i+)g=ui*ui+g; g=sqrt(g);for(i=0;i501;i+)yi=ui/g;x0=y0; /*求解x、y*/for(i=2;i=501;i+) B=0; for(t=1;t=i-1;t+)if(fabs(i-t)=1;i-) B=0; for(t=i+1;t=501;t+) if(fabs(i-t)3) B=B+Ai-t+2t-1*ut-1; ui-1=(xi-1-B)/A2i-1; B=0;for(i=0;i1e-12); B=1/B;return(B);最终主程序： #include#includevoid blank(void);/定义函数输出空格。void star1(void);/定义函数输出星号。void star2(void);/定义函数输出星号。void blank(void)/函数代码printf( );void star1(void)/函数代码。blank();for(int ii=0;ii=61;ii+)printf(*);printf(n);void star2(void)/函数代码。for(int ii=0;ii=17;ii+)printf(*);void main()extern void chushihua(double A5501); extern double mifa(double A5501); extern double fanmifa(double A5501);int i,j,k=0;double A5501,u501,R39,R_s,cond,B1,B2,det=1,uk39;chushihua(A);B1=mifa(A);for(i=0;i501;i+)A2i=A2i-B1;B2=mifa(A);B2=B2+B1;chushihua(A); /求条件数condR_s=fanmifa(A);cond=fabs(B1/R_s);for(i=0;i501;i+) /求行列式的值detdet=A2i*det;for(i=0;i39;i+)uki=B1+(i+1)*(B2-B1)/40;for(i=0;i39;i+) chushihua(A); for(j=0;jB1) /输出结果printf(1=%.11en501=%.11en,B1,B2);elseprintf(1=%.11en501=%.11en,B2,B1);printf(s=%.11en,R_s); printf(cond=%.11en,cond);printf(det=%.11en,det);for(i=0;i39;i+)printf(R%d=%.11en,i+1,Ri);3.计算结果4. 迭代初始向量的选取对迭代结果的影响在初始向量为u0=1 . . 1 1的情况下，幂法求得= -1.07001136150e+001 迭代343次；反幂法求=-5.55791079423e-003，迭代22次。当使用初始向量u0=1 . . 0 0时，=-2.08098108524e+000，幂法迭代160次；=2.66888692378e-002，反幂法迭代231次；相应的其他特征值和其迭代次数也发生了变化。因