Non-allocating matrix inversion

General Usage Performance
8 
using LinearAlgebra, BenchmarkTools, BenchmarkPlots, StatsPlots
bg = BenchmarkGroup()
bg[:inv] = @benchmarkable A*inv(B)*b setup=(A=rand(100,100);B=rand(100,100);b=rand(100))
bg[:backslash] = @benchmarkable A*(B\b) setup=(A=rand(100,100);B=rand(100,100);b=rand(100))
bg[:slash] = @benchmarkable (A/B)*b setup=(A=rand(100,100);B=rand(100,100);b=rand(100))
bg[:qr] = @benchmarkable A*(qr(B)\b) setup=(A=rand(100,100);B=rand(100,100);b=rand(100))
bg[:lu] = @benchmarkable A*(lu(B)\b) setup=(A=rand(100,100);B=rand(100,100);b=rand(100))
bg[:inv_pre] = @benchmarkable A*Q*b setup=(A=rand(100,100);B=rand(100,100);Q=inv(B);b=rand(100))
bg[:qr_pre] = @benchmarkable A*(Q\b) setup=(A=rand(100,100);B=rand(100,100);Q=qr(B);b=rand(100))
bg[:lu_pre] = @benchmarkable A*(Q\b) setup=(A=rand(100,100);B=rand(100,100);Q=lu(B);b=rand(100))
trial = run(bg)
plot(trial; yscale=:log10)

image

Even ignoring the numerical instability of the inversion process, you can see that A*(B\b) is faster than inversion, and that qr(B) is faster when used as a precomputed auxiliary.

Edit: Added in LU decomposition, because it’s the preferred algorithm.

6 Likes