浮点溢出

【问题标题】：Floating point overflow浮点溢出
【发布时间】：2016-12-05 22:07:23
【问题描述】：

我在这部分代码中遇到浮点溢出错误。哪位大神帮我看看原因。

do j=1,ny-1
  do i=1,nx-1

    sum = 0.0d0

    do k=0,1000
      n=2.0d0*dfloat(k)+ 1.0d0
      sum = sum + ((dsinh(n*pi*x(i))*dcos(n*pi*y(j)))/((n*n*pi*pi)*dsinh(2*n*pi)))
    end do

    ue(i,j)= (x(i)/(4.0d0))- 4.0d0*sum

  end do
end do

【问题讨论】：

一些 IMPLICIT NONE 和实际程序可能会有所帮助。我将从 sum = sum + 右边的东西的中间产品开始。然后，如果/当 k 变大时，您可能希望 sum 是 double (8) 的乘积是 float (*4)。由于 SUM 是一个内在函数，我会将“sum”重命名为 MySum，或者有一个数组维度（0:1000）并使用 SUM（MySum）。为什么要从 0 开始数组？ 1001分似乎很奇怪。 (npipi) 和 (2*npi) 和 (n*pi) 可以是预先计算的 A、B 和 C 或参数，然后它更快，看起来更干净.
欢迎来到 Stack Overflow。你应该展示一个完整的程序并描述你是如何编译它的，并展示你得到的实际输出。不要忘记输入数据。请务必阅读帮助页面stackoverflow.com/tour stackoverflow.com/help/how-to-ask
请努力正确格式化您的帖子和代码。我对您的代码应用了一些基本的缩进。

标签： fortran fortran90 fortran95

【解决方案1】：

问题是中间词dsinh(2*n*pi)。考虑k=1000。然后n=2001 所以我们需要评估dsinh(2001*pi) 大约是0.5*exp(6286) 或超过10^2700！这远远高于任何可以用双精度表示的数字。您需要重新评估计算总和的方式。 dsinh(n*pi*x(i)) 这个词也有问题。

我的猜测是，对于商 dsinh(n*pi*x(i))/dsinh(2*n*pi) 的稳健评估需要某种渐近扩展。对于0<x(i)<2，这个术语应该表现为exp(n*pi*(x(i)-2))，因为n 变大了。这将是很好的表现。

【讨论】：