Consider the following Python code that multiplies two complex numbers:
import numpy as np
a = np.matrix('28534314.10478439+28534314.10478436j').astype(np.complex128)
b = np.matrix('-1.39818115e+09+1.39818115e+09j').astype(np.complex128)
#Verify values
print(a)
print(b)
c=np.dot(a.getT(),b)
#Verify product
print(c)
Now the product should be -7.979228021897728000e+16 + 48j
which is correct when I run on Spyder. However, if I receive the values a
and b
from a sender to a receiver via MPI on an MPI4py program (I verify that they have been received correctly) the product is wrong and specifically -7.97922801e+16+28534416.j
. In both cases I am using numpy 1.14.3 and Python 2.7.14. The only difference in the latter case is that prior to receiving the values I initialize the matrices with:
a = np.empty_like(np.matrix([[0]*(1) for i in range(1)])).astype(np.complex128)
b = np.empty_like(np.matrix([[0]*(1) for i in range(1)])).astype(np.complex128)
and then the function MPI::Comm::Irecv()
gives them the correct values.
What is going wrong in the latter case if the a
and b
are correct but c
is wrong? Is numpy arbitrarily setting the imaginary part since it's quite smaller than the real part of the product?
First, this doesn't address the mp
stuff, but since it was raised in comments:
np.matrix
can takes a string argument, and produce a numeric matrix from it. Also notice that the shape is (1,1)
In [145]: a = np.matrix('28534314.10478439+28534314.10478436j')
In [146]: a
Out[146]: matrix([[28534314.10478439+28534314.10478436j]])
In [147]: a.dtype
Out[147]: dtype('complex128')
String input to np.array
produces a string:
In [148]: a = np.array('28534314.10478439+28534314.10478436j')
In [149]: a
Out[149]: array('28534314.10478439+28534314.10478436j', dtype='<U36')
But omit the quotes and we get the complex array, with shape () (0d):
In [151]: a = np.array(28534314.10478439+28534314.10478436j)
In [152]: a
Out[152]: array(28534314.10478439+28534314.10478436j)
In [153]: a.dtype
Out[153]: dtype('complex128')
And the product of these values:
In [154]: b = np.array(-1.39818115e+09+1.39818115e+09j)
In [155]: a*b # a.dot(b) same thing
Out[155]: (-7.979228021897728e+16+48j)
Without using mp
, I assume the initialization and setting is something like this:
In [179]: x=np.empty_like(np.matrix([[0]*(1) for i in range(1)])).astype(np.complex128)
In [180]: x[:]=a
In [181]: x
Out[181]: matrix([[28534314.10478439+28534314.10478436j]])
In [182]: y=np.empty_like(np.matrix([[0]*(1) for i in range(1)])).astype(np.complex128)
In [183]: y[:]=b
In [184]: y
Out[184]: matrix([[-1.39818115e+09+1.39818115e+09j]])
In [185]: x*y
Out[185]: matrix([[-7.97922802e+16+48.j]])
It may be worth trying np.zeros_like
instead of np.empty_like
. That will ensure the imaginary part is 0, instead of something random. Then if the mp
process is just setting the real part, you should get something different.