# Creating a numpy array of 3D coordinates from three 1D arrays

To use numpy mesh grid on the above example the following will work:

np.vstack(np.meshgrid(x_p,y_p,z_p)).reshape(3,-1).T


Numpy meshgrid for grids of more then two dimensions require numpy 1.7. To circumvent this and pulling the relevant data from the source code.

def ndmesh(*xi,**kwargs):
if len(xi) < 2:
msg = 'meshgrid() takes 2 or more arguments (%d given)' % int(len(xi) > 0)
raise ValueError(msg)

args = np.atleast_1d(*xi)
ndim = len(args)
copy_ = kwargs.get('copy', True)

s0 = (1,) * ndim
output = [x.reshape(s0[:i] + (-1,) + s0[i + 1::]) for i, x in enumerate(args)]

shape = [x.size for x in output]

# Return the full N-D matrix (not only the 1-D vector)
if copy_:
mult_fact = np.ones(shape, dtype=int)
return [x * mult_fact for x in output]
else:


Checking the result:

print np.vstack((ndmesh(x_p,y_p,z_p))).reshape(3,-1).T

[[ 1.  2.  8.]
[ 1.  2.  9.]
[ 1.  3.  8.]
....
[ 5.  3.  9.]
[ 5.  4.  8.]
[ 5.  4.  9.]]


For the above example:

%timeit sol2()
10000 loops, best of 3: 56.1 us per loop

%timeit np.vstack((ndmesh(x_p,y_p,z_p))).reshape(3,-1).T
10000 loops, best of 3: 55.1 us per loop


For when each dimension is 100:

%timeit sol2()
1 loops, best of 3: 655 ms per loop
In [10]:

%timeit points = np.vstack((ndmesh(x_p,y_p,z_p))).reshape(3,-1).T
10 loops, best of 3: 21.8 ms per loop


Depending on what you want to do with the data, you can return a view:

%timeit np.vstack((ndmesh(x_p,y_p,z_p,copy=False))).reshape(3,-1).T
100 loops, best of 3: 8.16 ms per loop