numpy：逐列复杂的img部分为零答案

【问题标题】：numpy: complex row by column img part is zeronumpy：逐列复杂的img部分为零
【发布时间】：2021-09-22 11:56:11
【问题描述】：

我有一行和一列 (1, 3) 和 (3,1) 用这些值

a = np.array([[-0.000147028898936696350574+0.003143749432638287544250j,
               -0.000184992561116814613342+0.001760595943778753280640j,
               -0.000147028898936696350574+0.003143749432638287544250j]], dtype=np.complex64)
b = np.array([[-147.028900146484375000000000-3143.749511718750000000000000j],
              [-184.992553710937500000000000-1760.595947265625000000000000j],
              [-147.028900146484375000000000-3143.749511718750000000000000j]], dtype=np.complex64)

我想把两者相乘

c = np.matmul(a, b, dtype=np.complex64)
print (c)

像这样设置打印选项

np.seterr(all='raise')
np.set_printoptions(precision=24)
np.set_printoptions(floatmode="fixed")
np.set_printoptions(threshold=sys.maxsize)

我明白了

>>> print (c)
[[22.943477630615234375000000+0.000000000000000000000000j]]

如果我像手动那样做乘法

d = a[0][0]*b[0][0] + a[0][1]*b[1][0] + a[0][2]*b[2][0]

我明白了

print (d)
(22.943478+5.9604645e-08j

与匹配的 C 示例一致

#include <complex.h>
#include <stdio.h>

int main()
{
  complex float a[] = { -0.000147028898936696350574+0.003143749432638287544250*I, 
    -0.000184992561116814613342+0.001760595943778753280640*I,
    -0.000147028898936696350574+0.003143749432638287544250*I};

  complex float b[] = {-147.028900146484375000000000-3143.749511718750000000000000*I,
    -184.992553710937500000000000-1760.595947265625000000000000*I,
    -147.028900146484375000000000-3143.749511718750000000000000*I };

  complex float  den  = a[0] * b[0] + a[1] * b[1] +  a[2] * b[2];

  printf("%+.30e  %+.30ej\n", crealf(den), cimagf(den) );
  return 0;
}

x@y:~$ gcc main.c && ./a.out 
+2.294347763061523437500000000000e+01  +5.960464477539062500000000000000e-08j

我如图实部匹配但虚部不匹配，为什么？

【问题讨论】：

我猜这个问题是由于浮点精度造成的。 Numpy 使用一些高度优化的算法来计算矩阵乘法，因此很难比较您的 C 代码和 python 代码。但是，如果您使用琐碎的np.sum(a*b.T)，您将获得（几乎）与C 代码相同的结果。
@obchardon 如果我做d = a[0][0]*b[0][0] + a[0][1]*b[1][0] + a[0][2]*b[2][0] print (d) 我得到>>> print (d) (22.943478+5.9604645e-08j) 但仍然0.000000000000000000000000j 与5.9604645e-08j 有很大不同
需要注意的是complex64不是一个复杂的值，每个部分都是float64（这是我最初假设的），它的总大小是64位，所以它由两个组成float32。您只有 23 位小数，因此您添加的数字之间的“细节”必须远大于数字大小的 ~1/8,000,000。

标签： python c numpy

【解决方案1】：

问题似乎在于精度。在计算 complex64 数字时，numpy 以某种方式将结果四舍五入为 0。如果您使用 dtype=complex128 或 dtype=complex256，虚部与 0 不同。但是将您提供的数字相乘会得到不同的虚部：

a = np.array([[-0.000147028898936696350574+0.003143749432638287544250j,
           -0.000184992561116814613342+0.001760595943778753280640j,
           -0.000147028898936696350574+0.003143749432638287544250j]], 
dtype=np.complex256)
b = np.array([[-147.028900146484375000000000-3143.749511718750000000000000j],
          [-184.992553710937500000000000-1760.595947265625000000000000j],
          [-147.028900146484375000000000-3143.749511718750000000000000j]], dtype=np.complex256)
np.matmul(a, b)
array([[22.94347681+2.93314883e-08j]], dtype=complex256)

更新：您的 C 代码输出和我的 python 的差异再次是由于精度，如果您在 C 代码中使用 double 值，那么输出将与我的相同。

【讨论】：