FORTRAN里TOL怎么定义,为什么我的计算不准
数值分析里的一道题目,在最后判断是否收敛时,如果使用下面这个语句则结果准确IF (DABS(DABS(M1)-DABS(M0))/DABS(M1)<TOL) EXIT
而使用
IF (DABS(DABS(M1)-DABS(M0))/DABS(M1)<1D-12) EXIT
则与正确答案不符合,这里不知道TOL,是什么东西,之前只是声明了变量,并没有赋值,是不是就是系统默认的最小值,
第二个问题,1D-12也很小了,为什么会不正确?
回答问题1:声明变量没有初始化或者其他赋值操作,为任意值(取决于系统当前运行状态),不知道是多少,也没有概率
回答问题2:误差限看需求以及变量精度而定,主要取决于需求,不一定越小越好
另外问题2中的不正确指的是什么? Fortran里面请参考epsilon函数
https://gcc.gnu.org/onlinedocs/gcc-7.2.0/gfortran/EPSILON.html#EPSILON
program test_epsilon
real :: x = 3.143
real(8) :: y = 2.33
print *, EPSILON(x)
print *, EPSILON(y)
end program test_epsilon 本帖最后由 shaoqikao100 于 2017-11-13 15:48 编辑
Jackdaw 发表于 2017-11-12 23:37
回答问题1:声明变量没有初始化或者其他赋值操作,为任意值(取决于系统当前运行状态),不知道是多少,也 ... SUBROUTINE DOOLITTLE(A,L,U,N)
!Doolittle分解
IMPLICIT REAL*8(A-Z)
INTEGER::N,I,K,R
REAL*8::A(N,N),L(N,N),U(N,N)
N = 501
!U的第一行
U(1,:)=A(1,:)
!L的第一列
L(:,1)=A(:,1)/U(1,1)
DO K=2,N
L(K,K)=1
DO J=K,N
S=0
DO M=1,K-1
S=S+L(K,M)*U(M,J)
END DO
U(K,J)=A(K,J)-S
END DO
DO I=K+1,N
S=0
DO M=1,K-1
S=S+L(I,M)*U(M,K)
END DO
L(I,K)=(A(I,K)-S)/U(K,K)
END DO
END DO
END SUBROUTINE DOOLITTLE
SUBROUTINE UPTRI(A,B,X,N)
!解上三角方程组
IMPLICIT REAL*8(A-Z)
INTEGER::I,J,K,N
REAL*8::A(N,N),B(N),X(N)
N = 501
X(N)=B(N)/A(N,N)
!回代部分
DO I=N-1,1,-1
X(I)=B(I)
DO J=I+1,N
X(I)=X(I)-A(I,J)*X(J)
END DO
X(I)=X(I)/A(I,I)
END DO
END SUBROUTINE UPTRI
SUBROUTINE DOWNTRI(A,B,X,N)
!解下三角方程组
IMPLICIT REAL*8(A-Z)
INTEGER::I,J,N
REAL*8::A(N,N),B(N),X(N)
N = 501
X(1)=B(1)/A(1,1)
DO K=2,N
X(K)=B(K)
DO I=1,K-1
X(K)=X(K)-A(K,I)*X(I)
END DO
X(K)=X(K)/A(K,K)
END DO
END SUBROUTINE DOWNTRI
SUBROUTINE UNMIFA(A,NAMDA4)
IMPLICIT NONE
INTEGER::N,I,K,L,S,J,M,K1
DIMENSION A(501,501),U(501),U0(501),V(501),U1(501,501),L1(501,501),Y(501)
DOUBLE PRECISION MA,A,NAMDA4,U,U0,V,M0,M1,TOL,U1,L1,Y
N = 501
DO I=1,N
U0(N)=10D0
END DO
!给出u0
U=U0
M0=0
!------------------------------------------------
CALL DOOLITTLE(A,L1,U1,N)
!------------------------------------------------
DO K1=1,500
CALL DOWNTRI(L1,U,Y,N)
CALL UPTRI(U1,Y,V,N)
!------------------------------------------------
!用解方程组的方法求矩阵v
!求最大值
M1=DABS(V(1))
DO I=2,N
IF (DABS(V(I))>M1) THEN
M1=DABS(V(I))
L = I !用k记录指标,但是ma不取r(i)
END IF
END DO
M1 = V(L)
!找出v中的最大值,返回到m1 SUBROUTINE MAX_ROU(R,N,MA)
U = V/M1
!IF (DABS(M1-M0)<1E-20) EXIT
IF (DABS(DABS(M1)-DABS(M0))/DABS(M1)<TOL) EXIT
M0=M1
END DO
NAMDA4=M1
NAMDA4=1/NAMDA4
END
PROGRAM MAIN
IMPLICIT NONE
INTEGER I,J,N
DIMENSION A(501,501)
DOUBLE PRECISION A,NAMDA4
!给a赋值!!!!!!!!!!!!!!!!!!!!!!!!!
DO I = 1,501
DO J = 1,501
IF(I == J) THEN
A(I,J) =(1.64-0.024*I)*SIN(0.2*I)-0.64*EXP(0.1/I)
ELSE IF(I == J+1) THEN
A(I,J) = 0.16
ELSE IF(I+1 == J) THEN
A(I,J) = 0.16
ELSE IF(I == J+2) THEN
A(I,J) = -0.064
ELSE IF(I+2 == J) THEN
A(I,J) = -0.064
!ELSE
! A(I,J) = 0
END IF
END DO
END DO
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
CALL UNMIFA(A,NAMDA4)
WRITE(*,*) NAMDA4
END
您好,非常感谢大神回答,小白一枚,这是看到书上的一些代码后自己写的一个用反幂法求特征值的程序,其中在unmifa的子程序里,如果使用tol判别,结果正确,而用我自己的则结果错误。能否看一下是不是代码本身有问题啊。
如果使用()<1D-12判别结果为 -0.178858282679491
使用()<TOL 判别结果为 -5.558469826156424E-003结果是对的,并且用于其它矩阵也可以。
请看一下,这是不是巧合
pasuka 发表于 2017-11-13 10:04
Fortran里面请参考epsilon函数
https://gcc.gnu.org/onlinedocs/gcc-7.2.0/gfortran/EPSILON.html#EPSILON
...
十分感谢,我的代码在上面的回复里贴出来了,这里如果要求的精度很明确的话怎么操作。我是看的白海波老师的那本fortran,他的程序里直接使用了tol这个参数,我不太懂 shaoqikao100 发表于 2017-11-13 15:41
十分感谢,我的代码在上面的回复里贴出来了,这里如果要求的精度很明确的话怎么操作。我是看的白海波老师 ...
IMPLICIT REAL*8(A-Z)
往后见着上面这行代码,请直接无视余下所有代码,莫要浪费青春
希望明天就能交作业,请出门左转去马云家
希望三个礼拜达到完成作业的Fortran水平,先把本站的教学视频看三遍,彭国伦的Fortran95程序设计看两遍
http://v.fcode.cn/
http://fcode.cn/resource_ebook-1-1.html pasuka 发表于 2017-11-13 16:15
往后见着上面这行代码,请直接无视余下所有代码,莫要浪费青春
希望明天就能交作业,请出门左转去马云家 ...
好的,正在学习。等看完看看能不能解决这个问题
页:
[1]