如何在 Julia 中有效地计算二次型? [英] How can I efficiently calculate a quadratic form in Julia?

问题描述

我想计算一个二次型:x' Q y in Julia.
对于这些情况，最有效的计算方法是什么:

I would like to calculate a quadratic form: x' Q y in Julia.
What would be the most efficient way to calculate this for the cases:

1. 没有假设.
2. Q 是对称的.
3. x 和 y 相同(x = y).
4. Q 都是对称的，x = y.

1. No assumption.
2. Q is symmetric.
3. x and y are the same (x = y).
4. Both Q is symmetric and x = y.

我知道 Julia 有 dot().但我想知道它是否比 BLAS 调用更快.

I know Julia has dot(). But I wonder if it is faster than BLAS call.

推荐答案

Julia 的 LinearAlgebra 标准库具有 3 参数 dot 的本机实现，以及专门用于对称/厄米矩阵的版本.您可以查看源这里和这里.

Julia's LinearAlgebra stdlib has native implementations of 3-argument dot, and also a version that is specialized for symmetric/hermitian matrices. You can view the source here and here.

您可以使用 BenchmarkTools.@btime 或 BenchmarkTools.@ballocated 确认它们没有分配(记得使用 $ 插入变量).矩阵的对称性被利用了，但是查看源代码，我看不出 x == y 是如何实现任何严重的加速的，除了可能节省一些数组查找.

You can confirm that they do not allocate using BenchmarkTools.@btime or BenchmarkTools.@ballocated (remember to interpolate variables using $). The symmetry of the matrix is exploited, but looking at the source, I don't see how x == y could enable any serious speedup, except perhaps saving a few array lookups.

比较BLAS版本和原生的执行速度，可以做

To compare the execution speed of the BLAS version and the native one, you can do

1.7.0> using BenchmarkTools, LinearAlgebra

1.7.0> X = rand(100,100); y=rand(100);

1.7.0> @btime $y' * $X * $y
  42.800 μs (1 allocation: 896 bytes)
1213.5489200642382

1.7.0> @btime dot($y, $X, $y)
  1.540 μs (0 allocations: 0 bytes)
1213.548920064238

这是原生版本的一大胜利.但是，对于更大的矩阵，情况会发生变化:

This is a big win for the native version. For bigger matrices, the picture changes, though:

1.7.0> X = rand(10000,10000); y=rand(10000);

1.7.0> @btime $y' * $X * $y
  33.790 ms (2 allocations: 78.17 KiB)
1.2507105095988091e7

1.7.0> @btime dot($y, $X, $y)
  44.849 ms (0 allocations: 0 bytes)
1.2507105095988117e7

可能是因为 BLAS 使用线程，而 dot 不是多线程的.还有一些浮点差异.

Possibly because BLAS uses threads, while dot is not multithreaded. There are also some floating point differences.

这篇关于如何在 Julia 中有效地计算二次型?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！