北航数值分析大作业第一题(fortran).docxVIP

“数值分析”计算实习报告——SY1415215孔维鹏一、算法说明1、采用幂法计算实对称矩阵A的按模最大特征值。由于是矩阵A的所有特征值当中模最大的一个，因此是或者。2、采用幂法计算矩阵A平移后矩阵的按模最大特征值。假设，那么所求得的平移后矩阵B的按模最大特征值为平移后的特征值；反之，，那么所求得的按模最大特征值为平移后的特征值。比较与，两者大的为，小的为。3、采用反幂法计算矩阵A的按模最小特征值。在反幂法中计算时，采用Doolittle分解法。4、非奇异的实对称矩阵A的条件数。5、采用Doolittle分解法将实对称矩阵A分解为一个上三角矩阵U，和单位下三角矩阵L，即。6、在存储矩阵A的元素时，设置一个数组C(m,n)来存放A的带内元素，其中，。数组C的第j列存放A的第j列元素，并使A得主对角元素存于Ｃ的第行中。在程序计算时，使用来代替。7、本程序中幂法采用范数，初始向量。二、源程序（FORTRAN）PROGRAM SY1415215DIMENSION XL(40),B(5,501),C(5,501),U(5,501)PARAMETER (EI=2.718281828459)REAL(8) X1,X2,X3,X,XM,X0!X为按模最大特征值，XM为带原点平移后的按模最小特征值REAL(8) Y,XL,B,C,UREAL(8) DETINTEGER RN=501R=2OPEN(1,FILE=计算实习大作业第一题计算结果.TXT)WRITE(1,(计算实习大作业第一题计算结果:,/))!存储矩阵A!存储矩阵A的所有非零元素至数组CDO I=1,5 IF((I==1).OR.(I==5)) THEN DO J=3,N C(I,J)=-0.064 ENDDO ELSE IF((I==2).OR.(I==4)) THEN DO J=2,N C(I,J)=0.16 ENDDO ELSE IF(I==3) THEN DO J=1,N C(I,J)=(1.64-0.024*J)*SIN(0.2*J)-0.64*EI**(0.1/J) ENDDO ENDIFENDDO!采用幂法计算A的按模最大特征值与特征向量CALL POVERMETHOD(C,N,X0)!采用带原点平移的幂法求矩阵A的最大与最小特征值DO I=1,5 DO J=1,N IF(I==3) THEN B(I,J)=C(I,J)-X0 ELSE B(I,J)=C(I,J) ENDIF ENDDOENDDO !求出平移后的矩阵CALL POVERMETHOD(B,N,XM)!运用幂法求出平移后矩阵的按模最大特征值X1=X0X2=XM+X0IF (X1X2) THENELSEX=X1X1=X2X2=XENDIFWRITE(1,(矩阵A的最小特征值为：,E20.13)) X1WRITE(1,(矩阵A的最大特征值为：,E20.13)) X2!运用反幂法求解矩阵A的按模最小特征值CALL INPOVERMETHOD(C,N,X3)WRITE(1,(矩阵A的按模最小特征值为：,E20.13)) X3WRITE(1,( ))WRITE(1,(与指定数最接近的特征值为：))WRITE(1,(3X,K,16X,μ,18X,λ))!求与指定数D最接近的特征值，特征值存储到XL中DO K=1,39 D=X1+K*(X2-X1)/40 !采用带原点平移的反幂法求与指定数D最接近的特征值 DO I=1,5 DO J=1,N IF(I==3) THEN B(I,J)=C(I,J)-D ELSE B(I,J)=C(I,J) ENDIF ENDDO ENDDO !求出平移后的矩阵 CALL INPOVERMETHOD(B,N,XL(K)) WRITE(1,(1X,I3,2X,E20.13,2X,E20.13)) K,D,XL(K)+D ENDDO!计算A的（谱范数）条件数cond（A）和行列式detACOND=ABS(X0/X3)WRITE(1,(A的（谱范数）条件数cond（A）为：,E20.13)) COND!doolittle分解得到的上三角阵乘积为detACALL DETA(C,N,U)DET=1DO I=1,N DET=DET*U(3,I) WRITE(1,(E20.13)) U(3,I)ENDDOWRITE(1,(A的行列式的值为：,E20.13)) DETCLOSE(1)END!幂法子程序SUBROUTINE POVERMETHOD(C1,N,X)DIMENSION U(N),Y(N),C(5,N),C1(5,501)REAL(8) U,Y,C,C1REAL(8) X,Z1,Z2,S1,S2,S3,B,EC=C1DO I=1,N U(I)=1ENDDO !任取非零向量UZ2=EIZ1=1.0