朱利亚 - 如何高效地将矩阵对角线上的元素变为零？

Question

在Julia中，将矩阵的对角线变为零的高效方法是什么？

- Davi Barreira

Answer 1

假设你有一个大小为 N x N 的矩阵 m，可以按如下方式完成：

setindex!.(Ref(m), 0.0, 1:N, 1:N)

另一个选项：

using LinearAlgebra
m[diagind(m)] .= 0.0

还有一些性能测试：

julia> using LinearAlgebra, BenchmarkTools

julia> m=rand(20,20);

julia> @btime setindex!.(Ref($m), 0.0, 1:20, 1:20);
  55.533 ns (1 allocation: 240 bytes)

julia> @btime $m[diagind($m)] .= 0.0;
  75.386 ns (2 allocations: 80 bytes)

- Przemyslaw Szufel

Answer 2

就性能而言，简单的循环更快（并且更明确，但这取决于个人口味）

julia> @btime foreach(i -> $m[i, i] = 0, 1:20)
  11.881 ns (0 allocations: 0 bytes)

julia> @btime setindex!.(Ref($m), 0.0, 1:20, 1:20);
  50.923 ns (1 allocation: 240 bytes)

它比diagind版本更快，但并没有快多少。

julia> m = rand(1000, 1000);

julia> @btime foreach(i -> $m[i, i] = 0.0, 1:1000)
  1.456 μs (0 allocations: 0 bytes)

julia> @btime foreach(i -> @inbounds($m[i, i] = 0.0), 1:1000)
  1.338 μs (0 allocations: 0 bytes)

julia> @btime $m[diagind($m)] .= 0.0;
  1.495 μs (2 allocations: 80 bytes)

- Andrej Oskin

Answer 3

以下是一种更通用的方法，如何通过访问自定义索引在数组上高效使用setindex!： 使用数组进行数组索引

线性索引是最佳性能的，这就是为什么diagind比笛卡尔索引运行得更好。

- Oskar

Answer 4

Przemyslaw Szufel的解决方案在1_000 x 1_000矩阵大小的基准测试中表明，diagind表现最佳：

julia> @btime setindex!.(Ref($m), 0.0, 1:1_000, 1:1_000);
  2.222 μs (1 allocation: 7.94 KiB)

julia> @btime $m[diagind($m)] .= 0.0;
  1.280 μs (2 allocations: 80 bytes)

- Rafael_Guerra