Meshgrid function in Julia

I am also confused. I use this for contour in PyPlot, where the help files tell me to use meshgrid. How can I avoid meshgrid without finding a different plotting package?

1 Like

You don’t need a meshgrid for contour plotting with PyPlot — you can pass 1d arrays for the axes and it knows how to broadcast them, and you write 2d functions of the axes variables by broadcast operations as well:

using PyPlot
x = range(0,2,length=100)'  # note ': this is a row vector
y = range(0,4,length=200)
z = @. sin(x) * cos(y)  # broadcasts to 2d array
contour(x, y, z)

image

14 Likes

Thanks. I’m even older than @jlchan and transcribed my stuff from (no joke) 20 year old matlab.

1 Like

This is the corresponding Matlab code, btw:

x = linspace(0, 2, 100);
y = linspace(0, 4, 200).';
z = sin(x) .* cos(y);  % broadcasts to 2d array
contour(x, y, z)

2 Likes

I have seen the error in my ways and will go forth and sin in a different way in the future.

9 Likes

Remember to cos as well!

8 Likes

I will do that as I sit in the sun and work on my tan.

8 Likes

Instead of sinning, sec and ye shall find!

11 Likes

That

is neat! I would have done .* but this trick is much better!

1 Like

Is the discussion in the link available in plain-text anywhere? I don’t seem to have access to the page.

I just found an easier way to generate the grid data:

x = 1:3;
y = 4:6;
z = 7:9;

xv = getindex.(Iterators.product(x, y, z), 1)  # first.(Iterators.product(x, y, z), 1) is also ok
yv = getindex.(Iterators.product(x, y, z), 2)
zv = getindex.(Iterators.product(x, y, z), 3) # last.(Iterators.product(x, y, z), 3) is also ok

And the result is

julia> getindex.(Iterators.product(x, y, z), 1)
3×3×3 Array{Int64, 3}:
[:, :, 1] =
 1  1  1
 2  2  2
 3  3  3

[:, :, 2] =
 1  1  1
 2  2  2
 3  3  3

[:, :, 3] =
 1  1  1
 2  2  2
 3  3  3

julia> getindex.(Iterators.product(x, y, z), 2)
3×3×3 Array{Int64, 3}:
[:, :, 1] =
 4  5  6
 4  5  6
 4  5  6

[:, :, 2] =
 4  5  6
 4  5  6
 4  5  6

[:, :, 3] =
 4  5  6
 4  5  6
 4  5  6

julia> getindex.(Iterators.product(x, y, z), 3)
3×3×3 Array{Int64, 3}:
[:, :, 1] =
 7  7  7
 7  7  7
 7  7  7

[:, :, 2] =
 8  8  8
 8  8  8
 8  8  8

[:, :, 3] =
 9  9  9
 9  9  9
 9  9  9

which is same as np.mgrid’s result. If you want np.meshgrid style, the first two dimension should be swapped.

3 Likes

The easiest way is not to do it at all :wink:

7 Likes

To get meshgrid style, x and y shall be swapped as

xv = getindex.(Iterators.product(y, x, z), 2)
yv = getindex.(Iterators.product(y, x, z), 1)

However, if you really need this, by far the easiest and fastest I’ve tested is list comprehension:

function meshgrid(x, y)
   X = [x for _ in y, x in x]
   Y = [y for y in y, _ in x]
   X, Y
end

function meshgrid2(x, y)
   X = getindex.(Iterators.product(y, x), 2)
   Y = getindex.(Iterators.product(y, x), 1)
   X, Y
end

x = 1:1000
y = 1001:2000
julia> @time meshgrid(x,y);
  0.007225 seconds (5 allocations: 15.259 MiB)

julia> @time meshgrid2(x,y);
  0.028644 seconds (9 allocations: 45.777 MiB, 18.41% gc time)

julia> @btime meshgrid($x,$y);
  1.033 ms (4 allocations: 15.26 MiB)

julia> @btime meshgrid2($x,$y);
  4.810 ms (8 allocations: 45.78 MiB)

There is now a package that implements such a lazy grid for ngrid (with lazy meshgrid as a one-line variation explained in the docs):

Hopefully this will help a few more people transition from Matlab to Julia :slight_smile:

9 Likes

function meshgrid(n,L) #n= number of grids, L= size of interval
a=transpose(range(0,L,length=n))
b=repeat(a,n,1)
return b
end

EDIT: I had a comparison with and without using Meshgrid, but I removed it because there was a mistake in the code and the results were erroneous (on top of using time in the global scope)

The part about LazyGrids being faster seemed to be correct though.

EDIT: using LazyGrids (https://github.com/JuliaArrays/LazyGrids.jl)
made it even faster
0.000021 seconds (15 allocations: 928 bytes)

Don’t time in the global scope. You’re measuring compilation time.

change your name from jlchan to jichan

1 Like

You mean Oji-chan?

2 Likes

:laughing: yup

1 Like