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)
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 
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 
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
yup
1 Like