Fortran COMPLEX 计算不同于 C++
Fortran COMPLEX calculates different from C++
我已经完成了从 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 专家吗?谢谢!
我建议将所有 Fortran 内容更改为 DP
1.23456789_8(或 1.23456789D00)等
并使用 DIMAG 代替 AIMAG
在 Fortran 代码中替换
c = (1.23456789, 3.45678901)
和
c = (1.23456789d0, 3.45678901d0)
如果没有 kind
,您在 rhs 上使用的真实文字很可能是 32 位实数,您可能需要 64 位实数。后缀 d0
使编译器创建最接近您提供的值的 64 位实数。我已经掩盖了其中的一些细节,还有其他(可能更好)的方法来指定实数文字的种类,但这种方法应该适用于任何当前的 Fortran 编译器。
我不是很懂C++,不知道C++代码是否有同样的问题。
如果我没看错你的问题,这两个代码对 8sf(单精度的极限)产生相同的答案。
至于 IEEE-754 合规性,据我所知,该标准并未涵盖复杂算术问题。我希望在大多数情况下,在幕后使用的 f-p 算术在预期误差范围内产生复数结果,但我不知道它们是有保证的,因为 f-p 算术的误差范围是。
我已经完成了从 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 专家吗?谢谢!
我建议将所有 Fortran 内容更改为 DP
1.23456789_8(或 1.23456789D00)等
并使用 DIMAG 代替 AIMAG
在 Fortran 代码中替换
c = (1.23456789, 3.45678901)
和
c = (1.23456789d0, 3.45678901d0)
如果没有 kind
,您在 rhs 上使用的真实文字很可能是 32 位实数,您可能需要 64 位实数。后缀 d0
使编译器创建最接近您提供的值的 64 位实数。我已经掩盖了其中的一些细节,还有其他(可能更好)的方法来指定实数文字的种类,但这种方法应该适用于任何当前的 Fortran 编译器。
我不是很懂C++,不知道C++代码是否有同样的问题。
如果我没看错你的问题,这两个代码对 8sf(单精度的极限)产生相同的答案。
至于 IEEE-754 合规性,据我所知,该标准并未涵盖复杂算术问题。我希望在大多数情况下,在幕后使用的 f-p 算术在预期误差范围内产生复数结果,但我不知道它们是有保证的,因为 f-p 算术的误差范围是。