Sympy: slicing ndarrays

[mgenet]

Hello,
I am confused by the behavior of ndArrays:

a = sympy.Array([1,2])
b = sympy.MutableDenseNDimArray.zeros(4)
b[0::2] # returns [0 0], OK
b[0::2] = a # fails with "object of type 'slice' has no len()"
b[0::2,0] = a # fails with "Index [2] out of border"

What would be the proper way to do that? In numpy this works

a = numpy.array([1,2])
b = numpy.zeros(4)
b[0::2] = a

Thanks a lot!
Martin

[naveensaigit]

The indexing in SymPy is slightly different when compared to standard python/numpy indexing. Other similar issues are - #18710, #18654, #18637, #18459, #18411
I think sympy indexing should be aligned to numpy indexing since most people are familiar with that. @oscarbenjamin What do you think?

[asmeurer]

It looks like multiple bugs. I don't know what "out of border" means. We should try to use error messages that are similar to NumPy, if not the exact same.

[asmeurer]

Is there a specific way that SymPy is different from NumPy, or is it just a buggy implementation of indexing? For NDimArray I think it should match NumPy exactly. For Matrix, I think we do something different, where a single integer index indexes across the whole matrix, rather than giving a single column. If we want to change this to match NumPy, it would be a backwards compatibility break. These are separate issues, though, as they are unrelated objects.

If we are doing manipulation of slices or other complicated index manipulation, we should consider using ndindex, as it has semantics on indexing that match NumPy exactly, and it is very well tested so it is basically guaranteed to be correct (n.b. ndindex is a library that I wrote).

[oscarbenjamin]

We should be careful to distinguish between slicing arrays and slicing matrices. Sympy's array types should probably match numpy's ndarray. Slicing a matrix always gives a (2D) matrix whereas slicing a 2D array gives a 1D array.

This issue just looks like a bug in the array class that can be fixed.

Changing the behaviour for matrix would not be backwards compatible although we should think about what the ideal behaviour should be so that it could be changed in future if a new matrix class is implemented or if there is a backwards incompatible change at some point.

[asmeurer]

If I were to implement matrix slicing from scratch, I would make it so that it matches NumPy semantics, with the exception that any index that would return a 1d array in NumPy would return a (1, n) or (n, 1) (as appropriate) Matrix. But I'm not sure if it's worth changing compatibility to fix this. It's unfortunate how this API was designed, but it's been this way for a long time. We can add a special attribute to matrix like M.idx that lets you index in a NumPy-like way, like M.idx[0] which would give the first row rather than the first element.