对于小型固定大小矩阵，特征值要快多少

Question

我现在正在使用 Julia，但我有一个性能关键函数，它需要在小型固定大小矩阵（3 维或 4 维）上进行大量重复矩阵运算。似乎Julia中的所有矩阵操作都由BLAS和LAPACK后端处理。在其中一些函数中似乎也发生了大量的内存分配。

有一个用于小型矩阵的 julia 库，该库拥有令人印象深刻的 3x3 矩阵加速，但它已经 3 年没有更新了。我正在考虑在本征中重写我的性能关键函数

我知道 Eigen

声称非常适合固定大小的矩阵，但我仍在尝试判断我是否应该用 Eigen 重写这个函数。性能基准适用于动态大小矩阵。有没有人有任何数据表明从固定大小的矩阵中获得多少性能？我正在做的运算类型是矩阵 x 矩阵、矩阵 x 向量、正定线性求解。

Answer 1

如果你想对小矩阵进行快速操作，我强烈推荐StaticArrays。例如（注意：这最初是在 BenchmarkTools 包之前编写的，现在推荐）：

using StaticArrays
using LinearAlgebra
function foo(A, b, n)
    s = 0.0
    for i = 1:n
        s += sum(A*b)
    end
    s
end
function foo2(A, b, n)
    c = A*b
    s = 0.0
    for i = 1:n
        mul!(c, A, b)
        s += sum(c)
    end
    s
end
A = rand(3,3)
b = rand(3)
Af = SMatrix{3,3}(A)
bf = SVector{3}(b)
foo(A, b, 1)
foo2(A, b, 1)
foo(Af, bf, 1)
@time foo(A, b, 10^6)
@time foo2(A, b, 10^6)
@time foo(Af, bf, 10^6)

结果：

julia> include("/tmp/foo.jl")
  0.080535 seconds (1.00 M allocations: 106.812 MiB, 14.86% gc time)
  0.064963 seconds (3 allocations: 144 bytes)
  0.001719 seconds (2 allocations: 32 bytes)

foo2试图变得聪明并避免内存分配，但它只是被使用StaticArrays时的幼稚实现所震撼。