On Fri, May 16, 2008 at 11:47 AM, Stuart Brorson <[EMAIL PROTECTED]> wrote: > Hi -- > > Sorry to be a pest with corner cases, but I found another one. > > In this case, if you try to take the arccos of numpy.inf in the > context of a complex array, you get a bogus return (IMO). Like this: > > In [147]: R = numpy.array([1, numpy.inf]) > > In [148]: numpy.arccos(R) > Warning: invalid value encountered in arccos > Out[148]: array([ 0., NaN]) > > In [149]: C = numpy.array([1+1j, numpy.inf]) > > In [150]: numpy.arccos(C) > Warning: invalid value encountered in arccos > Out[150]: array([ 0.90455689-1.06127506j, NaN -Infj]) > > The arccos(numpy.inf) in the context of a real array is OK, but taking > arcocs(numpy.inf) in the context of a complex array should return > NaN + NaNj, IMO. > > Thoughts?
Hmm, this works fine on OS X. This may be a problem with one of the lower-level math functions which we defer to the platform. In [1]: from numpy import * In [2]: arccos(nan+nan*1j) Out[2]: (nan+nanj) In [3]: arccos(nan+0j) Out[3]: (nan+nanj) In [4]: arccos(nan) Out[4]: nan In [5]: arccos([1.0+0j, nan]) Out[5]: array([ 0. -0.j, NaN NaNj]) The implementations of the complex versions are in umathmodule.c.src (for the expanded versions, see umathmodule.c after building); they are all prefixed with "nc_". E.g. the following are calling nc_acos() for doubles. Here is the source: static void nc_acos(cdouble *x, cdouble *r) { nc_prod(x,x,r); nc_diff(&nc_1, r, r); nc_sqrt(r, r); nc_prodi(r, r); nc_sum(x, r, r); nc_log(r, r); nc_prodi(r, r); nc_neg(r, r); return; /* return nc_neg(nc_prodi(nc_log(nc_sum(x,nc_prod(nc_i, nc_sqrt(nc_diff(nc_1,nc_prod(x,x)))))))); */ } I suspect the problem comes from the nc_log() which calls your platform's atan2() for the imaginary part of the result. -- Robert Kern "I have come to believe that the whole world is an enigma, a harmless enigma that is made terrible by our own mad attempt to interpret it as though it had an underlying truth." -- Umberto Eco _______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://projects.scipy.org/mailman/listinfo/numpy-discussion