Workspaces
Workspaces
represent the buffers and temporary storage that are used during the computations of LAPACK functions
. Upon initialization with a template matrix, work buffers will be allocated that are appropriate to be used during the factorization of matrices similar to the template, e.g. both Float64
and Float32
Matrices work, but also Complex
numbers are allowed when appropriate.
Workspace
The following convenience function is supplied in order to construct the correct Workspace
for a given LAPACK function
. This can then be used to perform the decompositions without extra allocations.
FastLapackInterface.Workspace
— TypeWorkspace(lapack_function, A)
Will create the correct Workspace
for the target lapack_function
and matrix A
.
Examples
julia> A = [1.2 2.3
6.2 3.3]
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
julia> ws = Workspace(LAPACK.geqrt!, A)
QRWYWs{Float64, Matrix{Float64}}
work: 4-element Vector{Float64}
T: 2×2 Matrix{Float64}
julia> LinearAlgebra.QRCompactWY(factorize!(ws, A)...)
LinearAlgebra.QRCompactWY{Float64, Matrix{Float64}, Matrix{Float64}}
Q factor: 2×2 LinearAlgebra.QRCompactWYQ{Float64, Matrix{Float64}, Matrix{Float64}}
R factor:
2×2 Matrix{Float64}:
-6.31506 -3.67692
0.0 -1.63102
Each Workspace
also has a function to resize!
to allow for its use with larger matrices or with more features (e.g. the computation of left eigenvectors and right eigenvectors using EigenWs
).
Base.resize!
— Methodresize!(ws, A; kwargs...)
Resizes the ws
to be appropriate for use with matrix A
. The kwargs
can be used to communicate which features should be supported by the Workspace
, such as left and right eigenvectors while using EigenWs
. This function is mainly used for automatic resizing inside LAPACK functions
.
QR
FastLapackInterface.QRWs
— TypeQRWs
Workspace for standard LinearAlgebra.QR
factorization using the LAPACK.geqrf!
function.
Examples
julia> A = [1.2 2.3
6.2 3.3]
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
julia> ws = QRWs(A)
QRWs{Float64}
work: 64-element Vector{Float64}
τ: 2-element Vector{Float64}
julia> t = QR(LAPACK.geqrf!(ws, A)...)
QR{Float64, Matrix{Float64}, Vector{Float64}}
Q factor: 2×2 LinearAlgebra.QRPackedQ{Float64, Matrix{Float64}, Vector{Float64}}
R factor:
2×2 Matrix{Float64}:
-6.31506 -3.67692
0.0 -1.63102
julia> Matrix(t)
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
FastLapackInterface.QRWYWs
— TypeQRWYWs
Workspace to be used with the LinearAlgebra.QRCompactWY
representation of the blocked QR factorization which uses the LAPACK.geqrt!
function. By default the blocksize for the algorithm is taken as min(36, min(size(template)))
, this can be overridden by using the blocksize
keyword of the constructor.
Examples
julia> A = [1.2 2.3
6.2 3.3]
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
julia> ws = QRWYWs(A)
QRWYWs{Float64, Matrix{Float64}}
work: 4-element Vector{Float64}
T: 2×2 Matrix{Float64}
julia> t = LinearAlgebra.QRCompactWY(LAPACK.geqrt!(ws, A)...)
LinearAlgebra.QRCompactWY{Float64, Matrix{Float64}, Matrix{Float64}}
Q factor: 2×2 LinearAlgebra.QRCompactWYQ{Float64, Matrix{Float64}, Matrix{Float64}}
R factor:
2×2 Matrix{Float64}:
-6.31506 -3.67692
0.0 -1.63102
julia> Matrix(t)
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
FastLapackInterface.QRPivotedWs
— TypeQRPivotedWs
Workspace to be used with the LinearAlgebra.QRPivoted
representation of the QR factorization which uses the LAPACK.geqp3!
function.
Examples
julia> A = [1.2 2.3
6.2 3.3]
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
julia> ws = QRPivotedWs(A)
QRPivotedWs{Float64, Float64}
work: 100-element Vector{Float64}
rwork: 0-element Vector{Float64}
τ: 2-element Vector{Float64}
jpvt: 2-element Vector{Int64}
julia> t = QRPivoted(LAPACK.geqp3!(ws, A)...)
QRPivoted{Float64, Matrix{Float64}, Vector{Float64}, Vector{Int64}}
Q factor: 2×2 LinearAlgebra.QRPackedQ{Float64, Matrix{Float64}, Vector{Float64}}
R factor:
2×2 Matrix{Float64}:
-6.31506 -3.67692
0.0 -1.63102
permutation:
2-element Vector{Int64}:
1
2
julia> Matrix(t)
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
FastLapackInterface.QROrmWs
— TypeQROrmWs
Workspace to be used with the LinearAlgebra.LAPACK.ormqr!
function. It requires the workspace of a QR
or a QRPivoted
previous factorization
Examples
```jldoctest julia> A = [1.2 2.3 6.2 3.3] 2×2 Matrix{Float64}: 1.2 2.3 6.2 3.3
julia> ws = QRPivotedWs(A) QRPivotedWs{Float64, Float64} work: 100-element Vector{Float64} rwork: 0-element Vector{Float64} τ: 2-element Vector{Float64} jpvt: 2-element Vector{Int64}
julia> C=[1 0.5; 2 1] 2×2 Matrix{Float64}: 1.0 0.5 2.0 1.0
julia> ormws = QROrmWs(ws, 'L', 'N', A, C) QROrmWs{Float64} work: 4224-element Vector{Float64} τ: 2-element Vector{Float64}
FastLapackInterface.QROrgWs
— TypeQROrgWs
Workspace to be used with the LinearAlgebra.LAPACK.orgqr!
function. It requires the workspace of a QR
or a QRPivoted
previous factorization
Examples
```jldoctest julia> A = [1.2 2.3 6.2 3.3] 2×2 Matrix{Float64}: 1.2 2.3 6.2 3.3
julia> ws = QRPivotedWs(A) QRPivotedWs{Float64, Float64} work: 100-element Vector{Float64} rwork: 0-element Vector{Float64} τ: 2-element Vector{Float64} jpvt: 2-element Vector{Int64}
julia> C=[1 0.5; 2 1] 2×2 Matrix{Float64}: 1.0 0.5 2.0 1.0
julia> orgws = QROrgWs(ws, 'L', 'N', A, C) QROrgWs{Float64} work: 4224-element Vector{Float64} τ: 2-element Vector{Float64}
Schur
FastLapackInterface.SchurWs
— TypeSchurWs
Workspace to be used with the LinearAlgebra.Schur
representation of the Schur decomposition which uses the LAPACK.gees!
function.
Examples
julia> A = [1.2 2.3
6.2 3.3]
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
julia> ws = SchurWs(A)
SchurWs{Float64}
work: 68-element Vector{Float64}
wr: 2-element Vector{Float64}
wi: 2-element Vector{Float64}
vs: 2×2 Matrix{Float64}
sdim: Base.RefValue{Int64}
bwork: 2-element Vector{Int64}
eigen_values: 2-element Vector{ComplexF64}
julia> t = Schur(LAPACK.gees!(ws, 'V', A)...)
Schur{Float64, Matrix{Float64}, Vector{Float64}}
T factor:
2×2 Matrix{Float64}:
-1.6695 -3.9
0.0 6.1695
Z factor:
2×2 Matrix{Float64}:
-0.625424 -0.780285
0.780285 -0.625424
eigenvalues:
2-element Vector{Float64}:
-1.6695025194532018
6.169502519453203
julia> Matrix(t)
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
FastLapackInterface.GeneralizedSchurWs
— TypeGeneralizedSchurWs
Workspace to be used with the LinearAlgebra.GeneralizedSchur
representation of the Generalized Schur decomposition which uses the LAPACK.gges!
function.
Examples
julia> A = [1.2 2.3
6.2 3.3]
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
julia> B = [8.2 0.3
1.7 4.3]
2×2 Matrix{Float64}:
8.2 0.3
1.7 4.3
julia> ws = GeneralizedSchurWs(A)
GeneralizedSchurWs{Float64}
work: 90-element Vector{Float64}
αr: 2-element Vector{Float64}
αi: 2-element Vector{Float64}
β: 2-element Vector{Float64}
vsl: 2×2 Matrix{Float64}
vsr: 2×2 Matrix{Float64}
sdim: Base.RefValue{Int64}
bwork: 2-element Vector{Int64}
eigen_values: 2-element Vector{ComplexF64}
julia> t = GeneralizedSchur(LAPACK.gges!(ws, 'V','V', A, B)...)
GeneralizedSchur{Float64, Matrix{Float64}, Vector{ComplexF64}, Vector{Float64}}
S factor:
2×2 Matrix{Float64}:
-1.43796 1.63843
0.0 7.16295
T factor:
2×2 Matrix{Float64}:
5.06887 -4.00221
0.0 6.85558
Q factor:
2×2 Matrix{Float64}:
-0.857329 0.514769
0.514769 0.857329
Z factor:
2×2 Matrix{Float64}:
-0.560266 0.828313
0.828313 0.560266
α:
2-element Vector{ComplexF64}:
-1.4379554610733563 + 0.0im
7.162947865097022 + 0.0im
β:
2-element Vector{Float64}:
5.068865029631368
6.855578082442485
LU
FastLapackInterface.LUWs
— TypeLUWs
Workspace to be used with the LinearAlgebra.LU
representation of the LU factorization which uses the LAPACK.getrf!
function.
Examples
julia> A = [1.2 2.3
6.2 3.3]
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
julia> ws = LUWs(A)
LUWs
ipiv: 2-element Vector{Int64}
julia> t = LU(LAPACK.getrf!(ws, A)...)
LU{Float64, Matrix{Float64}, Vector{Int64}}
L factor:
2×2 Matrix{Float64}:
1.0 0.0
0.193548 1.0
U factor:
2×2 Matrix{Float64}:
6.2 3.3
0.0 1.66129
Eigen
FastLapackInterface.EigenWs
— TypeEigenWs
Workspace for LinearAlgebra.Eigen
factorization using the LAPACK.geevx!
function.
Examples
julia> A = [1.2 2.3
6.2 3.3]
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
julia> ws = EigenWs(A, rvecs=true)
EigenWs{Float64, Matrix{Float64}, Float64}
work: 260-element Vector{Float64}
rwork: 2-element Vector{Float64}
VL: 0×2 Matrix{Float64}
VR: 2×2 Matrix{Float64}
W: 2-element Vector{Float64}
scale: 2-element Vector{Float64}
iwork: 0-element Vector{Int64}
rconde: 0-element Vector{Float64}
rcondv: 0-element Vector{Float64}
julia> t = LAPACK.geevx!(ws, 'N', 'N', 'V', 'N', A);
julia> LinearAlgebra.Eigen(t[2], t[5])
Eigen{Float64, Float64, Matrix{Float64}, Vector{Float64}}
values:
2-element Vector{Float64}:
-1.6695025194532018
6.169502519453203
vectors:
2×2 Matrix{Float64}:
-0.625424 -0.420019
0.780285 -0.907515
FastLapackInterface.HermitianEigenWs
— TypeHermitianEigenWs
Workspace to be used with Hermitian diagonalization using the LAPACK.syevr!
function. Supports both Real
and Complex
Hermitian matrices.
Examples
julia> A = [1.2 2.3
6.2 3.3]
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
julia> ws = HermitianEigenWs(A, vecs=true)
HermitianEigenWs{Float64, Matrix{Float64}, Float64}
work: 66-element Vector{Float64}
rwork: 0-element Vector{Float64}
iwork: 20-element Vector{Int64}
w: 2-element Vector{Float64}
Z: 2×2 Matrix{Float64}
isuppz: 4-element Vector{Int64}
julia> LinearAlgebra.Eigen(LAPACK.syevr!(ws, 'V', 'A', 'U', A, 0.0, 0.0, 0, 0, 1e-6)...)
Eigen{Float64, Float64, Matrix{Float64}, Vector{Float64}}
values:
2-element Vector{Float64}:
-0.2783393759541063
4.778339375954106
vectors:
2×2 Matrix{Float64}:
-0.841217 0.540698
0.540698 0.841217
FastLapackInterface.GeneralizedEigenWs
— TypeGeneralizedEigenWs
Workspace that can be used for LinearAlgebra.GeneralizedEigen
factorization using LAPACK.ggev!
. Supports Real
and Complex
matrices.
Examples
julia> A = [1.2 2.3
6.2 3.3]
2×2 Matrix{Float64}:
1.2 2.3
6.2 3.3
julia> B = [8.2 1.7
5.9 2.1]
2×2 Matrix{Float64}:
8.2 1.7
5.9 2.1
julia> ws = GeneralizedEigenWs(A, rvecs=true)
GeneralizedEigenWs{Float64, Matrix{Float64}, Float64}
work: 78-element Vector{Float64}
vl: 0×2 Matrix{Float64}
vr: 2×2 Matrix{Float64}
αr: 2-element Vector{Float64}
αi: 2-element Vector{Float64}
β: 2-element Vector{Float64}
julia> αr, αi, β, _, vr = LAPACK.ggev!(ws, 'N', 'V', A, B);
julia> LinearAlgebra.GeneralizedEigen(αr ./ β, vr)
GeneralizedEigen{Float64, Float64, Matrix{Float64}, Vector{Float64}}
values:
2-element Vector{Float64}:
-0.8754932558185097
1.6362721153456299
vectors:
2×2 Matrix{Float64}:
-0.452121 -0.0394242
1.0 1.0
BunchKaufman
FastLapackInterface.BunchKaufmanWs
— TypeBunchKaufmanWs
Workspace for LinearAlgebra.BunchKaufman
factorization using the LAPACK.sytrf!
or LAPACK.sytrf_rook!
functions for symmetric matrices, and LAPACK.hetrf!
or LAPACK.hetrf_rook!
functions for hermitian matrices (e.g. with ComplexF64
or ComplexF32
elements).
Examples
julia> A = [1.2 7.8
7.8 3.3]
2×2 Matrix{Float64}:
1.2 7.8
7.8 3.3
julia> ws = BunchKaufmanWs(A)
BunchKaufmanWs{Float64}
work: 128-element Vector{Float64}
ipiv: 2-element Vector{Int64}
julia> A, ipiv, info = LAPACK.sytrf!(ws, 'U', A)
([1.2 7.8; 7.8 3.3], [-1, -1], 0)
julia> t = LinearAlgebra.BunchKaufman(A, ipiv,'U', true, false, info)
BunchKaufman{Float64, Matrix{Float64}, Vector{Int64}}
D factor:
2×2 Tridiagonal{Float64, Vector{Float64}}:
1.2 7.8
7.8 3.3
U factor:
2×2 UnitUpperTriangular{Float64, Matrix{Float64}}:
1.0 0.0
⋅ 1.0
permutation:
2-element Vector{Int64}:
1
2
julia> A = [1.2 7.8
7.8 3.3]
2×2 Matrix{Float64}:
1.2 7.8
7.8 3.3
julia> ws = BunchKaufmanWs(A)
BunchKaufmanWs{Float64}
work: 128-element Vector{Float64}
ipiv: 2-element Vector{Int64}
julia> A, ipiv, info = LAPACK.sytrf_rook!(ws, 'U', A)
([1.2 7.8; 7.8 3.3], [-1, -2], 0)
julia> t = LinearAlgebra.BunchKaufman(A, ipiv,'U', true, true, info)
BunchKaufman{Float64, Matrix{Float64}, Vector{Int64}}
D factor:
2×2 Tridiagonal{Float64, Vector{Float64}}:
1.2 7.8
7.8 3.3
U factor:
2×2 UnitUpperTriangular{Float64, Matrix{Float64}}:
1.0 0.0
⋅ 1.0
permutation:
2-element Vector{Int64}:
1
2
Cholesky
FastLapackInterface.CholeskyPivotedWs
— TypeCholeskyPivotedWs
Workspace for LinearAlgebra.CholeskyPivoted
factorization using the LAPACK.pstrf!
function. The standard LinearAlgebra.Cholesky
uses LAPACK.potrf!
which is non-allocating and does not require a separate Workspace
.
Examples
julia> A = [1.2 7.8
7.8 3.3]
2×2 Matrix{Float64}:
1.2 7.8
7.8 3.3
julia> ws = CholeskyPivotedWs(A)
CholeskyPivotedWs{Float64}
work: 4-element Vector{Float64}
piv: 2-element Vector{Int64}
julia> AA, piv, rank, info = LAPACK.pstrf!(ws, 'U', A, 1e-6)
([1.816590212458495 4.293758683992806; 7.8 -17.236363636363635], [2, 1], 1, 1)
julia> CholeskyPivoted(AA, 'U', piv, rank, 1e-6, info)
CholeskyPivoted{Float64, Matrix{Float64}, Vector{Int64}}
U factor with rank 1:
2×2 UpperTriangular{Float64, Matrix{Float64}}:
1.81659 4.29376
⋅ -17.2364
permutation:
2-element Vector{Int64}:
2
1
LSE
FastLapackInterface.LSEWs
— TypeLSEWs
Workspace for the least squares solving function LAPACK.geqrf!
.
Examples
julia> A = [1.2 2.3 6.2
6.2 3.3 8.8
9.1 2.1 5.5]
3×3 Matrix{Float64}:
1.2 2.3 6.2
6.2 3.3 8.8
9.1 2.1 5.5
julia> B = [2.7 3.1 7.7
4.1 8.1 1.8]
2×3 Matrix{Float64}:
2.7 3.1 7.7
4.1 8.1 1.8
julia> c = [0.2, 7.2, 2.9]
3-element Vector{Float64}:
0.2
7.2
2.9
julia> d = [3.9, 2.1]
2-element Vector{Float64}:
3.9
2.1
julia> ws = LSEWs(A, B)
LSEWs{Float64}
work: 101-element Vector{Float64}
X: 3-element Vector{Float64}
julia> LAPACK.gglse!(ws, A, c, B, d)
([0.19723156207005318, 0.0683561362406917, 0.40981438442398854], 13.750943845251626)