北航数值分析大作业第二题.docxVIP

数值分析第二次大作业 史立峰 SY1505327  方案 （1）利用循环结构将(i,j=1,2,……,10)进行赋值，得到需要变换的矩阵A； （2）然后，对矩阵A利用Householder矩阵进行相似变换，把A化为上三角矩阵A(n-1)。 对A拟上三角化，得到拟上三角矩阵A(n-1)，具体算法如下： 记A(1)=A，并记A(r)的第r列至第n列的元素为。 对于执行 若全为零，则令A(r+1) =A(r),转5；否则转2。 计算 令。 计算 继续。 （3）使用带双步位移的QR方法计算矩阵A(n-1)的全部特征值，也是A的全部特征值，具体算法如下： 给定精度水平和迭代最大次数。 记，令。 如果，则得到的一个特征值，置（降阶），转4；否则转5。 如果，则得到的一个特征值，转11；如果，则转3。 求2阶子阵 的两个特征值和，即计算二次方程 的两个根和。 如果，则得到的两个特征值和，转11；否则转7。 如果，则得到的两个特征值和，置（降阶），转4；否则转8 如果，则计算终止，未得到的全部特征值；否则转9。 记，计算 置，转3。 的全部特征值已计算完毕，停止计算。 其中，的分解与的计算用下列算法实现： 记。 对于执行 若全为零，则令,转5；否则转2。 计算 令。 计算 继续。 此算法执行完后，就得到。 （4）计算Q,R，一边求R*Q矩阵。 记 对于执行 1.若全为零，则令转5；否则转 2. 2.计算 3.令 。 4.计算 5.继续 当此算法执行完毕后，就得到正交矩阵和上三角矩阵 6.然后计算出矩阵 （5）用列主元素Gauss消去法计算矩阵对应于实特征值的特征向量，具体算法如下： 记 消元过程 对于执行 选行号，使。 交换与所含的数值。 对于计算 回代过程 最终得到的向量的即为对应于实特征值的特征向量。  二、源程序 #includeiostream.h #includemath.h #includeiomanip.h #define n 10 #define E 1.0e-12 void Householder(double a[n][n]);//拟上三角化函数 double sgn(double a);//符号函数 void QRfenjie(double a[n][n],double Q[n][n],double R[n][n]); void QR(double a[n][n],double L[n][2]);//带双步位移的QR分解 void MxM(double M[n][n],double A[n][n],double B[n][n],int m);//矩阵相乘 void solve(double a[n][n],double s1[2],double s2[2],int m);//解方程函数 void Gauss(double a[n][n], double x[n]);//定义列主元高斯消去法函数 void main() { int i,j,k; double a[n][n],qr[n][n],q[n][n],r[n][n],L[n][2],x[n]; for(i=0;in;i++)//录入要进行变换的矩a，并对q初始赋值。 { for(j=0;jn;j++) { if(i==j) a[i][j]=1.5*cos(i+1+1.2*(j+1)); else a[i][j]=sin(0.5*(i+1)+0.2*(j+1)); } } // coutendl; Householder(a);//调用拟上三角化函数 coutendlendl对矩阵A拟上三角化前七列的结果:endl; for(i=0;in;i++) { //k=0;//为了显示方便，每行显示三个元素，使用k来实现 //cout矩阵A的第i+1行元素为：endl; for(j=0;j7;j++) { if (fabs(a[i][j])E) a[i][j]=0; coutsetprecision(12)setiosflags(ios::scientific)a[i][j] ; // k++; //coutendl; //if(k%3==0)//判断每行是否到达三个元素 } coutendl; } cout对矩阵A拟上三角化后三列的结果:endl; for(i=0;in;i++) { //k=0;//为了显示方便，每行显示三个元素，使用k来实现 //cout矩阵A的第i+1行元素为：endl; for(j=7;jn;j++) { if (fabs(a[i][j])E) a[i][j]=0; c