array of multivectors fails

[arsenovic]

@rkern or @moble, do either of you know why this happens?

import numpy as np
from clifford import Cl 
layout,blades = Cl(2) # note i changed firstIdx=1 by default in Cl
locals().update(blades)

np.array([e1,e2],dtype=np.object) 

yields

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

perhaps there are some special array-like methods that Multivector implements? is there a way to have an array of multivectors, instead of an array of their values?

[rkern]

Multivector objects are also sequences. np.array() consumes sequences as if they were lists. To create a dtype=object array of sequences, the magic of np.array() that automatically determines the shape of the output is interfering. The usual workaround is to create an empty object array of the desired shape, then assign the elements.

vs = [e1, e2]
a = np.empty(len(vs), dtype=object)
a[:] = vs

[arsenovic]

sweet thanks a lot.

fyi, i am currently working on a Frame object. its a light object, and i was looking for simple element-wise arithmetic to work, like i could implement A*B with

[a*b for a,b in zip(A,B)]

but then it seemed that i should just to re-use numpy's array for this.

[arsenovic]

i am re-opening. because clifford is numerical package, we need to make working with arrays easier, so that broadcasting works cleaner.

ie i want something like this to work.

from clifford.g3 import * 
x = linspace(0,1,101)*e3

[arsenovic]

i have tried a few things with __array_ufunc__, and __array_wrap__ , but the difficulty of this problem doesnt make sense to me. all that is needed is for ndarray's to see MultiVectors as numbers instead of sequences, and i think all the broadcasting and ufuncs(that are defined) should work?

is there a way to prevent numpy from seeing the sequence character of MultiVector? MultiVector.__numpy_ignore_mysequenceness_for_broadcasting__==True

[hugohadfield]

@eric-wieser you probably know all about this stuff, any chance of a little help :) ? We are basically aiming to be able to write something like:

(np.random.rand(4)*e1 )+e12

or

e12+(e1*np.random.rand(4))

We currently have an array wrap setup here: https://github.com/pygae/clifford/blob/09945d8ab36934f726386ec189f58ec122263b6a/clifford/init.py#L981

And a numeric coercing thing here: https://github.com/pygae/clifford/blob/19bfa2d42d81df89b57c0b4bbd307b2bcb70c27b/clifford/init.py#L1004

Currently numpy tries to read the MultiVector class as a sequence and broadcast those values, but we would like it to broadcast the whole MultiVector object

[hugohadfield]

Here is my first stab at it: https://github.com/pygae/clifford/pull/89

[eric-wieser]

You should leave array_wrap behind and implement array_ufunc instead

[hugohadfield]

Oof, it's true I saw in the numpy docs it's getting deprecated, when is the array_wrap being phased out? Is array_ufunc basically stable between numpy versions?

[eric-wieser]

__array_wrap__ will stick around forever, but it's not powerful enough to solve the problems you'll be wanting to solve.

[arsenovic]

so, hugo's #89 implementation works great, and i dont want to beat a dead horse. but i dont understand why there is not an easy way to prevent numpy from consuming a custom sequence-like object. it seems like it should be a one-line thing to implement this feature.

[eric-wieser]

I think this will do the trick:

def __array__(self):
    # we are a scalar, and the only appropriate dtype is an object array
    arr = np.empty((), dtype=object)
    arr[()] = self
    return arr

[eric-wieser]

The big problem with __array_wrap__ is it runs after numpy has already done the computation. If you're going to just do your own computation, you should use __array_ufunc__.

But if you're not trying to be an array at all, then telling numpy how to turn yourself into one is sufficient.

[arsenovic]

eric this works great, thanks a lot. one last question. if i want the ndarray of objects to have attributes access of the elements, like

a = np.rand(3)*e1 # returns ndarray of MultiVectors
b = np.array([x.dual() for x in a]) # this is clunky
#  it would be sweet to write 
b = a.dual()

so it would work like ndarray.max(). possible?

[eric-wieser]

No, that's much harder. You should at least use:

dual_array = np.vectorize(MultiVector.dual)
b = dual_array(a)

To actually make those methods, you'd probably have to make a new MultiVectorArray class subclassing ndarray, and make array return one of those.

[arsenovic]

thanks for the definitive reply!