Thank you very much, it was really helpful!
On Tuesday, March 8, 2016 at 1:00:53 AM UTC-6, Tomas Lycken wrote: > > The thing with meshgrid, is that it’s terribly inefficient compared to > the equivalent idioms in Julia. Instead of implementing it and keeping with > inefficient habits, embrace the fact that there are better ways to do > things :) > > Since this question is recurring, I think the problem might rather be > documentation. As a user looking for a meshgrid function, where did you > look? What did you look for, more specifically? Would it have helped to > have some examples along the lines of “when you want to do this, using > meshgrid, in MATLAB, here’s how to do the same thing without meshgrid in > Julia”? I could turn the following writeup into a FAQ entry or something, > but it would be nice to know where to place it to actually have it make a > difference. > Making do without meshgrid (and getting out ahead!) > > For starters, if what you are after is a straight-up port, meshgrid is > *very* easy to implement for yourself, so that you don’t have to think > about the rest. *Mind you, this is not the most performant way to do > this!* > > meshgrid(xs, ys) = [xs[i] for i in 1:length(xs), j in length(ys)], [ys[j] for > i in 1:length(xs), j in 1:length(ys)] > > Try it out with X, Y = meshgrid(1:15, 1:8). (If you’re not happy with the > orientation of the axes, just switch i and j in the meshgrid > implementation.) > > However, most of the uses of meshgrid in MATLAB is for stuff where you > want to do a bunch of calculations for each point on a grid, and where each > calculation is dependent on the coordinates of just one point. Usually > something along the lines of > > X, Y = meshgrid(xs, ys) > Z = sin(X) + 2 * Y.^2 ./ cos(X + Y) > > In Julia, this is better done with a list comprehension altogether, > eliminating the need for allocating X and Y: > > Z = [sin(x) + 2 * y^2 / cos(x+y) for x in xs, y in ys] > > If you think it’s now become too difficult to read what the actual > transformation is, just refactor it out into a function: > > foo(x,y) = sin(x) + 2 * y^2 / cos(x+y) > Z = [foo(x,y) for x in xs, y in ys] > > If you do this a lot, you might want to avoid writing the list > comprehension altogether; that’s easy too. Just create a function that does > that for you. With a little clever interface design, you’ll see that this > is even more useful than what you could do with meshgrid. > > inflate(f, xs, ys) = [f(x,y) for x in xs, y in ys] > > # call version 1: with a named function > inflate(foo, xs, ys) > > # call version 2: anonymous function syntax > inflate((x,y) -> sin(x) + 2 * y^2 / cos(x,y), xs, ys) > > # call version 3: with a do block > inflate(xs, ys) do x, y > sin(x) + 2 * y^2 / cos(x,y) > end > > Note that the last version of this call is very flexible; in fact, it’s > just as flexible as having a named function. For example, you can do > everything in multiple steps, allowing for intermediate results that don’t > have to allocate large arrays (as they would have, had you used meshgrid), > making it easier to be explicit about whatever algorithm you’re using. > > However, the two last versions of these calls are not the most performant > as of yet, since they’re using anonymous functions. This is changing, > though, and in 0.5 this won’t be an issue anymore. Until then, for large > arrays I’d use the first call version, with a named function, to squeeze as > much as possible out of this. > > I did a small benchmark of this using a straight list comprehension (i.e. > without the inflate function) vs using meshgrid. Here are the results: > > julia> foo(x,y) = sin(x) + 2y.^2 ./ cos(x+y); > > julia> function bench_meshgrid(N) > X, Y = meshgrid(1:N, 1:N) > foo(X, Y) > end; > > julia> function bench_comprehension(N) > [foo(x,y) for x in 1:N, y in 1:N] > end; > > # warmup omitted > > julia> @time bench_comprehension(1_000); > 0.033284 seconds (6 allocations: 7.630 MB) > > julia> @time bench_meshgrid(1_000); > 0.071993 seconds (70 allocations: 68.666 MB, 5.12% gc time) > > julia> @time bench_comprehension(10_000); > 3.547387 seconds (6 allocations: 762.940 MB, 0.06% gc time) > > julia> @time bench_meshgrid(10_000); > ERROR: OutOfMemoryError() > in .* at arraymath.jl:118 > in foo at none:1 > in bench_meshgrid at none:3 > > Note the memory usage: for the list comprehension, it’s entirely dominated > by the result (10_000^2 * 8 / 1024^2 ~= 762.939 MB), while meshgrid eats > almost ten times as much memory (the exact ratio will vary very much > depending on the function foo and how many intermediate arrays it > creates). For large inputs, meshgrid doesn't even make it. > > // T > > On Monday, March 7, 2016 at 5:35:13 PM UTC+1, Christoph Ortner wrote: > > Have a look at >> http://matlabcompat.github.io >> >> If it doesn't have a mesh grid function, then maybe mesh grid from the >> examples folder could be added to this package. >> >> Christoph >> > >
