Array and __rmul__ operator in Python Numpy

Question

In a project, I created a class, and I needed an operation between this new class and a real matrix, so I overloaded the __rmul__ function like this

class foo(object):

    aarg = 0

    def __init__(self):
        self.aarg = 1


    def __rmul__(self,A):
        print(A)
        return 0

    def __mul__(self,A):
        print(A)
        return 0

but when I called it, the result wasn't what I expected

A = [[i*j for i in np.arange(2)  ] for j in np.arange(3)]
A = np.array(A)
R = foo()
C =  A * R

Output:

0
0
0
1
0
2

It seems that the function is called 6 times, once for each elements.

Instead, the __mul__ function works greatly

C = R * A

Output:

[[0 0]
 [0 1]
 [0 2]]

If A isn't an array, but only a list of lists, both work fine

A = [[i*j for i in np.arange(2)  ] for j in np.arange(3)]
R = foo()
C =  A * R
C = R * A

Output

[[0, 0], [0, 1], [0, 2]]
[[0, 0], [0, 1], [0, 2]]

I'd really want for my __rmul__ function to work also on arrays (my original multiplication function isn't commutative). How can I solve it?

Answer 1

The behaviour is expected.

First of all you have to understand how an operation like x*y is actually executed. The python interpreter will first try to compute x.__mul__(y). If this call returns NotImplemented it will then try to compute y.__rmul__(x). Except when y is a proper subclass of the type of x, in this case the interpreter will first consider y.__rmul__(x) and then x.__mul__(y).

Now what happens is that numpy treats arguments differently depending on whether or not he thinks the argument are scalar or arrays.

When dealing with arrays * does element-by-element multiplication, while scalar multiplication multiplies all the entry of an array by the given scalar.

In your case foo() is considered as a scalar by numpy, and thus numpy multiplies all elements of the array by foo. Moreover, since numpy doesn't know about the type foo it returns an array with dtype=object, so the object returned is:

array([[0, 0],
       [0, 0],
       [0, 0]], dtype=object)

Note: numpy's array does not return NotImplemented when you try to compute the product, so the interpreter calls numpy's array __mul__ method, which performs scalar multiplication as we said. At this point numpy will try to multiply each entry of the array by your "scalar" foo(), and here's is where your __rmul__ method gets called, because the numbers in the array return NotImplemented when their __mul__ is called with a foo argument.

Obviously if you change the order of the arguments to the initial multiplication your __mul__ method gets called immediately and you don't have any trouble.

So, to answer your question, one way to handle this is to have foo inherit from ndarray, so that the second rule applies:

class foo(np.ndarray):
    def __new__(cls):
       # you must implement __new__
    # code as before

Warning however that subclassing ndarray isn't straightforward. Moreover you might have other side effects, since now your class is an ndarray.

Answer 2

You can define __numpy_ufunc__ function in your class. It works even without subclassing the np.ndarray. You can find the documentation here.

Here is an example based on your case:

class foo(object):

    aarg = 0

    def __init__(self):
        self.aarg = 1

    def __numpy_ufunc__(self, *args):
        pass

    def __rmul__(self,A):
        print(A)
        return 0

    def __mul__(self,A):
        print(A)
        return 0

And if we try it,

A = [[i*j for i in np.arange(2)  ] for j in np.arange(3)]
A = np.array(A)
R = foo()
C =  A * R

Output:

[[0 0]
 [0 1]
 [0 2]]

It works!