Fir*_*man 55 c++ python numpy infinity complex-numbers
我一直在 C++ 和 numpy 中发现关于函数处理复杂无限数log的行为的奇怪行为。log具体来说,log(inf + inf * 1j)等于(inf + 0.785398j)我期望的值(inf + nan * 1j)。
当取复数的对数时,实部是输入的绝对值的对数,虚部是输入的相位。返回 0.785398 作为 的虚部log(inf + inf * 1j)意味着它假设inf实部和虚部中的 s 具有相同的长度。这个假设似乎与其他计算不一致,例如 ,inf - inf == nan它inf / inf == nan假设 2 infs 不一定具有相同的值。
为什么假设log(inf + inf * 1j)不同?
重现 C++ 代码:
#include <complex>
#include <limits>
#include <iostream>
int main() {
double inf = std::numeric_limits<double>::infinity();
std::complex<double> b(inf, inf);
std::complex<double> c = std::log(b);
std::cout << c << "\n";
}
Run Code Online (Sandbox Code Playgroud)
重现 Python 代码 (numpy):
import numpy as np
a = complex(float('inf'), float('inf'))
print(np.log(a))
Run Code Online (Sandbox Code Playgroud)
编辑:感谢所有参与历史原因和数学原因讨论的人。你们所有人都把这个天真的问题变成了一个非常有趣的讨论。提供的答案都是高质量的,我希望我能接受 1 个以上的答案。然而,我决定接受@simon的答案,因为它更详细地解释了数学原因,并提供了解释逻辑的文档的链接(尽管我无法完全理解它)。
And*_*eak 51
C99 规范的免费最终草案第 491 页上说明
\n\n\n\nclog(+\xe2\x88\x9e, +i\xe2\x88\x9e) 返回 +\xe2\x88\x9e + i\xcf\x80/4。
\n
\n\n该函数的语义旨在与 C 函数保持一致
\nclog。
我同意,从数学的角度来看,这种行为令人困惑,并且可以说与其他inf语义不一致,正如您所指出的。但实际上,它是 C 标准的一部分,这使其成为 C++ 标准的一部分,并且由于 NumPy 通常依赖于 C 行为(即使在令人困惑的情况下),因此在 Python 示例中继承了这一点。
标准库cmath.log()函数具有相同的行为(如果你测试正确......):
>>> import cmath\n\n>>> cmath.log(complex(float(\'inf\'), float(\'inf\')))\n(inf+0.7853981633974483j)\nRun Code Online (Sandbox Code Playgroud)\n我无法调查 C 标准的基本原理。我认为在考虑这些复杂功能如何相互作用时,可能会做出务实的选择。
\nsim*_*mon 39
0.785398 的值(实际上是 pi/4)至少与其他一些函数一致:正如您所说,复数的对数的虚部与该数的相位角相同。这可以重新表述为一个自己的问题: 的相位角是多少inf + j * inf?
z我们可以通过计算复数的相位角atan2(Im(z), Re(z))。对于给定的数字,这可以归结为计算atan2(inf, inf),对于 Numpy 和 C/C++ 来说,这也是 0.785398(或 pi/4)。所以现在可以问类似的问题:为什么是atan2(inf, inf) == 0.785398?
我对后者没有答案(除了“C/C++ 规范是这么说的”,正如其他人已经回答的那样),我只有一个猜测:至于atan2(y, x) == atan(y / x),x > 0可能有人在这种情况下决定不解释为inf / inf“未定义”,而是“一个非常大的数字除以同一个非常大的数字”。该比率的结果将为 1,并且atan(1) == pi/4根据 的数学定义atan。
也许这不是一个令人满意的答案,但至少我可以希望表明给log定边缘情况下的定义与相关函数定义的类似边缘情况并不完全不一致。
编辑:正如我所说,与其他一些功能一致:np.angle(complex(np.inf, np.inf)) == 0.785398例如,它也与 一致。
编辑2:查看实际atan2实现的源代码会出现以下代码注释:
请注意,不明显的情况是 y 和 x 均为无穷大或均为零。有关详细信息,请参阅W. Kahan 的《复杂初等函数的分支切割》或《无事生非的符号位》
我找到了参考文档,您可以在这里找到副本。在本参考文献的第 8 章“复数零和无穷大”中,威廉·卡汉(William Kahan,数学家和计算机科学家,根据维基百科,他是“浮点之父”)介绍了复数的零和无穷大边缘情况,以及到达 pi/4 以输入inf + j * inf函数arg(arg是计算复数相位角的函数,就像np.angle上面一样)。您可以在链接的 PDF 的第 17 页上找到此结果。我不是一个足够的数学家,无法总结卡汉的基本原理(也就是说:我并不真正理解它),但也许其他人可以。
小智 5
如果我们从纯数学的角度考虑这一点,那么我们可以根据极限来看待运算,例如当 x 趋向无穷大时,1/x 趋于 0(表示为 lim(x => inf) 1/x = 0),这就是我们通过浮点观察到的情况。
对于 2 个无穷大的运算,我们单独考虑每个无穷大。因此:
lim(x => inf) x/1 = inf
lim(x => inf) 1/x = 0
Run Code Online (Sandbox Code Playgroud)
一般来说,我们说 inf/x = inf,并且 x/inf = 0。因此:
lim(x => inf) inf/x = inf
lim(x => inf) x/inf = 0
Run Code Online (Sandbox Code Playgroud)
我们应该选择这两个中的哪一个?float 的规范通过将其声明为 NaN 来回避。
然而,对于复杂的日志,我们观察到:
lim(x=>inf) log(x + 0j) = inf + 0j
lim(x=>inf) log(0 + xj) = inf + pi/2j
lim(x=>inf) log(inf + xj) = inf + 0j
lim(x=>inf) log(x + infj) = inf + pi/2j
Run Code Online (Sandbox Code Playgroud)
仍然存在矛盾,但不是在 0 和 inf 之间,而是在 0 和 pi/2 之间,因此规范作者选择将差异分开。我不能说为什么他们做出这个选择,但浮点无穷大不是数学无穷大,而是代表“这个数字太大而无法表示”。鉴于 log(complex) 的使用可能比减法和除法更纯粹的数学,作者可能认为保留恒等式 im(log(x+xj)) == pi/4 是有用的。