Hi, Juan. Meshgrid can actually give what you want, but you must use the options: copy=False and indexing=’ij’.
In [7]: %timeit np.meshgrid(np.arange(512), np.arange(512)) 1000 loops, best of 3: 1.24 ms per loop In [8]: %timeit np.meshgrid(np.arange(512), np.arange(512), copy=False) 10000 loops, best of 3: 27 µs per loop In [9]: %timeit np.meshgrid(np.arange(512), np.arange(512), copy=False, indexing='ij') 10000 loops, best of 3: 23 µs per loop Best regards Per A. Brodtkorb From: NumPy-Discussion [mailto:numpy-discussion-boun...@scipy.org] On Behalf Of Juan Nunez-Iglesias Sent: 9. mars 2017 04:20 To: Discussion of Numerical Python Subject: Re: [Numpy-discussion] Why do mgrid and meshgrid not return broadcast arrays? Hi Warren, ogrid doesn’t solve my problem. Note that my code returns arrays that would evaluate as equal to the mgrid output. It’s just that they are copied in mgrid into a giant array, instead of broadcast: In [176]: a0, b0 = np.mgrid[:5, :5] In [177]: a1, b1 = th.broadcast_mgrid((np.arange(5), np.arange(5))) In [178]: a0 Out[178]: array([[0, 0, 0, 0, 0], [1, 1, 1, 1, 1], [2, 2, 2, 2, 2], [3, 3, 3, 3, 3], [4, 4, 4, 4, 4]]) In [179]: a1 Out[179]: array([[0, 0, 0, 0, 0], [1, 1, 1, 1, 1], [2, 2, 2, 2, 2], [3, 3, 3, 3, 3], [4, 4, 4, 4, 4]]) In [180]: a0.strides Out[180]: (40, 8) In [181]: a1.strides Out[181]: (8, 0) On 9 Mar 2017, 2:05 PM +1100, Warren Weckesser <warren.weckes...@gmail.com<mailto:warren.weckes...@gmail.com>>, wrote: On Wed, Mar 8, 2017 at 9:48 PM, Juan Nunez-Iglesias <jni.s...@gmail.com<mailto:jni.s...@gmail.com>> wrote: I was a bit surprised to discover that both meshgrid nor mgrid return fully instantiated arrays, when simple broadcasting (ie with stride=0 for other axes) is functionally identical and happens much, much faster. Take a look at ogrid: https://docs.scipy.org/doc/numpy/reference/generated/numpy.ogrid.html Warren I wrote my own function to do this: def broadcast_mgrid(arrays): shape = tuple(map(len, arrays)) ndim = len(shape) result = [] for i, arr in enumerate(arrays, start=1): reshaped = np.broadcast_to(arr[[...] + [np.newaxis] * (ndim - i)], shape) result.append(reshaped) return result For even a modest-sized 512 x 512 grid, this version is close to 100x faster: In [154]: %timeit th.broadcast_mgrid((np.arange(512), np.arange(512))) 10000 loops, best of 3: 25.9 µs per loop In [156]: %timeit np.meshgrid(np.arange(512), np.arange(512)) 100 loops, best of 3: 2.02 ms per loop In [157]: %timeit np.mgrid[:512, :512] 100 loops, best of 3: 4.84 ms per loop Is there a conscious design decision as to why this isn’t what meshgrid/mgrid do already? Or would a PR be welcome to do this? Thanks, Juan. _______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org<mailto:NumPy-Discussion@scipy.org> https://mail.scipy.org/mailman/listinfo/numpy-discussion _______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org<mailto:NumPy-Discussion@scipy.org> https://mail.scipy.org/mailman/listinfo/numpy-discussion
_______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org https://mail.scipy.org/mailman/listinfo/numpy-discussion
