fortran - Fortran 中快速向上/向下取整 double ？

Question

在 Fortran 中是否有快速向上/向下舍入的方法？

由于正双数的位表示的线性顺序，可以实现如下舍入。

pinf 和 ninf 是全局常量，分别为+/-无穷大

function roundup(x)
    double precision ,intent(in) :: x
    double precision :: roundup

    if (isnan(x))then 
        roundup = pinf
        return
    end if 
    if (x==pinf)then
        roundup = pinf
        return 
    end if
    if (x==ninf)then
        roundup = ninf
        return 
    end if
    if (x>0)then
        roundup = transfer((transfer(x,1_8)+1_8),1d0)
    else if (x<0) then
        roundup = transfer((transfer(x,1_8)-1_8),1d0)
    else
        if (transfer(x,1_8)==Z'0000000000000000')then
            roundup = transfer((transfer(x,1_8)+1_8),1d0)
        else
            roundup = transfer((transfer(-x,1_8)+1_8),1d0)
        end if 
    end if 
end function roundup

我觉得这不是最好的方法，因为它很慢，但它几乎只使用位操作。

另一种方法是使用乘法和一些 epsiloneps = epsilon (1d0)

function roundup2(x)
    double precision ,intent(in) :: x
    double precision :: roundup2
    if (isnan(x)) then 
        roundup2 = pinf
        return 
    else if (x>=eps) then 
        roundup2 = x*(1d0+eps)
    else if (x<=-eps) then 
        roundup2 = x*(1d0-eps)
    else 
        roundup2 = eps
    end if 
end function roundup2

对于某些x，两个函数返回相同的结果(1d0、158d0)，对于某些则不然(0.1d0、15d0)。

第一个函数更准确，但比第二个函数慢约 3.6 倍(10^9 轮测试中 11.1 秒与 3.0 秒)

    print * ,x,y,abs(x-y)
    do i = 1, 1000000000
        x = roundup(x)
        !y = roundup2(y)
    end do 
    print * ,x,y,abs(x-y)

在不检查 NaN/Infinities 的情况下，第一个功能测试需要 8.5 秒 (-20%)。

我真的很难使用round函数，并且在程序配置文件中花费了很多时间。是否有跨平台的方法可以更快地舍入而不损失精度？

更新

该问题怀疑当时对舍入和舍入的调用，无法对它们重新排序。为了简短起见，我没有提到舍入。

提示:第一个函数使用两个 transfer 函数和一个添加函数。在第二种情况下，它比一次乘法和一次加法慢。当传输对数字位没有任何作用时，为什么传输成本如此之高？是否可以用更快的函数取代传输或完全避免添加调用？

Answer 1

我建议您查看 Fortran 标准 IEEE 浮点内在模块(IEEE_ARITHMETIC、IEEE_FEATURES、IEEE_EXCEPTIONS)。它们提供了 IEEE_SET_ROUNDING_MODE，您可以在其中设置后续操作的舍入模式。理想情况下，您可以使用 IEEE_GET_ROUNDING_MODE 获取当前模式并保存它，设置新模式，执行操作，然后恢复模式。

一些警告 - 更改处理器舍入模式本身就是一个缓慢的操作，但如果您执行一次然后进行多次舍入，那将是一个胜利。并非所有当前的 Fortran 编译器都支持 IEEE 内在模块，但大多数合理的编译器都应该支持。您可能需要告诉编译器您正在使用 IEEE 环境 - 对于 Intel Fortran，请使用“-fp-model strict”。