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
>
