浮点数在 Fortran 数组构造函数中造成问题
Floats causing trouble in Fortran array constructor
我遇到了一些与 Fortran 数组构造函数相关的错误。也许你能帮我理解错误？我有两个问题（见下文）。
我的最小工作示例是：
program debug
use ieee_arithmetic
implicit none
integer :: ii
real, parameter :: KZERO = 1
real, parameter :: LAMBDA = 1.3
integer, parameter :: pad = 5
integer, parameter :: nsh = 60 + 2*pad
integer, parameter :: top=nsh-pad, bot=pad+1
real(kind=8), parameter :: dt2 = 1.0D-5/2
real(kind=8), parameter :: VK_SS = 6.0D-10
real(kind=8), parameter :: VM_SS = 6.0D-05
real(kind=8), parameter,dimension(nsh) :: k = [(KZERO*LAMBDA**(ii-pad), ii=1, nsh)]
real(kind=8), parameter,dimension(nsh) :: NUK_K = [(-dt2*VK_SS*k(ii)**2, ii=1, nsh)]
real(kind=8), parameter,dimension(nsh) :: NUK_M = [(-dt2*VM_SS*k(ii)**2, ii=1, nsh)]
real(kind=8), parameter, dimension(nsh) :: A = [(.0,ii=1,pad), ( REAL(exp(NUK_K(ii))), ii=bot, top), (.0,ii=1,pad)]
real(kind=8), parameter, dimension(nsh) :: B = [(.0,ii=1,pad), ( REAL( NUK_M(ii) ), ii=bot, top), (.0,ii=1,pad)]
!real(kind=8), parameter, dimension(nsh) :: C = [(.0,ii=1,pad), ( REAL(exp(NUK_M(ii))), ii=bot, top), (.0,ii=1,pad)]
print *,A
print *,'-------------------------'
print *,B
print *,'-------------------------'
!print *,C
end program debug
第一个问题： 为什么我需要在定义了 A、B 和 C 的构造函数中使用 REAL() 进行类型转换？如果我不这样做，我会得到这个错误：
/opt/intel/composer_xe_2011_sp1.9.293/bin/intel64/ifort -o debug debug.f90
debug.f90(23): error #7113: Each ac-value expression in an array-constructor must have the same type and type parameters. [EXP]
real(kind=8), parameter, dimension(nsh) :: A = [(.0,ii=1,pad), ( exp(NUK_K(ii)), ii=bot, top), (.0,ii=1,pad)]
--------------------------------------------------------------------^
...但是 do-loop 元素已经是 REAL kind=8 因为 NUK_K 是?
第二个问题： 当我取消注释定义 C 的行（也是通过数组构造函数）时，出现以下错误：
/opt/intel/composer_xe_2011_sp1.9.293/bin/intel64/ifort -o debug debug.f90
debug.f90(25): warning #7919: The value was too small when converting to REAL(KIND=4); the result is zero.
real(kind=8), parameter, dimension(nsh) :: C = [(.0,ii=1,pad), ( REAL(exp(NUK_M(ii))), ii=bot, top), (.0,ii=1,pad)]
--------------------------------------------------------------------^
debug.f90(25): error #7768: This operation on this data type is currently inaccurate.
real(kind=8), parameter, dimension(nsh) :: C = [(.0,ii=1,pad), ( REAL(exp(NUK_M(ii))), ii=bot, top), (.0,ii=1,pad)]
--------------------------------------------------------------------^
..警告没问题，因为我猜这意味着编译器建议 exp(-[large number]) -> 0，对吗？
但是我该如何处理错误 "This operation on this data type is currently inaccurate." ?
我已经解决这个问题很长时间了，希望你能帮助我！
编辑
首先非常感谢您的热心解答！我很感激。然而，不准确的问题仍然存在。有任何想法吗？我最初的猜测是使用四边形而不是双精度，然后将结果转换为双精度（我的 CPU 本身不支持四边形）。这也不起作用。我的新最小工作示例是（注意我删除了零填充以使示例更简单）：
program debug
use ieee_arithmetic
use ISO_Fortran_env
implicit none
integer, parameter :: rd = real64
integer, parameter :: rq = real128
integer, parameter :: df = rq ! Default float kind
real(df), parameter :: KZERO = 1
real(df), parameter :: LAMBDA = 1.3
integer, parameter :: pad = 5
integer, parameter :: nsh = 60 + 2*pad
integer, parameter :: top=nsh-pad, bot=pad+1
real(df), parameter :: dt2 = 1.0D-5/2
real(df), parameter :: VK_SS = 6.0D-10
real(df), parameter :: VM_SS = 6.0D-05
integer :: ii
real(df), parameter,dimension(nsh) :: k = [(KZERO*LAMBDA**(ii-pad), ii=1, nsh)]
real(df), parameter,dimension(nsh) :: NUK_K = [(-dt2*VK_SS*k(ii)**2, ii=1, nsh)]
real(df), parameter,dimension(nsh) :: NUK_M = [(-dt2*VM_SS*k(ii)**2, ii=1, nsh)]
real(df), parameter, dimension(nsh) :: A = [ ( exp(NUK_K(ii)) , ii=1, nsh) ]
real(df), parameter, dimension(nsh) :: B = [ ( NUK_M(ii) , ii=1, nsh) ]
real(df), parameter, dimension(nsh) :: C = [ ( exp(NUK_M(ii)) , ii=1, nsh) ]
print *,A
print *,'-------------------------'
print *,B
print *,'-------------------------'
print *,C
end program debug
但这仍然报错
/opt/intel/composer_xe_2011_sp1.9.293/bin/intel64/ifort -o debug debug.f90
debug.f90(30): error #7768: This operation on this data type is currently inaccurate.
real(df), parameter, dimension(nsh) :: C = [ ( exp(NUK_M(ii)) , ii=1, nsh) ]
--------------------------------------------------^
第一个问题：
您正在尝试将默认类型 real(4) 的数组 (.0,ii=1,pad) 与 ( REAL(exp(NUK_K(ii))), ii=bot, top) 数组 real(8) 连接起来。最简单的解决方案是将第一个数组也变成 real(8) 数组，例如将浮点数 .0 替换为 0.0d0（双精度数的科学计数法）。
示例：
real(kind=8), parameter, dimension(nsh) :: A = [(0.0d0,ii=1,pad), ( exp(NUK_K(ii)), ii=bot, top), (0.0d0,ii=1,pad)]
real(kind=8), parameter, dimension(nsh) :: B = [(0.0d0,ii=1,pad), ( NUK_M(ii) , ii=bot, top), (0.0d0,ii=1,pad)]
real(kind=8), parameter, dimension(nsh) :: C = [(0.0d0,ii=1,pad), ( exp(NUK_M(ii)), ii=bot, top), (0.0d0,ii=1,pad)]
第二个问题：
默认情况下，real(num) 不会将数字 num 类型转换为 real(8)，而是类型转换为默认类型 real(4)。如果您想将其类型转换为 real(8)，您必须使用符号 real(num,kind=8) 明确地这样做。
所以换句话说，您一直在做的是显式地将 ( exp(NUK_M(ii)), ii=bot, top) 从双精度类型转换为单精度浮点数，然后将其与单精度零数组连接起来，然后将结果再次为双精度数组。这就是编译器抱怨精度损失的原因。
可能感兴趣的人的跟进：
我决定做一个 hacky 解决方案来解决有关内部函数 exp() 精度的问题。
问题是我希望计算 exp(-[非常大的数字])，它变得太接近于零，以至于偶数四边形浮点数无法表示大（负）指数。问题是，我不能简单地将这些条目设置为零，因为这会扰乱我的完整问题的数字。
因此，我所做的是用 exp(-7.0D2) 替换 exp(-[number lager than 7.0D2])，这接近于使用双精度表示的极限。
我遇到了一些与 Fortran 数组构造函数相关的错误。也许你能帮我理解错误？我有两个问题（见下文）。
我的最小工作示例是：
program debug
use ieee_arithmetic
implicit none
integer :: ii
real, parameter :: KZERO = 1
real, parameter :: LAMBDA = 1.3
integer, parameter :: pad = 5
integer, parameter :: nsh = 60 + 2*pad
integer, parameter :: top=nsh-pad, bot=pad+1
real(kind=8), parameter :: dt2 = 1.0D-5/2
real(kind=8), parameter :: VK_SS = 6.0D-10
real(kind=8), parameter :: VM_SS = 6.0D-05
real(kind=8), parameter,dimension(nsh) :: k = [(KZERO*LAMBDA**(ii-pad), ii=1, nsh)]
real(kind=8), parameter,dimension(nsh) :: NUK_K = [(-dt2*VK_SS*k(ii)**2, ii=1, nsh)]
real(kind=8), parameter,dimension(nsh) :: NUK_M = [(-dt2*VM_SS*k(ii)**2, ii=1, nsh)]
real(kind=8), parameter, dimension(nsh) :: A = [(.0,ii=1,pad), ( REAL(exp(NUK_K(ii))), ii=bot, top), (.0,ii=1,pad)]
real(kind=8), parameter, dimension(nsh) :: B = [(.0,ii=1,pad), ( REAL( NUK_M(ii) ), ii=bot, top), (.0,ii=1,pad)]
!real(kind=8), parameter, dimension(nsh) :: C = [(.0,ii=1,pad), ( REAL(exp(NUK_M(ii))), ii=bot, top), (.0,ii=1,pad)]
print *,A
print *,'-------------------------'
print *,B
print *,'-------------------------'
!print *,C
end program debug
第一个问题： 为什么我需要在定义了 A、B 和 C 的构造函数中使用 REAL() 进行类型转换？如果我不这样做，我会得到这个错误：
/opt/intel/composer_xe_2011_sp1.9.293/bin/intel64/ifort -o debug debug.f90
debug.f90(23): error #7113: Each ac-value expression in an array-constructor must have the same type and type parameters. [EXP]
real(kind=8), parameter, dimension(nsh) :: A = [(.0,ii=1,pad), ( exp(NUK_K(ii)), ii=bot, top), (.0,ii=1,pad)]
--------------------------------------------------------------------^
...但是 do-loop 元素已经是 REAL kind=8 因为 NUK_K 是?
第二个问题： 当我取消注释定义 C 的行（也是通过数组构造函数）时，出现以下错误：
/opt/intel/composer_xe_2011_sp1.9.293/bin/intel64/ifort -o debug debug.f90
debug.f90(25): warning #7919: The value was too small when converting to REAL(KIND=4); the result is zero.
real(kind=8), parameter, dimension(nsh) :: C = [(.0,ii=1,pad), ( REAL(exp(NUK_M(ii))), ii=bot, top), (.0,ii=1,pad)]
--------------------------------------------------------------------^
debug.f90(25): error #7768: This operation on this data type is currently inaccurate.
real(kind=8), parameter, dimension(nsh) :: C = [(.0,ii=1,pad), ( REAL(exp(NUK_M(ii))), ii=bot, top), (.0,ii=1,pad)]
--------------------------------------------------------------------^
..警告没问题，因为我猜这意味着编译器建议 exp(-[large number]) -> 0，对吗？ 但是我该如何处理错误 "This operation on this data type is currently inaccurate." ?
我已经解决这个问题很长时间了，希望你能帮助我！
编辑
首先非常感谢您的热心解答！我很感激。然而，不准确的问题仍然存在。有任何想法吗？我最初的猜测是使用四边形而不是双精度，然后将结果转换为双精度（我的 CPU 本身不支持四边形）。这也不起作用。我的新最小工作示例是（注意我删除了零填充以使示例更简单）：
program debug
use ieee_arithmetic
use ISO_Fortran_env
implicit none
integer, parameter :: rd = real64
integer, parameter :: rq = real128
integer, parameter :: df = rq ! Default float kind
real(df), parameter :: KZERO = 1
real(df), parameter :: LAMBDA = 1.3
integer, parameter :: pad = 5
integer, parameter :: nsh = 60 + 2*pad
integer, parameter :: top=nsh-pad, bot=pad+1
real(df), parameter :: dt2 = 1.0D-5/2
real(df), parameter :: VK_SS = 6.0D-10
real(df), parameter :: VM_SS = 6.0D-05
integer :: ii
real(df), parameter,dimension(nsh) :: k = [(KZERO*LAMBDA**(ii-pad), ii=1, nsh)]
real(df), parameter,dimension(nsh) :: NUK_K = [(-dt2*VK_SS*k(ii)**2, ii=1, nsh)]
real(df), parameter,dimension(nsh) :: NUK_M = [(-dt2*VM_SS*k(ii)**2, ii=1, nsh)]
real(df), parameter, dimension(nsh) :: A = [ ( exp(NUK_K(ii)) , ii=1, nsh) ]
real(df), parameter, dimension(nsh) :: B = [ ( NUK_M(ii) , ii=1, nsh) ]
real(df), parameter, dimension(nsh) :: C = [ ( exp(NUK_M(ii)) , ii=1, nsh) ]
print *,A
print *,'-------------------------'
print *,B
print *,'-------------------------'
print *,C
end program debug
但这仍然报错
/opt/intel/composer_xe_2011_sp1.9.293/bin/intel64/ifort -o debug debug.f90
debug.f90(30): error #7768: This operation on this data type is currently inaccurate.
real(df), parameter, dimension(nsh) :: C = [ ( exp(NUK_M(ii)) , ii=1, nsh) ]
--------------------------------------------------^
第一个问题：
您正在尝试将默认类型 real(4) 的数组 (.0,ii=1,pad) 与 ( REAL(exp(NUK_K(ii))), ii=bot, top) 数组 real(8) 连接起来。最简单的解决方案是将第一个数组也变成 real(8) 数组，例如将浮点数 .0 替换为 0.0d0（双精度数的科学计数法）。
示例：
real(kind=8), parameter, dimension(nsh) :: A = [(0.0d0,ii=1,pad), ( exp(NUK_K(ii)), ii=bot, top), (0.0d0,ii=1,pad)]
real(kind=8), parameter, dimension(nsh) :: B = [(0.0d0,ii=1,pad), ( NUK_M(ii) , ii=bot, top), (0.0d0,ii=1,pad)]
real(kind=8), parameter, dimension(nsh) :: C = [(0.0d0,ii=1,pad), ( exp(NUK_M(ii)), ii=bot, top), (0.0d0,ii=1,pad)]
第二个问题：
默认情况下，real(num) 不会将数字 num 类型转换为 real(8)，而是类型转换为默认类型 real(4)。如果您想将其类型转换为 real(8)，您必须使用符号 real(num,kind=8) 明确地这样做。
所以换句话说，您一直在做的是显式地将 ( exp(NUK_M(ii)), ii=bot, top) 从双精度类型转换为单精度浮点数，然后将其与单精度零数组连接起来，然后将结果再次为双精度数组。这就是编译器抱怨精度损失的原因。
可能感兴趣的人的跟进：
我决定做一个 hacky 解决方案来解决有关内部函数 exp() 精度的问题。 问题是我希望计算 exp(-[非常大的数字])，它变得太接近于零，以至于偶数四边形浮点数无法表示大（负）指数。问题是，我不能简单地将这些条目设置为零，因为这会扰乱我的完整问题的数字。 因此，我所做的是用 exp(-7.0D2) 替换 exp(-[number lager than 7.0D2])，这接近于使用双精度表示的极限。