北航数值分析计算实习第一题打印.docVIP

北航数值作业第一题：； 2、求A的最大最小特征值，可以对A用幂法求出模最大的特征值λm1，然后以λm1为平移量对A进行平移后再用幂法求出λm，令λm2=λm1+λm，可以得知λm1，λm2为A的最大最小特征值，从而可得到A的最大特征值λ501，最小特征值λ1； 3、因为在用反幂法的过程中要解方程组,还要求A的行列式，所以对A进行LU分解，并将分解后的元素也存在C中，以便在求A的模最小的特征值和行列式的时候调用； 4、用反幂法求出A的模最小的特征值λs； 5、求与数μk最接近的特征值： 依次k=1到39执行：对A按μk进行平移得B=A-μkI，并将B存在C中，再对B用反幂法求出B的特征值λ，则离μk最近的特征值为λik=λ+μk； 5、A为实对称矩阵，所以A的条件数,其中λm1，λs分别为A的模最大和模最小的特征值； 6、调用A的LU分解程序，则分解后的U矩阵的对角线相乘即为矩阵的。 #includestdio.h #includemath.h #includeconio.h //定义A中元素 double C[5][501]; double a[501]; double b; double c; //声明所有函数 void YaSuoJZ(double C[5][501],double a[501],double b,double c) ;//压缩矩阵函数 double mifa(double C[5][501]); //幂法函数 void daizhuangLU(double A[5][501]); //带状矩阵的LU分解 double fanmifa(double C[5][501]);//反幂法函数 //最值函数 int max2(int x,int y); int max3(int x,int y,int z); int min(int x,int y); //最值函数 int max2(int x,int y) //求2个数的最大值 { int z; z=xy?x:y; return(z); } int max3(int x,int y,int z) //求3个数的最大值 { int w; w = z max2(x,y)? z:max2(x,y); return(w); } int min(int x,int y) //求2个数的最小值 { int z; z=xy?y:x; return(z); } //将矩阵A压缩存储在矩阵C中 void YaSuoJZ(double C[5][501],double a[501],double b,double c) { int i; for(i=0;i=500;i++) { if(i=2) C[0][i]=c; else C[0][i]=0; if(i=1) C[1][i]=b; else C[1][i]=0; if(i=499) C[3][i]=b; else C[3][i]=0; if(i=498) C[4][i]=c; else C[4][i]=0; C[2][i]=a[i]; } } //幂法函数：用幂法求矩阵模最大的特征值 double mifa(double C[5][501]) { double u[501]; double y[501]={0},η=0; double β,βk=0; double ε=1; // ε为精度 double sumu=0,sumAY=0; int i,j,k=1; //k为循环次数 for (i=0;i=500;i++) //取任一非零向量u0 u[i] = 1.0; while(ε=1e-12) { for(i=0;i=500;i++) //求u(k-1)的2范数η sumu=sumu+u[i]*u[i]; η=sqrt(sumu); sumu=0; for(i=0;i=500;i++) //求y(k-1) y[i]=u[i]/η; for(i=0;i=500;i++) //求u(k)的各分量u[i] { for(j=max2(0,i-2);j=min(i+2,500);j++) sumAY=sumAY+C[i-j+2][j]*y[j]; u[i]=sumAY; sumAY=0; } //求幂法中的βk