SimpleArray should handle negative NumPy strides and Python slice semantics consistently

[ThreeMonth03]

SimpleArray does not consistently follow Python/NumPy slicing semantics when negative bounds, negative steps, slice offsets, or negative NumPy strides are involved. Therefore, SimpleArray only works reliably when bounds, steps, and strides are positive and the array view starts from the beginning of the underlying buffer.

Affected APIs include:

• SimpleArray.__getitem__, implemented by SimpleArray::value_type const & at(ssize_t it)

  • Example: sarr[-8,-1,-2] = rhs

• SimpleArray.__setitem__, implemented by ArrayPropertyHelper::setitem_parser()

  • Example: sarr[-8:-1:2, 8:1:-2] = rhs

• SimpleArray(array=ndarray) constructor

  • Example: SimpleArrayFloat64(array=ndarr[::-1])

Tasks

This issues could be devided to a lot of tasks:

• Support negetive index for SimpleArray.__getitem__.
• Remove the usage of SimpleArray::value_type const & at(size_t it).
• Change stride storage from small_vector<size_t> to
  small_vector<ssize_t>.
• Change shape and indices storage from small_vector<size_t> to
  small_vector<ssize_t>.
• Support negative strides for as_mdspan.
• Add a local workaround for SimpleArray.setitem.
• Verify that the related issues continue to behave correctly.

Related Issues

• issue [Bug] Wrong validate_range() in SimpleArray #580
• issue Bug: F-contigous check in SimpleArray breaks #603

How to reproduce the error

$ python3 - <<'PY'
import numpy as np
import modmesh

ndarr_shape = (8, 9, 10, 11, 12)
ndarr = np.arange(np.prod(ndarr_shape), dtype='float64')
ndarr = ndarr.reshape(ndarr_shape)
expected = ndarr.copy()
sarr = modmesh.SimpleArrayFloat64(array=ndarr)

key = np.index_exp[-8:-1:2, 8:1:-2, 1:-1:3, ..., -2:-8:-2]
value_shape = expected[key].shape

value = np.empty(value_shape, dtype='float64')
next_value = 1000.0
for i0 in range(value_shape[0]):
    for i1 in range(value_shape[1]):
        for i2 in range(value_shape[2]):
            for i3 in range(value_shape[3]):
                for i4 in range(value_shape[4]):
                    value[i0, i1, i2, i3, i4] = next_value
                    next_value += 1

sarr[key] = value
expected[key] = value
np.testing.assert_array_equal(sarr.ndarray, expected)
PY
Traceback (most recent call last):
  File "<stdin>", line 23, in <module>
RuntimeError: Broadcast input array from shape(4, 4, 3, 11, 3) into shape(4, 4, 1, 11, 3)

$ python3 - <<'PY'
import numpy as np
import modmesh

ndarr_shape = (16, 18, 20, 22, 24)
ndarr = np.arange(np.prod(ndarr_shape), dtype='float64')
ndarr = ndarr.reshape(ndarr_shape)

view = ndarr[15::-2, 1::2, 19::-2, 1::2, 23::-2]
expected = np.array(view, copy=True)
sarr = modmesh.SimpleArrayFloat64(array=view)

key = np.index_exp[::-2, 1:-1:2, -8:-1:2, 10:2:-3, ...]
value_shape = expected[key].shape

value = np.empty(value_shape, dtype='float64')
next_value = 1000.0
for i0 in range(value_shape[0]):
    for i1 in range(value_shape[1]):
        for i2 in range(value_shape[2]):
            for i3 in range(value_shape[3]):
                for i4 in range(value_shape[4]):
                    value[i0, i1, i2, i3, i4] = next_value
                    next_value += 1

sarr[key] = value
expected[key] = value
np.testing.assert_array_equal(view, expected)
PY
Traceback (most recent call last):
  File "<stdin>", line 24, in <module>
RuntimeError: Broadcast input array from shape(4, 4, 4, 3, 12) into shape(18446744073709551612, 18446744073709551615, 4, 3, 12)

[yungyuc]

The leading dimension in SimpleArray supports ghost index, which treats negative indices differently from Python/numpy. Please take this into consideration.

The ghost indices can complicate the support for the Python/numpy negative index value.

[ThreeMonth03]

The leading dimension in SimpleArray supports ghost index, which treats negative indices differently from Python/numpy. Please take this into consideration.

The ghost indices can complicate the support for the Python/numpy negative index value.

I plan to correct SimpleArray.setitem at first, but I'm not familiar with the application of ghost index, so I needs some clarifications.

1. For SimpleArray, should apis __setitem__ and __getitem__ in pybind11 wrapper works like Python/numpy if we don't consider the first dimension?
   I find a counterexample in unittest. The unittest means that SimpleArray didn't support negetive index before. I'm not sure how we treat negative start and negative end in slicing nowadays. By the way, I think it is acceptable for negative step in slice.
   

   modmesh/tests/test_buffer.py

   Lines 625 to 629 in 01f5fba

   	with self.assertRaisesRegex(
   	IndexError,
   	r"SimpleArray: dim 1 in \[0, -1, 0\] < 0"
   	):
   	sarr[0, -1, 0] = 1

2. For SimpleArray, if we just only consider the first dimension, the start and end value in slice should be the exact index among the range(-ngost, nbody).
   If start and end value are not specified, they will be deduced to -ngost or nbody - 1 according to the sign of step.
   Here is what I think. I just don't know how we could align the second point with the first point.

[yungyuc]

For a ghost array (with non-zero nghost), the range of array index is extended from 0 to negative at the beginning. The index range is (-nghost)..(nbody-1) (i.e., a close-open interval [-nghost, nbody)). For a normal array, the index range is 0..(nbody-1) (close-open interval [0, nbody)).

As such, a normal array is also a ghost array whose nghost == 0. Python/numpy negative indexing should work with ghost index in the following way:

sa = SimpleArray(8, dtype='int8')
sa.ndarray[:] = np.arange(8)
sa.nghost = 3
assert sa.nbody == 5
assert sa[-4] == 7  # rewind from back.
assert sa[-3] == 0  # first element in the buffer.

Ghost index does not affect negative stride. It should work the same with or without ghost index.

[ThreeMonth03]

For a ghost array (with non-zero nghost), the range of array index is extended from 0 to negative at the beginning. The index range is (-nghost)..(nbody-1) (i.e., a close-open interval [-nghost, nbody)). For a normal array, the index range is 0..(nbody-1) (close-open interval [0, nbody)).

As such, a normal array is also a ghost array whose nghost == 0. Python/numpy negative indexing should work with ghost index in the following way:

sa = SimpleArray(8, dtype='int8')
sa.ndarray[:] = np.arange(8)
sa.nghost = 3
assert sa.nbody == 5
assert sa[-4] == 7 # rewind from back.
assert sa[-3] == 0 # first element in the buffer.
Ghost index does not affect negative stride. It should work the same with or without ghost index.

Thanks, the explainaion is clear now. After discussion, I found that we may align and correct the behavior of at(ssize_t idx), so I decide to create and adjust the order of the tasks.

modmesh/cpp/modmesh/buffer/SimpleArray.hpp

Lines 1791 to 1813 in ee8174b

	value_type const & at(size_t it) const
	{
	validate_range(it);
	return data(it);
	}
	value_type & at(size_t it)
	{
	validate_range(it);
	return data(it);
	}
	value_type const & at(ssize_t it) const
	{
	validate_range(it);
	it += m_nghost;
	return data(it);
	}
	value_type & at(ssize_t it)
	{
	validate_range(it);
	it += m_nghost;
	return data(it);
	}