[Numpy-discussion] Apparently non-deterministic behaviour of complex array multiplication

Olivier Delalleau shish at keba.be
Wed Nov 30 20:57:11 EST 2011


I guess it's just a typo on your part, but just to make sure, you are using
.transpose(), not .transpose, correct?

-=- Olivier

2011/11/30 Karl Kappler <magnetotellurics at gmail.com>

> Hello,
> I am somewhat new to scipy/numpy so please point me in the right direction
> if I am posting to an incorrect forum.
>
> The experience which has prompted my post is the following:
> I have a numpy array Y where the elements of Y are
> type(Y[0,0])
> Out[709]: <type 'numpy.complex128'>
>
> The absolute values of the real and complex values do not far exceed say
> 1e-10.  The shape of Y is (24, 49218).
> When I perform the operation: C = dot(Y,Y.conj().transpose), i.e. I form
> the covariance matrix by multiplying T by its conjugate transpose, I
> sometimes get NaN in the array C.
>
> I can imagine some reasons why this may happen, but what is truly puzzling
> to me is that I will be working in ipython and will execute for example:
> find(isnan(C)) and will be returned an list of elements of C which are
> NaN,
> fine, but then I recalculate C, and repeat the find(isnan(C)) command and
> I get a different answer.
>
> I type:
> find(isnan(dot(Y,Y.conj().transpose)))
> and an empty array is returned.  Repeated calls of the same command
> however result in a non-empty array.  In fact, the sequence of arrays
> returned from several consecutive calls varies. Sometimes there are tens of
> NaNs, sometimes none.
>
> I have been performing a collection of experiments for some hours and
> cannot get to the bottom of this;
> Some things I have tried:
> 1. Cast the array Y as a matrix X and calculate X*X.H --- in this case i
> get the same result in that sometimes I have NaN and sometimes I do not.
> 2. set A=X.H and calculate X*A --- same results*
> 3. set b=A.copy() and calc X*b --- same results*.
> 4. find(isnan(real(X*X.H))) --- same results*
> 5. find(isnan(real(X)*real(X.H))) - no NaN appear
>
> *N.B. "Same results" does not mean that the same indices were going NaN,
> simply that I was getting back a different result if I ran the command say
> a dozen times.
>
> So it would seem that it has something to do with the complex
> multiplication.   I am wondering if there is too much dynamic range being
> used in the calculation?  It absolutely amazes me that I can perform the
> same complex-arithmetic operation sitting at the command line and obtain
> different results each time.  In one case I ran a for loop where I
> performed the multiplication 1000 times and found that 694 trials had no
> NaN and 306 trials had NaN.
>
> Saving X to file and then reloading it in a new ipython interpreter
> typically resulted in no NaN.
>
> For a fixed interpretter and instance of X or Y, the indices which go NaN
> (when they do) sometimes repeat many times and sometimes they vary
> apparently at random.
>
> Also note that I have had a similar problem with much smaller arrays, say
> 24 x 3076
>
> I have also tried 'upping' the numpy array to complex256, I have like 12GB
> of RAM...
>
> This happens both in ipython and when I call my function from the command
> line.
>
> Does this sound familiar to anyone?  Is my machine possessed?
>
>
> _______________________________________________
> NumPy-Discussion mailing list
> NumPy-Discussion at scipy.org
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20111130/c58b756f/attachment.html>


More information about the NumPy-Discussion mailing list