using numpy functions on an array of objects
Hey, What is the best solution to get this code working? Anyone a good idea? ------------------ test.py ----------------------------------------------- import numpy import numpy.linalg class afloat: def __init__(self,x): self.x = x def __add__(self,rhs): return self.x + rhs.x def sin(self): return numpy.sin(self.x) def inv(self): return numpy.linalg.inv(self.x) def trace(self): return 0 y = afloat(numpy.eye(3)) z = afloat(numpy.ones(3)) print y + z # works print numpy.sin(y) # works print numpy.trace(y) # doesn't work...??? print numpy.linalg.inv(y) # doesn't work ...??? -------------------- end test.py -------------------------------- === Explanation why I need that === I have the following problem. I need to do numerical calculations on generalized versions of real numbers. In particular with truncated Taylor polynomials. I've implemented that as a class that I called TC. To define what the multiplication of two Taylor polynomials is I use operator overloading. Additionally, I need to compute sine, cosine, exp, etc. of Taylor polynomials. For that, I can use the numpy functions. Numpy is apparently smart enough to call the member function sin(self) of my class afloat when it realizes that the argument of numpy.sin is not a known type. This is really great. However, some functions are not as smart: Among them trace inv dot As a workaround, i could this: def inv(X): if X.__class__ == afloat: return X.inv() else: return numpy.linalg.inv(X) This is somewhat OK, but I'd like to use already existing Python code that uses Numpy internally. So I have to rely that numpy.inv(X) calls X.inv() when X is object of my class. best regards, Sebastian
I think you want to subclass an ndarray here. It's a bit tricky to so, but if you look in the wiki and these mailing list archives, you'll find advise on how to do it. -Chris -- Christopher Barker, Ph.D. Oceanographer Emergency Response Division NOAA/NOS/OR&R (206) 526-6959 voice 7600 Sand Point Way NE (206) 526-6329 fax Seattle, WA 98115 (206) 526-6317 main reception Chris.Barker@noaa.gov
On Fri, Jan 30, 2009 at 13:18, Christopher Barker
I think you want to subclass an ndarray here. It's a bit tricky to so, but if you look in the wiki and these mailing list archives, you'll find advise on how to do it.
That still won't work. numpy.linalg.inv() simply does a particular algorithm on float and complex arrays and nothing else. -- 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
Wouldn't it be nice to have numpy a little more generic?
All that would be needed was a little check of the arguments.
If I do:
numpy.trace(4)
shouldn't numpy be smart enough to regard the 4 as a 1x1 array?
numpy.sin(4) works!
and if
x = my_class(4)
wouldn't it be nice if
numpy.trace(x)
would call
x.trace() ?
numpy.sin(my_class(4)) works!
Wouldn't it be nice if numpy worked a little more consistent.
Is this worth a ticket? Or am I missing something here?
On Fri, Jan 30, 2009 at 10:05 PM, Robert Kern
On Fri, Jan 30, 2009 at 13:18, Christopher Barker
wrote: I think you want to subclass an ndarray here. It's a bit tricky to so, but if you look in the wiki and these mailing list archives, you'll find advise on how to do it.
That still won't work. numpy.linalg.inv() simply does a particular algorithm on float and complex arrays and nothing else.
-- 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
On Sat, Jan 31, 2009 at 10:30, Sebastian Walter
Wouldn't it be nice to have numpy a little more generic? All that would be needed was a little check of the arguments.
If I do: numpy.trace(4) shouldn't numpy be smart enough to regard the 4 as a 1x1 array?
Why? It's not a 1x1 array. It's a scalar. If you want a 1x1 array, give it a 1x1 array.
numpy.sin(4) works!
Yes, numpy.sin() operates on scalars in addition to arrays.
and if x = my_class(4)
wouldn't it be nice if
numpy.trace(x) would call x.trace() ?
numpy.sin(my_class(4)) works!
Wouldn't it be nice if numpy worked a little more consistent. Is this worth a ticket? Or am I missing something here?
numpy.sin() is a ufunc. Unary ufuncs will call the method of the same name on objects in an object array (or the scalar itself if given an object scalar). For example: In [8]: class MyClass(object): ...: def __init__(self, x): ...: self.x = x ...: def __repr__(self): ...: return 'MyClass(%r)' % (self.x,) ...: def sin(self): ...: return MyClass(self.x+1) ...: ...: In [9]: sin(MyClass(4)) Out[9]: MyClass(5) In [10]: sin([MyClass(4), MyClass(5)]) Out[10]: array([MyClass(5), MyClass(6)], dtype=object) You'll notice that numpy.sin() does not try to call the list.sin() method when given the list. It interprets it as an object array, and calls the MyClass.sin() method on each of the elements. numpy.trace() is not an unary ufunc. It's just a function that operates on (N>=2)-D arrays. You simply couldn't apply the same rules as numpy.sin(). Otherwise, it would try to call the .trace() method on each of the objects in your container, and obviously you can't implement trace that way. Having numpy.trace(x) simply call x.trace() would not be making numpy more consistent. Now, that said, the implementation of numpy.trace(x, *args) is actually simply asarray(x).trace(*args). That should probably be asanyarray(x) in order to allow ndarray subclasses. But this only works because ndarray.trace() already exists. Making every function in numpy check for a method first is just not going to happen. -- 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
On Sun, Feb 1, 2009 at 12:24 AM, Robert Kern
On Sat, Jan 31, 2009 at 10:30, Sebastian Walter
wrote: Wouldn't it be nice to have numpy a little more generic? All that would be needed was a little check of the arguments.
If I do: numpy.trace(4) shouldn't numpy be smart enough to regard the 4 as a 1x1 array?
Why? It's not a 1x1 array. It's a scalar. If you want a 1x1 array, give it a 1x1 array.
numpy.sin(4) works!
Yes, numpy.sin() operates on scalars in addition to arrays.
and if x = my_class(4)
wouldn't it be nice if
numpy.trace(x) would call x.trace() ?
numpy.sin(my_class(4)) works!
Wouldn't it be nice if numpy worked a little more consistent. Is this worth a ticket? Or am I missing something here?
numpy.sin() is a ufunc. Unary ufuncs will call the method of the same name on objects in an object array (or the scalar itself if given an object scalar). For example:
In [8]: class MyClass(object): ...: def __init__(self, x): ...: self.x = x ...: def __repr__(self): ...: return 'MyClass(%r)' % (self.x,) ...: def sin(self): ...: return MyClass(self.x+1) ...: ...:
In [9]: sin(MyClass(4)) Out[9]: MyClass(5)
In [10]: sin([MyClass(4), MyClass(5)]) Out[10]: array([MyClass(5), MyClass(6)], dtype=object)
You'll notice that numpy.sin() does not try to call the list.sin() method when given the list. It interprets it as an object array, and calls the MyClass.sin() method on each of the elements.
numpy.trace() is not an unary ufunc. It's just a function that operates on (N>=2)-D arrays. You simply couldn't apply the same rules as numpy.sin(). Otherwise, it would try to call the .trace() method on each of the objects in your container, and obviously you can't implement trace that way.
Having numpy.trace(x) simply call x.trace() would not be making numpy more consistent.
Now, that said, the implementation of numpy.trace(x, *args) is actually simply asarray(x).trace(*args). That should probably be asanyarray(x) in order to allow ndarray subclasses. But this only works because ndarray.trace() already exists. Making every function in numpy check for a method first is just not going to happen.
Ok I see. I understand your reasoning. Nonetheless, I didn't suggest that trace() and sin() are the same, because they are not. I just wanted to express that they should act the same if the object is of *unknown type*. I mean, numpy.sin(MyClass(3)) works. In the worst of all possible worlds, numpy would raise an exception because MyClass is *not an array or a scalar*. But it doesn't. And that is really cool! It's awesome that numpy works on arbitrary type for the most part. In contrast, if trace(X) encounters an unknown type it simply raises an exception. It could as well try *in the very end* to call X.trace(). I.e. *not* "numpy check for a method first" but *numpy check for a method as last resort*. That wouldn't do any harm, would it? And it is not a major effort to add those simple checks to dot(), trace(), inv() . I could provide a patch. But if that is deemed "not going to happen" for whatever reasons. Is there a good workaround? I.e. if I do import numpy import mypackage can I overwrite the functions of numpy? I mean, that is quite a hack, but is the next *best* option. The reason for that need is, that I am writing a Python module to compute higher order derivatives of functions that are given as an algorithm on matrix operations: e.g. we want to compute the Hessian of the function def f(X,Y,Z): """ X is (N,M) array, Y is (M,K) array, Z is (K, L) array""" V = numpy.dot(X,Y) W = numpy.dot(Y,Z) return numpy.trace(V*W) To do that, I generalized the real numbers to the field of truncated Taylor polynomials of scalars and real matrices to truncated Taylor polynomials of matrices. The theory is explained on http://en.wikipedia.org/wiki/Automatic_differentiation You can have a look at the unit test, e.g. at def test_2x2dot2x2_reverse() on http://github.com/b45ch1/algopy/blob/bd7154e2e7a7e6e1931addc0d9ec0604d488d73... Mtc is a class of Matrix Taylor Coefficients. The class Function is used to build the computational graph of function nodes. I think it would be a nice addition to everyone who is doing scientific programming with numpy because derivatives are often required and divided differences suck for anything but first order on small algorithms and especially if one wants to differentiate solutions of ODEs or PDEs. If trace, dot, inv, etc. don't work that way, ppl would have to define two versions of the function f to make it work. best regards, Sebastian Walter
-- 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
participants (3)
-
Christopher Barker
-
Robert Kern
-
Sebastian Walter