0

I am trying to implement a custom datatype with built-in support for numpy using Numpy-C API. My code is based on the quaterion implementation by https://github.com/moble/quaternion. Currently that project has supported functions like 'add', 'multiply', 'divide', 'power'.

>>> a = np.random.rand(2, 2).astype(np.quaternion)
>>> a
array([[quaternion(0.135325974754691, 0, 0, 0),
        quaternion(0.536005227432872, 0, 0, 0)],
       [quaternion(0.86691238892986, 0, 0, 0),
        quaternion(0.306838076780839, 0, 0, 0)]], dtype=quaternion)
>>> a * a # multiply per element
array([[quaternion(0.0183131194433073, 0, 0, 0),
        quaternion(0.287301603835365, 0, 0, 0)],
       [quaternion(0.751537090080076, 0, 0, 0),
        quaternion(0.0941496053625638, 0, 0, 0)]], dtype=quaternion)

But found it really difficult to support matmul (@) using Numpy C API. And the quaternion project does not support the matmul ufunc:

>>> a @ a
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
TypeError: ufunc 'matmul' not supported for the input types, and the inputs could not be safely coerced to any supported types according to the casting rule ''safe''

Is there any Demo code or Document that can help me to write support for np.matmul?

Shen Messy
  • 11
  • 1
  • The `@` delegates to `arr.__matmul__`. All operators delegate to a `.__....` method. I don't know exactly what the relation is between the `ufunc` `np.matmul` and this method. Both are compiled. But for quaternion what do you expect `matmul` to do. I'm a little rusty on those, what's quaternion "matrix product"? You are aware, aren't you that `(x * y).sum(axis=1)` does the matrix product for suitably aligned arrays? – hpaulj May 11 '20 at 07:06
  • Thanks @hpaulj. I am trying to implement my own datatype which accounts on matmul pretty much, and my code is based on the quaternion implementation. I am not sure what you mean by suitably aligned but I believe matmul is something like ``[i,j] = sum(x[i, :] * y[:, j])``. Is there a easier way to replace matmul with multiply and sum without looping on ``i`` and ``j``? – Shen Messy May 11 '20 at 07:53

0 Answers0