def indices(dimensions, dtype=int):
dimensions = tuple(dimensions)<EOL>N = len(dimensions)<EOL>if N == <NUM_LIT:0>:<EOL><INDENT>return array([], dtype=dtype)<EOL><DEDENT>res = empty((N,)+dimensions, dtype=dtype)<EOL>for i, dim in enumerate(dimensions):<EOL><INDENT>tmp = arange(dim, dtype=dtype)<EOL>tmp.shape = (<NUM_LIT:1>,)*i + (dim,)+(<NUM_LIT:1>,)*(N-i-<NUM_LIT:1>)<EOL>newdim = dimensions[:i] + (<NUM_LIT:1>,)+ dimensions[i+<NUM_LIT:1>:]<EOL>val = zeros(newdim, dtype)<EOL>add(tmp, val, res[i])<EOL><DEDENT>return res<EOL>
Return an array representing the indices of a grid.

Compute an array where the subarrays contain index values 0,1,...
varying only along the corresponding axis.

Parameters
----------
dimensions : sequence of ints
    The shape of the grid.
dtype : dtype, optional
    Data type of the result.

Returns
-------
grid : ndarray
    The array of grid indices,
    ``grid.shape = (len(dimensions),) + tuple(dimensions)``.

See Also
--------
mgrid, meshgrid

Notes
-----
The output shape is obtained by prepending the number of dimensions
in front of the tuple of dimensions, i.e. if `dimensions` is a tuple
``(r0, ..., rN-1)`` of length ``N``, the output shape is
``(N,r0,...,rN-1)``.

The subarrays ``grid[k]`` contains the N-D array of indices along the
``k-th`` axis. Explicitly::

    grid[k,i0,i1,...,iN-1] = ik

Examples
--------
>>> grid = np.indices((2, 3))
>>> grid.shape
(2, 2, 3)
>>> grid[0]        # row indices
array([[0, 0, 0],
       [1, 1, 1]])
>>> grid[1]        # column indices
array([[0, 1, 2],
       [0, 1, 2]])

The indices can be used as an index into an array.

>>> x = np.arange(20).reshape(5, 4)
>>> row, col = np.indices((2, 3))
>>> x[row, col]
array([[0, 1, 2],
       [4, 5, 6]])

Note that it would be more straightforward in the above example to
extract the required elements directly with ``x[:2, :3]``.