在 Julia 中利用 ndgrid/meshgrid 功能

Question

我正在尝试在 Julia 中找到类似于 MATLAB 的

meshgrid

或

ndgrid

的功能。我知道 Julia 在

示例

中定义了 ndgrid，但是当我尝试使用它时，出现以下错误。

UndefVarError：ndgrid 未定义

有人知道如何让内置

ndgrid

函数工作，或者可能是我还没有找到的另一个函数或提供这些方法的库（内置函数是首选）？在这种情况下我宁愿不写自己的。

谢谢！

Answer 1

我们更愿意避免使用这些函数，因为它们分配通常不需要的数组。这些数组中的值具有规则的结构，因此不需要存储；它们只能在迭代期间计算。例如，一种替代方法是编写数组理解：

julia> [ 10i + j for i=1:5, j=1:5 ]
5×5 Array{Int64,2}:
 11  12  13  14  15
 21  22  23  24  25
 31  32  33  34  35
 41  42  43  44  45
 51  52  53  54  55

或者，您可以编写

for

循环，或迭代

product

迭代器：

julia> collect(Iterators.product(1:2, 3:4))
2×2 Array{Tuple{Int64,Int64},2}:
 (1, 3)  (1, 4)
 (2, 3)  (2, 4)

Answer 2

我确实发现有时在 numpy 中使用像

meshgrid

这样的函数很方便。通过列表理解很容易做到：

function meshgrid(x, y)
    X = [i for i in x, j in 1:length(y)]
    Y = [j for i in 1:length(x), j in y]
    return X, Y
end

例如

x = 1:4
y = 1:3
X, Y = meshgrid(x, y)

现在

julia> X
4×3 Array{Int64,2}:
 1  1  1
 2  2  2
 3  3  3
 4  4  4
julia> Y
4×3 Array{Int64,2}:
 1  2  3
 1  2  3
 1  2  3
 1  2  3

但是，我没有发现这使得代码比使用迭代运行得更快。这就是我的意思：

定义后

x = 1:1000
y = x
X, Y = meshgrid(x, y)

我对以下两个函数进行了基准测试

using Statistics

function fun1()
    return mean(sqrt.(X.*X + Y.*Y))
end

function fun2()
    sum = 0.0
    for i in 1:1000
        for j in 1:1000
            sum += sqrt(i*i + j*j)
        end
    end
    return sum / (1000*1000)
end

以下是基准测试结果：

julia> @btime fun1()
  8.310 ms (19 allocations: 30.52 MiB)
julia> @btime run2()
  1.671 ms (0 allocations: 0 bytes)

meshgrid

方法明显较慢并且占用更多内存。任何 Julia 专家都知道为什么吗？我知道 Julia 是一种与 Python 不同的编译语言，因此迭代不会比向量化慢，但我不明白为什么向量（数组）计算比迭代慢很多倍。 （对于更大的 N，这种差异甚至更大。）

编辑：阅读完这篇文章后，我有了以下“meshgrid”方法的更新版本。这个想法是不预先创建网格，而是通过 Julia 强大的元素数组操作在计算中完成它：

x = collect(1:1000)
y = x'

function fun1v2()
    mean(sqrt.(x .* x .+ y .* y))
end

这里的技巧是大小为 M 的列数组和大小为 N 的行数组之间的

.+

，后者返回一个 M×N 数组。它为你做“网格”。此功能比

fun1

快近 3 倍，尽管没有

fun2

快。

julia> @btime fun1v2()
  3.189 ms (24 allocations: 7.63 MiB)
765.8435104896155

Answer 3

在上面，@ChrisRackauckas 建议执行此操作的“正确方法”是使用懒惰的操作员，但他还没有抽出时间来执行此操作。

现在有一个注册包，里面有lazy

ndgrid

： https://github.com/JuliaArrays/LazyGrids.jl

比中的版本更通用 矢量化例程.jl 因为它可以处理不同类型的向量，例如，

ndgrid(1:3, Float16[0:2], ["x", "y", "z"])

。

文档中有

Literate.jl

示例，表明惰性性能非常好。

当然懒

meshgrid

就差一步了：

meshgrid(y,x) = (ndgrid_lazy(x,y)[[2,1]]...,)

Answer 4

给定两个轴

x = 1:4
y = 1:3

可以使用 生成 meshgrid

X = x' .+ 0y
Y = 0*x' .+ y

[将

x

与

0

相乘需要星号，因为

0x

试图解析为十六进制数（例如

0x0A

）并会产生

invalid numeric constant

]
这段代码产生

3×4 Matrix{Int64}:
 1  2  3  4
 1  2  3  4
 1  2  3  4
3×4 Matrix{Int64}:
 1  1  1  1
 2  2  2  2
 3  3  3  3

要生成扁平网格，可以使用

X = repeat(x, length(y))
Y = repeat(y, inner=length(x))

# [1, 2, 3, 4, 1, 2, 3, 4, 1, 2, 3, 4]
# [1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3]

或成对

julia> [(xi, yi) for yi=y for xi=x]
12-element Vector{Tuple{Int64, Int64}}:
 (1, 1)
 (2, 1)
 (3, 1)
 (4, 1)
 (1, 2)
 (2, 2)
 (3, 2)
 (4, 2)
 (1, 3)
 (2, 3)
 (3, 3)
 (4, 3)