2

我已经完成了从 Fortran 到 C++ 的移植，但发现了 COMPLEX 类型的一些差异。考虑以下代码：

PROGRAM CMPLX
    COMPLEX*16 c
    REAL*8 a
    c = (1.23456789, 3.45678901)
    a = AIMAG(1.0 / c)
    WRITE (*, *) a
END

和 C++：

#include <complex>
#include <iostream>
#include <iomanip>

int main()
{
    std::complex<double> c(1.23456789, 3.45678901);
    double a = (1.0 / c).imag();

    std::cout << std::setprecision(15) << " " << a << std::endl;
}

用 clang++ 或 g++ 编译 C++ 版本，我得到输出：-0.256561150444368 但是编译 Fortran 版本给了我：-0.25656115049876993

我的意思是，这两种语言不是都遵循 IEEE 754 吗？如果我在 Octave (Matlab) 中运行以下命令：

octave:1> c=1.23456789+ 3.45678901i
c =  1.2346 + 3.4568i
octave:2> c
c =  1.2346 + 3.4568i
octave:3> output_precision(15)
octave:4> c
c =  1.23456789000000e+00 + 3.45678901000000e+00i
octave:5> 1 / c
ans =  9.16290109820952e-02 - 2.56561150444368e-01i

我得到的和 C++ 版本一样。Fortran COMPLEX 类型是怎么回事？我错过了一些编译器标志吗？-ffast-math 不会改变任何东西。我想在 C++ 和 Fortran 中生成完全相同的 15 位小数，因此我更容易发现移植差异。

周围有 Fortran 大师吗？谢谢！

4

2 回答 2

5

在 Fortran 代码中替换

c = (1.23456789, 3.45678901)

和

c = (1.23456789d0, 3.45678901d0)

如果没有 akind您在 rhs 上使用的真实文字，很可能是 32 位实数，而您可能需要 64 位实数。后缀d0使编译器创建最接近您提供的值的 64 位实数。我已经忽略了其中的一些细节，还有其他（可能更好）方法来指定实数文字的种类，但是这种方法在任何当前的 Fortran 编译器上都可以正常工作。

我对C++不是很了解，我不确定C++代码是否有同样的问题。

如果我正确阅读了您的问题，这两个代码对 8sf（单精度限制）产生相同的答案。

至于 IEEE-754 合规性，据我所知，该标准并未涵盖复杂算术的问题。我希望在大多数情况下，幕后使用的 fp 算术会在预期的误差范围内产生复数的结果，但我不知道它们是否得到保证，因为 fp 算术的误差范围是。

于 2015-01-12T19:36:00.570 回答
0

我建议将所有 Fortran contants 更改为 DP

1.23456789_8（或 1.23456789D00）等

并使用 DIMAG 而不是 AIMAG

于 2015-01-12T19:34:45.790 回答