[chiragjangiir]

Problem

Vectorized (paired) indexing on arrays that implement both __array_function__ and __array_namespace__ — such as sparse arrays — raises TypeError: Vectorized indexing is not supported on every release since 2024.10.0. It worked in 2024.09.0 (reported in #9694):

da = xr.DataArray.from_series(s, sparse=True)
da.sel(x=x_idx, y=y_idx)  # TypeError: Vectorized indexing is not supported

Root cause

#9530 made as_indexable prefer __array_namespace__ over __array_function__ for arrays implementing both protocols, so sparse arrays are now wrapped in ArrayApiIndexingAdapter. Its _vindex_get unconditionally raised TypeError, even though the same array fully supports paired indexing through the numpy protocol it also implements.

Solution

When the wrapped array also implements __array_function__, delegate _vindex_get to NdArrayLikeIndexingAdapter — the same "use __array_function__ as a fallback for missing Array API functionality" approach #9530 established elsewhere. Arrays that only implement the Array API keep raising the existing TypeError (paired indexing is genuinely not part of the standard; that case is tracked separately). The result stays a COO array, so sparsity is preserved exactly as before 2024.10.0. No new API, 6 lines in xarray/core/indexing.py.

Regression coverage

Three tests added to xarray/tests/test_indexing.py:

• test_array_api_adapter_vindex_falls_back_to_array_function — adapter-level: a dual-protocol array supports paired indexing (failed with the reported TypeError before the fix).
• test_vectorized_indexing_sparse_array — end-to-end: .sel with DataArray indexers on a sparse-backed DataArray (the exact reported case); asserts the values and that the result is still sparse.COO (failed with the reported TypeError before the fix).
• test_array_api_adapter_vindex_without_array_function_raises — guards the remaining error path: array-API-only arrays still raise TypeError: Vectorized indexing is not supported.

The issue's MVCE now reproduces the exact 2024.09.0 output (dense array([9., 3., 5.]) and <COO: shape=(3,), nnz=3>).

Validation

Run locally on macOS/arm64, Python 3.12, against this branch:

• python -m pytest xarray/tests/ -q -n 4 → 10816 passed, 7020 skipped, 143 xfailed, 25 xpassed, 0 failed, 0 errors (skips are optional backends not installed locally: netCDF4, zarr, cupy, …)
• python -m pytest xarray/tests/test_indexing.py → 183 passed
• ruff check + ruff format --check on the changed files → clean
• pre-commit run --files … → all hooks pass (ruff, blackdoc, rst, typos, …)
• mypy xarray/core/indexing.py → no errors in the changed file (remaining diagnostics are missing optional-dependency stubs for libraries not installed locally, e.g. matplotlib)

Related issue

Closes #9694.

[headtr1ck]

Thanks @chiragjangiir! Falling back to __array_function__ looks right for sparse, but could we also support arrays that only implement the Array API instead of raising?

Since the 2024.12 revision, the standard supports integer array indexing (data-apis/array-api#669), but only for keys that contain nothing but integers and integer arrays. Slices can be converted to aranges that broadcast along their own axis. The result then already has the vindex dimension order, so no moveaxis is needed:

xp = self.array.__array_namespace__()
device = self.array.device
key = indexer.tuple
shape = np.broadcast_shapes(*(k.shape for k in key if not isinstance(k, slice)))
n_slices = sum(isinstance(k, slice) for k in key)
new_key = []
i_slice = 0
for size, k in zip(self.array.shape, key, strict=True):
    if isinstance(k, slice):
        k = xp.arange(*k.indices(size), device=device)
        new_shape = (
            (1,) * (len(shape) + i_slice) + k.shape + (1,) * (n_slices - i_slice - 1)
        )
        i_slice += 1
    else:
        k = xp.broadcast_to(xp.asarray(k, device=device), shape)
        new_shape = shape + (1,) * n_slices
    new_key.append(xp.reshape(k, new_shape))
return self.array[tuple(new_key)]

With array-api-strict, this gives the same results as NumpyVIndexAdapter, including 2D indexers and slices between the indexed axes. _vindex_set would still have to raise, since the standard doesn't cover assignment with integer arrays.

This doesn't help sparse, though. It claims __array_api_version__ = "2025.12", but doesn't implement integer array indexing, and its namespace has no arange:

import numpy as np
import sparse

x = sparse.COO.from_numpy(np.arange(6.0).reshape(2, 3))
x[np.array([[0], [1]]), np.array([[0, 2]])]
# IndexError: Only one-dimensional iterable indices supported.

So the __array_function__ fallback should stay first, and we can't decide based on __array_api_version__. The same limitation means vectorized indexing on sparse with multi-dimensional indexers still fails with this PR, as it did in 2024.09.

@shoyer, since you brought up data-apis/array-api#669 back then: does this split sound good to you (__array_function__ if available, standard integer array indexing otherwise)? Or would you rather go through the Array API for sparse too, and get it fixed there?

[dcherian]

We have some code to translate indexers that is used for backends. I think that could be used effectively here.

EDIT: I see that's what you have done but it's a bit indirect.

I propose that we update _vindex_get to translate to integer indexers for array-api implementers of the right minimum version.

[chiragjangiir]

Thanks @headtr1ck! I've updated ArrayApiIndexingAdapter._vindex_get to implement integer array indexing for arrays adhering to the Array API standard:

1. Dual-protocol precedence: The __array_function__ fallback remains first so dual-protocol arrays like sparse continue to route through their NumPy-compatible indexing path.
2. Standard integer array indexing: For Array-API-only arrays, slice keys are converted to xp.arange arrays that broadcast along their respective axes, and integer array keys are broadcasted and reshaped, yielding the vectorized dimension order directly without needing moveaxis.
3. Assignment: _vindex_set continues to raise TypeError("Vectorized indexing is not supported") because the standard does not specify integer array assignment.
4. Validation: Added test coverage in xarray/tests/test_indexing.py tested against array-api-strict for 1D and 2D indexers, interleaved slices between indexed axes, and unsupported assignment.

[dcherian]

I thought we did this already? Exmaple #10311 (comment)