浮点溢出

我得到这部分代码中的浮点溢出错误。你们能帮助我找出原因吗？浮点溢出

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

来源

2016-12-05 Slow learner

一些IMPLICIT NONE和一个实际的程序可能会有所帮助。我向dstart和sum = sum +右边的东西的中间产物跑去。然后，如果/当k变大时，您可能希望将sum作为产品的双倍（* 8）是浮点数（* 4）。由于SUM是一个内部函数，我会将'sum'重命名为MySum，或者使用一个数组维（0：1000）并使用SUM（MySum）。为什么你要在0开始数组？ 1001分似乎很奇怪。（n * pi * pi）和（2 * n * pi）和（n * pi）可以预先计算A，B和C或参数，然后它更快，看起来更干净。 – Holmz

欢迎来到Stack Overflow。你应该展示一个完整的程序，并描述你如何编译它并显示你得到的实际输出。不要忘记输入数据。请务必阅读帮助页面http://stackoverflow.com/tour http://stackoverflow.com/help/how-to-ask –

也请花些功夫来正确设置您的文章和代码的格式。我将一些基本缩进应用于您的代码。 –

问题是中期术语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会变大。这将是很好的表现。

来源

2016-12-06 00:33:34 RussF

回答

相关问题