# Should numpy.sqrt(-1) return 1j rather than nan?

Christopher Barker Chris.Barker at noaa.gov
Thu Oct 12 12:49:25 CDT 2006

```As someone on this list (sorry to quick with the delete button) said:

"numpy is closer to the metal than MATLAB"

MATLAB aside, numpy is somewhat "close to the metal". I think this is
clearly part of its design philosophy, and also something I personally

Getting back to MATLAB (at least version 4, the last one I used
regularly) -- everything in MATLAB is a matrix of double precision
Floating points. (or complex - that does get weird!). This makes things
simpler for a lot of use, but restricts the power and flexibility
substantially.

Everything in numpy is an array of arbitrary shape and many possible
data types -- this is more to think about, but give you a lot of power
and flexibility. In any case, that is what numpy is. Period.

Given that, arrays of different types behave differently -- that's
something you'd darn well better understand pretty soon in your numpy
career:

>>> a = N.array((1,2,3), dtype=N.uint8)
>>> a *= 200
>>> a
array([200, 144,  88], dtype=uint8)

Oh my gosh! 2 * 200 isn't 144! what happened?

This isn't any different than sqrt(-1) resulting in NaN or an error. In
fact, there are all sorts of places where upcasting is not happening
automagically -- in fact, I think that's become more consistent in the
new numpy. So, in numpy you need to choose the datatype appropriate for
your problem at hand. I know I always know whether I want to work in the
complex plane or not.

Another change recently made is to make doubles the default in more
places -- I like that change.

So, given the entire philosophy and current status of how numpy works,
the ONLY question that is legitimate here is:

Should the "default" datatype be complex or double? I vote double.

Travis Oliphant wrote:
> Now, it would be possible to give ufuncs a dtype keyword argument that
> allowed you to specify which underlying loop was to be used for the
> calculation.

> sqrt(a, dtype=complex).

I like this OK, and isn't is similar to the dtype argument in the
accumulating methods? (sum, etc.). However, it's probably not worth the
C-api change -- is it that different than calling:

sqrt(a.astype(complex)) ?

As for scipy, I heartily agree that scipy should be a superset of numpy,
and NOT change the behavior of numpy functions and methods. We're taking
the opportunity at this point for numpy to break backward compatibility
with Numeric/numarray -- this is probably the best time to   do the same
for scipy.

couldn't we at least introduce new names?

numpy.scimath.csqrt()

for instance?

-Chris

--
Christopher Barker, Ph.D.
Oceanographer

NOAA/OR&R/HAZMAT         (206) 526-6959   voice
7600 Sand Point Way NE   (206) 526-6329   fax
Seattle, WA  98115       (206) 526-6317   main reception

Chris.Barker at noaa.gov

-------------------------------------------------------------------------
Using Tomcat but need to do more? Need to support web services, security?
Get stuff done quickly with pre-integrated technology to make your job easier