[Numpy-tickets] [NumPy] #968: Enhance meshgrid to generate 3D grids + add option for sparse grids

NumPy numpy-tickets@scipy....
Mon Dec 8 20:43:09 CST 2008

```#968: Enhance meshgrid  to generate 3D grids + add option for sparse grids
-------------------------+--------------------------------------------------
Reporter:  pbrod        |       Owner:  somebody
Type:  enhancement  |      Status:  new
Priority:  normal       |   Milestone:  1.3.0
Component:  numpy.lib    |     Version:  devel
Severity:  normal       |    Keywords:  meshgrid
-------------------------+--------------------------------------------------
I would suggest to replace numpy's meshgrid with the meshgrid found in
SciTools at:

for the following reasons:

The `meshgrid` function in NumPy is limited to two dimensions only, while
the SciTools version can also work with 3D grids. In addition,
the NumPy version of `meshgrid` has no option for generating sparse
grids to conserve memory, like it is in SciTools by specifying the
`sparse` argument:
{{{
>>> xv, yv = meshgrid(linspace(-2,2,5), linspace(-1,1,3), sparse=True)
>>> xv
array([[-2., -1.,  0.,  1.,  2.]])
>>> yv
array([[-1.],
[ 0.],
[ 1.]])
}}}

Actually, this is the default behavior for the `meshgrid` function in
SciTools. In NumPy, however, we will in this case get a full 2D grid:
{{{
>>> xv, yv = numpy.meshgrid(linspace(-2,2,5), linspace(-1,1,3))
>>> xv
array([[-2., -1.,  0.,  1.,  2.],
[-2., -1.,  0.,  1.,  2.],
[-2., -1.,  0.,  1.,  2.]])
>>> yv
array([[-1., -1., -1., -1., -1.],
[ 0.,  0.,  0.,  0.,  0.],
[ 1.,  1.,  1.,  1.,  1.]])
}}}

This is the same result we get by setting `sparse=False` in `meshgrid`
in SciTools.

The NumPy functions `mgrid` and `ogrid` does provide support for,
respectively, full and sparse n-dimensional meshgrids, however,
these functions uses slices to generate the meshgrids rather than
one-dimensional coordinate arrays such as in Matlab. With slices, the
user does not have the option to generate meshgrid with, e.g.,
irregular spacings, like
{{{
>>> x = array([-1,-0.5,1,4,5], float)
>>> y = array([0,-2,-5], float)
>>> xv, yv = meshgrid(x, y, sparse=False)
>>> xv
array([[-1. , -0.5,  1. ,  4. ,  5. ],
[-1. , -0.5,  1. ,  4. ,  5. ],
[-1. , -0.5,  1. ,  4. ,  5. ]])
>>> yv
array([[ 0.,  0.,  0.,  0.,  0.],
[-2., -2., -2., -2., -2.],
[-5., -5., -5., -5., -5.]])
}}}

In addition to the reasons mentioned above, the meshgrid function in
NumPy supports only Cartesian indexing, i.e., x and y, not matrix
indexing, i.e., rows and columns (`mgrid` and `ogrid` supports only
matrix indexing). The `meshgrid` function in SciTools supports both
indexing conventions through the `indexing` keyword argument. Giving
the string `'ij'` returns a meshgrid with matrix indexing, while
`'xy'` returns a meshgrid with Cartesian indexing. The difference is
illustrated by the following code snippet:
{{{
nx = 10
ny = 15

x = linspace(-2,2,nx)
y = linspace(-2,2,ny)

xv, yv = meshgrid(x, y, sparse=False, indexing='ij')
for i in range(nx):
for j in range(ny):
# treat xv[i,j], yv[i,j]

xv, yv = meshgrid(x, y, sparse=False, indexing='xy')
for i in range(nx):
for j in range(ny):
# treat xv[j,i], yv[j,i]
}}}

It is not entirely true that matrix indexing is not supported by the
`meshgrid` function in NumPy because we can just switch the order of
the first two input and output arguments:
{{{
yv, xv = numpy.meshgrid(y, x)
}}}
is the same as
{{{
xv, yv = meshgrid(x, y, sparse=False, indexing='ij')
}}}
However, we think it is clearer to have the logical "x, y"

--
Ticket URL: <http://scipy.org/scipy/numpy/ticket/968>
NumPy <http://projects.scipy.org/scipy/numpy>
The fundamental package needed for scientific computing with Python.
```