# A Perspective on Numpy

In a former life, before joining Enthought, I was a researcher in operator algebras, a field of mathematics which can be briefly summed up as the study of certain sets of “infinite-dimensional matrices.” My first reaction when seeing numpy (actually, numerical python, this was a while back) was “A matrix module! In Python! Show me more!”

But then I looked a little deeper, and my second reaction was “Why on earth would you not have the product of two arrays be the matrix product?” And then: “And why is function application element-wise?”

From the point of view of an operator algebraist everything interesting comes out of non-commutative multiplication, and clearly the “right” way to apply a function to a matrix is

f(A) = U f(D) U *

where A = U D U * is a unitary diagonalization of A, and f is applied just to entries on the diagonal of D. The fact that this only works for unitarily diagonalizable matrices was a minor concern for me at the time… after all I was a pure mathematician!

Since then, particularly in a more recent life as a scientist at a data visualization company, I’ve come to appreciate the way that numpy slices and dices arrays of numbers. And numpy now meets me half-way with a first-class matrix type, so I can at least write A*B and have it mean what I want it to.

But for those who still think that numpy made the wrong choice for applying functions to arrays…

```from numpy import *
def mfunc(func):
def new_func(x):
x = mat(x)
if not allclose(x*x.H, x.H*x):
# matrix is not normal-ish
raise ValueError, "matrix must be normal"
# matrix is normal-ish
e, U = linalg.eig(x)
return U*mat(diag(func(e)))*U.H
return new_func```

Which you can then use as:

`mexp = mfunc(exp)`

for any ufunc, or as a decorator on your own functions:

```@mfunc
def f(x):
return cos(x) + 1.0j*sin(x)```

When I showed this to Travis Oliphant, he pointed out that scipy already has a function linalg.funm which does something similar to the above, but not as a decorator. So if you have scipy, you can instead use:

```import scipy.linalg
def mfunc(func):
def new_func(x):
return scipy.linalg.funm(x, func)
return new_func```
This entry was posted in Enthought Tool Suite, NumPy, Open Source, Python, SciPy on .