Bh107 wishlist

[dionhaefner]

This is a growing list of what I am personally missing from the new bh107 interface.

I will probably add some of these in a PR, but it is unlikely that I'll find the time to do them all.

Showstoppers

(critical features I would need to adopt bh107 in Veros)

• [x] Cumulative functions are missing, e.g. add.accumulate / cumsum

• [x] Element-wise comparisons are missing:

  >>> bh_array == 0
  False

• [ ] Boolean masking is missing:

  >>> b[bh107.ones(b.shape, dtype=bool)]
  IndexError: Only integers, slices (`:`), ellipsis (`...`), np.newaxis (`None`) and integer or boolean arrays are valid indices

• [x] Creating empty axes with newaxis / None doesn't work:

  >>> b[np.newaxis, :] * b
  TypeError: can only concatenate list (not "tuple") to list

• [ ] Missing ufuncs: mean, std

Quality of Life

• [ ] A readable __repr__ on BhArray
• [x] Conversions to NumPy by calling np.array(bh_array)
• [ ] Fancy indexing support
• [ ] Method versions of ufuncs: bh_array.sum()
• [ ] Abbreviations for ufuncs: sum instead of add.reduce, abs instead of absolute

[dionhaefner]

I am planning to contribute an implementation leveraging __array_ufunc__ and __array_function__ for bh107. The following minimal example works:

import os
os.environ['NUMPY_EXPERIMENTAL_ARRAY_FUNCTION'] = '1'

import numpy as np

class MyArray:
    def __init__(self, array):
        self.array = array

    def __array_function__(self, func, types, args, kwargs):
        cls = self.__class__
        args = (a.array if isinstance(a, cls) else a for a in args)
        return cls(func(*args, **kwargs))

    def __array_ufunc__(self, ufunc, method, *inputs, **kwargs):
        cls = self.__class__
        inputs = (a.array if isinstance(a, cls) else a for a in inputs)
        return cls(getattr(ufunc, method)(*inputs, **kwargs))

    def __repr__(self):
        if np.isscalar(self.array):
            return repr(self.array)
        return self.__class__.__name__ + repr(self.array)[5:]

Then, you can use it as

>>> a = MyArray(np.random.rand(10))
>>> np.divide(a, np.ones(10))
MyArray([0.22199372, 0.14071593, 0.45037375, 0.2900408 , 0.28482036,
       0.1315491 , 0.65760868, 0.29693582, 0.780299  , 0.57618724])
>>> np.add.accumulate(a)
MyArray([0.22199372, 0.36270965, 0.8130834 , 1.1031242 , 1.38794456,
       1.51949366, 2.17710234, 2.47403816, 3.25433716, 3.8305244 ])
>>> np.argsort(a)
MyArray([5, 1, 0, 4, 3, 7, 2, 9, 6, 8])
>>> np.allclose(a, np.ones(100))
False

That way, we would get sane default implementations for virtually all NumPy functions without -m bohrium or import bohrium as np magic, with the possibility to implement high-performance drop-ins for everything we want to support.

[madsbk]

That would be great!