【发布时间】:2013-12-07 22:12:12
【问题描述】:
我有以下问题:
我在Fortran 77 上有一个子例程,它必须返回一个一维数组 D_Y。所有需要的值都从主程序正确传输,但是矩阵 D_Y (D_Y(16-18)) 的三个元素设置为 0.0,我不知道为什么。
NDIM=18,P1和P2是对的,它们在方程中的所有参数也是对的,但是它们的总和总是设置为0!
SUBROUTINE DY(Y,T,PARAM,NDIM,D_Y)
IMPLICIT NONE
real(8),intent(in)::Y(*)
real(8),intent(in)::PARAM(*),T
real(8),intent(out)::D_Y(NDIM)
integer::D
integer,intent(in)::NDIM
real(8)::R,R1,R2,P1,P2
open(18,file='kinetic2.txt',status='unknown')
R=((Y(1)-Y(4))**2.D0+(Y(2)-Y(5))**2.D0+(Y(3)-Y(6))**2.D0)**0.5D0
R1=((Y(1)-Y(13))**2.D0+(Y(2)-Y(14))**2.D0+(Y(3)-Y(15))**2.D0)**0.5D0
R2=((Y(4)-Y(13))**2.D0+(Y(5)-Y(14))**2.D0+(Y(6)-Y(15))**2.D0)**0.5D0
DO D=1,6
D_Y(D)=Y(D+6)
END DO
DO D=7,9
D_Y(D)=-(PARAM(1)*PARAM(3))*(Y(D-6)-Y(D-3))/((R)**3.D0)
END DO
DO D=10,12
D_Y(D)=-(PARAM(1)*PARAM(2))*(Y(D-6)-Y(D-9))/((R)**3.D0)
END DO
DO D=13,15
D_Y(D)=Y(D+3)
END DO
DO D=16,18
P1=-(PARAM(1)*PARAM(3))*(Y(D-3)-Y(D-12))/((R2)**3.D0)
P2=-(PARAM(1)*PARAM(2))*(Y(D-3)-Y(D-15))/((R1)**3.D0)
D_Y(D)=P1+P2
write(18,*) P1,P2,D_Y(D)
END DO
RETURN
END
【问题讨论】:
-
那行真的是这样缩进的吗?
-
次要 - 你没有关闭通道 18。为什么不在循环中打印 param, y 和 r1, r2 而不是破坏你的大脑试图找出问题所在。从输出中可能会很明显。
-
另外,与在循环中计算 R**3.0D0 不同,在循环外计算它们会更有效。没有发现您的错误,但需要考虑的事情更少。这些数字真的很大吗?
-
你能发布一个显示错误的最小工作程序吗?这将使调试代码变得更加简单;-) 顺便说一下,不要使用
**2.D0- 这对编译器来说很难优化!最好使用**2甚至x*x。 -
Y,R1和R2是对的,我已经查过了。程序正在运行 :) 我没有收到任何错误消息。无论如何,谢谢大家的建议。