我一直在运行与调试模式下,英特尔的编译器版本13.1.3.192一个巨大的Fortran代码(与-O0 -g -traceback -fpe3
旗导通)。 它给了我下面的输出信息:
... ...
forrtl: warning (402): fort: (1): In call to MPI_ALLGATHER, an array temporary was created for argument #1
forrtl: error (65): floating invalid
Image PC Routine Line Source
arts 00000000016521D9 pentadiagonal_ser 385 pentadiagonal.f90
arts 0000000001644862 pentadiagonal_ 62 pentadiagonal.f90
arts 00000000004DF167 implicit_solve_x_ 1201 implicit.f90
arts 0000000000538079 scalar_bquick_inv 383 scalar_bquick.f90
arts 00000000004EFEAC scalar_step_ 190 scalar.f90
arts 0000000000401744 simulation_run_ 182 simulation.f90
arts 0000000000401271 MAIN__ 10 main.f90
arts 0000000000400FAC Unknown Unknown Unknown
arts 000000000420E444 Unknown Unknown Unknown
arts 0000000000400E81 Unknown Unknown Unknown
和错误的来源是从子程序pentadiagonal_serial,这是解决一个五对角矩阵:
subroutine pentadiagonal_serial(A,B,C,D,E,R,n,lot)
use precision
implicit none
integer, intent(in) :: n,lot
real(WP), dimension(lot,n) :: A ! LOWER-2
real(WP), dimension(lot,n) :: B ! LOWER-1
real(WP), dimension(lot,n) :: C ! DIAGONAL
real(WP), dimension(lot,n) :: D ! UPPER+1
real(WP), dimension(lot,n) :: E ! UPPER+2
real(WP), dimension(lot,n) :: R ! RHS - RESULT
real(WP), dimension(lot) :: const
integer :: i
if (n .eq. 1) then
! Solve 1x1 system
R(:,1) = R(:,1)/C(:,1)
return
else if (n .eq. 2) then
! Solve 2x2 system
const(:) = B(:,2)/C(:,1)
C(:,2) = C(:,2) - D(:,1)*const(:)
R(:,2) = R(:,2) - R(:,1)*const(:)
R(:,2) = R(:,2)/C(:,2)
R(:,1) = (R(:,1) - D(:,1)*R(:,2))/C(:,1)
return
end if
! Forward elimination
do i=1,n-2
! Eliminate A(2,i+1)
const(:) = B(:,i+1)/(C(:,i)+tiny(1.0_WP))
C(:,i+1) = C(:,i+1) - D(:,i)*const(:)
D(:,i+1) = D(:,i+1) - E(:,i)*const(:)
R(:,i+1) = R(:,i+1) - R(:,i)*const(:)
其中线
const(:) = B(:,i+1)/(C(:,i)+tiny(1.0_WP))
导致错误。 我试着打印出的值const(:)
,并发现有确实存在Infinity
价值。 不过,我不明白为什么它可以产生无穷大。 据我所见,避免零为分母,将tiny(1.0_WP)
被添加到C(:,i)
,现在它几乎是不可能的分母为零......我也查了,当这个子程序被调用,一切都被初始化或给予声明之后的值。 所以我想不通的地方可能出问题。