Sparse Linear Solvers

CHOLMOD.Factor{Tv,Ti} <: Factorization{Tv}

The Cholesky (LL') or LDL' factorization of a sparse symmetric or Hermitian matrix, returned by cholesky and ldlt. It wraps a pointer to a cholmod_factor struct, which holds the symbolic analysis, the fill-reducing permutation and the numeric factor. Tv is Float64, Float32, ComplexF64 or ComplexF32, and Ti is Int32 or Int64 (only Int32 on 32-bit systems).

With P the permutation matrix of F.p, the factorization is A == P'*L*L'*P for cholesky and A == P'*L*D*L'*P for ldlt. The properties of F are:

Apart from F.p, these are lazy components for use with \, each applying one step of a solve. Only sparse(F.L), and for ldlt sparse(F.LD) (the unit L with D on its diagonal), can be materialized. sparse(F) reconstructs the factorized matrix.

F supports \, ldiv!, det, logdet, logabsdet, diag, issuccess, nnz, copy, CHOLMOD.rcond, refactorization with cholesky! and ldlt!, and the low-rank modifications lowrankdowndate and lowrankupdate. copy(F) and deepcopy(F) return an independent factorization: nothing done to one affects the other.

CHOLMOD owns the memory, which is released by a finalizer. The pointer is null after deserialization, and using such a factorization throws an ArgumentError. ldiv! takes an optional CHOLMOD.CholmodWS to avoid allocating. A solve (\ or ldiv!) only reads F and takes no lock, so several tasks can solve with one F at the same time, provided that none of them refactorizes or updates it meanwhile; those calls take a lock internal to F.

Examples

julia> F = ldlt(sparse([4.0 2.0; 2.0 -3.0]));

julia> propertynames(F)
(:L, :U, :PtL, :UP, :D, :LD, :DU, :PtLD, :DUP, :p, :ptr)

julia> F \ [6.0, -1.0] ≈ [1.0, 1.0]
true

CHOLMOD.Sparse{Tv,Ti} <: AbstractSparseMatrix{Tv,Ti}

A sparse matrix in compressed column form stored in memory allocated by CHOLMOD, wrapping a pointer to a cholmod_sparse struct. Tv is Float64, Float32, ComplexF64 or ComplexF32, and Ti is Int32 or Int64 (only Int32 on 32-bit systems).

Sparse(A) copies a SparseMatrixCSC, or a Symmetric or Hermitian view of one, converting element types CHOLMOD does not support to a floating-point type it does. CHOLMOD records in the struct's stype whether the whole matrix or only one triangle of a symmetric or Hermitian matrix is stored; a plain SparseMatrixCSC that is Hermitian is stored as one triangle. sparse(S) copies back to a SparseMatrixCSC, wrapped in Symmetric or Hermitian when only a triangle is stored, so it is not type stable. Sparse(F) returns the L (or LD) factor of a CHOLMOD.Factor F.

The memory is released by a finalizer. The pointer is null after deserialization, and using such an object throws an ArgumentError.

CHOLMOD.Dense{Tv} <: DenseMatrix{Tv}

A dense matrix stored in memory allocated by CHOLMOD, wrapping a pointer to a cholmod_dense struct. Tv is Float64, Float32, ComplexF64 or ComplexF32. Dense(A) copies a strided vector or matrix A into CHOLMOD storage, and Matrix, Vector and copyto! copy it back into a Julia array. It is the dense right-hand side and solution type of the CHOLMOD solve routines; F \ b with a Julia array b does not need it.

The memory is released by a finalizer. The pointer is null after deserialization, and using such an object throws an ArgumentError.

cholesky(A::SparseMatrixCSC; shift = 0.0, check = true, perm = nothing) -> CHOLMOD.Factor

Compute the Cholesky factorization of a sparse positive definite matrix A. A must be a SparseMatrixCSC or a Symmetric/Hermitian view of a SparseMatrixCSC. Note that if A doesn't have the type tag, it must itself be symmetric or Hermitian. If perm is not given, a fill-reducing permutation is used. F = cholesky(A) is most frequently used to solve systems of equations with F\b, but also the methods diag, det, and logdet are defined for F. You can also extract individual factors from F, using F.L. However, since pivoting is on by default, the factorization is internally represented as A == P'*L*L'*P with a permutation matrix P; using just L without accounting for P will give incorrect answers. To include the effects of permutation, it's typically preferable to extract "combined" factors like PtL = F.PtL (the equivalent of P'*L) and LtP = F.UP (the equivalent of L'*P). The complete list of supported factors is :L, :PtL, :UP, :U. The permutation vector is available as F.p, defined such that L*L' == A[p, p],

The L component can be materialized as a sparse matrix using sparse(F.L). Other components cannot be materialized directly, but can be reconstructed from sparse(F.L) and F.p if needed.

When check = true, an error is thrown if the decomposition fails. When check = false, responsibility for checking the decomposition's validity (via issuccess) lies with the user.

Setting the optional shift keyword argument computes the factorization of A+shift*I instead of A. If the perm argument is provided, it should be a permutation of 1:size(A,1) giving the ordering to use (instead of CHOLMOD's default AMD ordering).

See also ldlt for a similar factorization that does not require positive definiteness, but can be significantly slower than cholesky.

Examples

In the following example, the fill-reducing permutation used is [3, 2, 1]. If perm is set to 1:3 to enforce no permutation, the number of nonzero elements in the factor is 6.

julia> A = [2 1 1; 1 2 0; 1 0 2]
3×3 Matrix{Int64}:
 2  1  1
 1  2  0
 1  0  2

julia> C = cholesky(sparse(A))
SparseArrays.CHOLMOD.Factor{Float64, Int64}
type:    LLt
method:  simplicial
maxnnz:  5
nnz:     5
success: true

julia> C.p
3-element Vector{Int64}:
 3
 2
 1

julia> L = sparse(C.L);

julia> Matrix(L)
3×3 Matrix{Float64}:
 1.41421   0.0       0.0
 0.0       1.41421   0.0
 0.707107  0.707107  1.0

julia> L * L' ≈ A[C.p, C.p]
true

julia> P = sparse(1:3, C.p, ones(3))
3×3 SparseMatrixCSC{Float64, Int64} with 3 stored entries:
  ⋅    ⋅   1.0
  ⋅   1.0   ⋅
 1.0   ⋅    ⋅

julia> P' * L * L' * P ≈ A
true

julia> C = cholesky(sparse(A), perm=1:3)
SparseArrays.CHOLMOD.Factor{Float64, Int64}
type:    LLt
method:  simplicial
maxnnz:  6
nnz:     6
success: true

julia> L = sparse(C.L);

julia> Matrix(L)
3×3 Matrix{Float64}:
 1.41421    0.0       0.0
 0.707107   1.22474   0.0
 0.707107  -0.408248  1.1547

julia> L * L' ≈ A
true

cholesky!(F::CHOLMOD.Factor, A::SparseMatrixCSC; shift = 0.0, check = true) -> CHOLMOD.Factor

Compute the Cholesky ($LL'$) factorization of A, reusing the symbolic factorization F. A must be a SparseMatrixCSC or a Symmetric/ Hermitian view of a SparseMatrixCSC. Note that if A doesn't have the type tag, it must itself be symmetric or Hermitian.

A must have the sparsity pattern of the matrix that F was computed from, since only the numerical factorization is redone.

See also cholesky.

Examples

julia> A = sparse([4.0 1 0; 1 4 1; 0 1 4]); b = [1.0, 2.0, 3.0];

julia> F = cholesky(A);

julia> cholesky!(F, 2A);

julia> F \ b ≈ (2A) \ b
true

ldlt(A::SparseMatrixCSC; shift = 0.0, check = true, perm=nothing) -> CHOLMOD.Factor

Compute the $LDL'$ factorization of a sparse matrix A. A must be a SparseMatrixCSC or a Symmetric/Hermitian view of a SparseMatrixCSC. Note that if A doesn't have the type tag, it must itself be symmetric or Hermitian. A fill-reducing permutation is used. F = ldlt(A) is most frequently used to solve systems of equations A*x = b with F\b. The returned factorization object F also supports the methods diag, det, logdet, and inv. You can extract individual factors from F using F.L. However, since pivoting is on by default, the factorization is internally represented as A == P'*L*D*L'*P with a permutation matrix P; using just L without accounting for P will give incorrect answers. To include the effects of permutation, it is typically preferable to extract "combined" factors like PtL = F.PtL (the equivalent of P'*L) and LtP = F.UP (the equivalent of L'*P). The complete list of supported factors is :L, :PtL, :D, :UP, :U, :LD, :DU, :PtLD, :DUP. Each one acts as the matrix its name spells out, so that for instance F.PtL \ b solves with P'*L and F.LD \ b solves with the product L*D. The permutation vector is available as F.p, defined such that L*D*L' == A[p, p].

Of these, only LD can be materialized, with sparse(F.LD). Beware that the matrix it returns is not the product L*D: it is CHOLMOD's packed $LDL'$ factor, which stores L with its unit diagonal overwritten by the diagonal of D. Solving with it is therefore not the same as solving with F.LD. Unpack it as

julia> A = sparse([4.0 1 0; 1 -4 1; 0 1 4]); F = ldlt(A);

julia> LD = sparse(F.LD);

julia> D = Diagonal(diag(LD));  # equivalently, Diagonal(diag(F))

julia> L = tril(LD, -1) + I;    # unit lower triangular

julia> L * D * L' ≈ A[F.p, F.p]
true

so that L*D*L' == A[F.p, F.p]. The remaining components cannot be materialized directly, but can be reconstructed from L, D and F.p.

Unlike the related Cholesky factorization, the $LDL'$ factorization does not require A to be positive definite. However, it still requires all leading principal minors to be well-conditioned and will fail if this is not satisfied.

When check = true, an error is thrown if the decomposition fails. When check = false, responsibility for checking the decomposition's validity (via issuccess) lies with the user.

Setting the optional shift keyword argument computes the factorization of A+shift*I instead of A. If the perm argument is provided, it should be a permutation of 1:size(A,1) giving the ordering to use (instead of CHOLMOD's default AMD ordering).

See also cholesky for a factorization that can be significantly faster than ldlt, but requires A to be positive definite.

ldlt!(F::CHOLMOD.Factor, A::SparseMatrixCSC; shift = 0.0, check = true) -> CHOLMOD.Factor

Compute the $LDL'$ factorization of A, reusing the symbolic factorization F. A must be a SparseMatrixCSC or a Symmetric/Hermitian view of a SparseMatrixCSC. Note that if A doesn't have the type tag, it must itself be symmetric or Hermitian.

A must have the sparsity pattern of the matrix that F was computed from, since only the numerical factorization is redone.

See also ldlt.

Examples

julia> A = sparse([4.0 1 0; 1 -4 1; 0 1 4]); b = [1.0, 2.0, 3.0];

julia> F = ldlt(A);

julia> ldlt!(F, 2A);

julia> F \ b ≈ (2A) \ b
true

lowrankupdate(F::CHOLMOD.Factor, C::AbstractArray) -> FF::CHOLMOD.Factor

Get an LDLt Factorization of A + C*C' given an LDLt or LLt factorization F of A.

The returned factor is always an LDLt factorization.

Only real factorizations are supported; CHOLMOD cannot update or downdate a complex one.

See also lowrankupdate!, lowrankdowndate, lowrankdowndate!.

Examples

julia> A = sparse([4.0 1 0; 1 4 1; 0 1 4]); c = [1.0, 0.0, 2.0]; b = [1.0, 2.0, 3.0];

julia> F = cholesky(A);

julia> G = lowrankupdate(F, c);

julia> G \ b ≈ (A + c * c') \ b, lowrankdowndate(G, c) \ b ≈ A \ b
(true, true)

lowrankupdate!(F::CHOLMOD.Factor, C::AbstractArray)

Update an LDLt or LLt Factorization F of A to a factorization of A + C*C'.

LLt factorizations are converted to LDLt.

Only real factorizations are supported; CHOLMOD cannot update or downdate a complex one.

See also lowrankupdate, lowrankdowndate, lowrankdowndate!.

Examples

julia> A = sparse([4.0 1 0; 1 4 1; 0 1 4]); c = [1.0, 0.0, 2.0]; b = [1.0, 2.0, 3.0];

julia> F = cholesky(A);

julia> lowrankupdate!(F, c);

julia> F \ b ≈ (A + c * c') \ b
true

lowrankdowndate(F::CHOLMOD.Factor, C::AbstractArray) -> FF::CHOLMOD.Factor

Get an LDLt Factorization of A - C*C' given an LDLt or LLt factorization F of A.

The returned factor is always an LDLt factorization.

Only real factorizations are supported; CHOLMOD cannot update or downdate a complex one.

See also lowrankdowndate!, lowrankupdate, lowrankupdate!.

lowrankdowndate!(F::CHOLMOD.Factor, C::AbstractArray)

Update an LDLt or LLt Factorization F of A to a factorization of A - C*C'.

LLt factorizations are converted to LDLt.

Only real factorizations are supported; CHOLMOD cannot update or downdate a complex one.

See also lowrankdowndate, lowrankupdate, lowrankupdate!.

rcond(F::CHOLMOD.Factor) -> Float64

Return CHOLMOD's rough estimate of the reciprocal condition number of the factorized matrix, computed from the diagonal of the factor alone: the smallest entry of abs.(diag(F)) divided by the largest, squared when F is an LL' factorization so that the result estimates the reciprocal condition number of the factorized matrix rather than of its factor.

This is much cheaper than a norm-based estimate such as cond(A, 1), but also much cruder. For positive definite A it is exact when A is diagonal, and otherwise an upper bound on 1 / cond(A, 2), so it can report a matrix as far better conditioned than it is. Use it to detect a badly conditioned or singular factorization, not to measure conditioning accurately. The LU counterpart is UMFPACK.rcond.

Returns 0 if the matrix is singular or the factor has a zero or NaN on its diagonal, and 1 if the matrix is 1-by-1. NaN is never returned.

Examples

julia> A = sparse(Diagonal([1.0, 2.0, 4.0]));

julia> SparseArrays.CHOLMOD.rcond(cholesky(A))
0.25

julia> SparseArrays.CHOLMOD.rcond(ldlt(A))
0.25

julia> SparseArrays.CHOLMOD.rcond(cholesky(sparse(Diagonal([1.0, 0.0])); check=false))
0.0

CHOLMOD.CholmodWS(F::CHOLMOD.Factor)

Scratch space for ldiv!(x, F, b; workspace), which makes repeated solves allocation-free. Without it, ldiv! allocates its scratch space on each call. A workspace can be reused with any Factor of the same index type, but not by two calls at once. Its memory is released by a finalizer or by CHOLMOD.free!.

UMFPACK.UmfpackLU{Tv,Ti} <: Factorization{Tv}

The LU factorization of a sparse matrix computed by UMFPACK, returned by lu. Tv is Float64 or ComplexF64, and Ti is Int32 or Int64 (only Int32 on 32-bit systems). F holds a zero-based copy of the factorized matrix together with UMFPACK's opaque symbolic and numeric objects, which are released by finalizers.

The factors are copied out of UMFPACK on each property access:

They satisfy F.L * F.U == (F.Rs .* A)[F.p, F.q].

F supports \, ldiv!, det, logabsdet, issuccess, nnz, adjoint, transpose, UMFPACK.rcond and refactorization with lu!. A serialized UmfpackLU carries the matrix rather than the factors, which are recomputed on first use after deserialization.

ldiv! takes an optional UMFPACK.UmfpackWS to avoid allocating. Solves (\ or ldiv!, with F, F' or transpose(F)), det, logabsdet and the factor properties only read F and take no lock, so several tasks can use one F at the same time, provided that none of them refactorizes it with lu! meanwhile; a task that does needs its own copy(F), which is independent of F and duplicates the factors. lu! and copy take an internal lock.

Examples

julia> A = sparse([4.0 1.0 0.0; 1.0 4.0 1.0; 0.0 1.0 4.0]);

julia> F = lu(A);

julia> F.L * F.U ≈ (F.Rs .* A)[F.p, F.q]
true

julia> L, U, p, q, Rs = F.:(:);

julia> L * U ≈ (Rs .* A)[p, q]
true

lu(A::AbstractSparseMatrixCSC; check = true, q = nothing, control = get_umfpack_control()) -> F::UmfpackLU

Compute the LU factorization of a sparse matrix A.

For sparse A with real or complex element type, the return type of F is UmfpackLU{Tv, Ti}, with Tv = Float64 or ComplexF64 respectively and Ti is an integer type (Int32 or Int64).

When check = true, an error is thrown if the decomposition fails. When check = false, responsibility for checking the decomposition's validity (via issuccess) lies with the user.

The column permutation q can either be a permutation vector or nothing. If no permutation vector is provided or q is nothing, UMFPACK's default is used. If the permutation is not zero-based, a zero-based copy is made.

The control vector defaults to the Julia SparseArrays package's default configuration for UMFPACK (NB: this is modified from the UMFPACK defaults to disable iterative refinement), but can be changed by passing a vector of length UMFPACK_CONTROL, see the UMFPACK manual for possible configurations. For example to reenable iterative refinement:

control = SparseArrays.UMFPACK.get_umfpack_control(Float64, Int64)  # the defaults for a Float64 matrix
control[SparseArrays.UMFPACK.JL_UMFPACK_IRSTEP] = 2.0  # UMFPACK's default number of refinement steps
F = lu(A; control)
x = F \ b   # solves with iterative refinement

The individual components of the factorization F can be accessed by indexing:

The relation between F and A is

F.L*F.U == (F.Rs .* A)[F.p, F.q]

F further supports the following functions:

• \

• det

See also lu!

Examples

julia> A = sparse([2.0 1 0; 0 3 1; 1 0 4]); b = [1.0, 2.0, 3.0];

julia> F = lu(A);

julia> F \ b ≈ A \ b
true

julia> F.L * F.U ≈ (F.Rs .* A)[F.p, F.q]
true

lu!(F::UmfpackLU, A::AbstractSparseMatrixCSC; check=true, reuse_symbolic=true, q=nothing) -> F::UmfpackLU

Compute the LU factorization of a sparse matrix A, reusing the symbolic factorization of an already existing LU factorization stored in F. Unless reuse_symbolic is set to false, the sparse matrix A must have an identical nonzero pattern as the matrix used to create the LU factorization F, otherwise an error is thrown. If the size of A and F differ, all vectors are resized accordingly and the symbolic factorization is computed anew.

When check = true, an error is thrown if the decomposition fails. When check = false, responsibility for checking the decomposition's validity (via issuccess) lies with the user.

The column permutation q can either be a permutation vector or nothing. If no permutation vector is provided or q is nothing, UMFPACK's default is used. If the permutation is not zero based, a zero based copy is made.

See also lu

Examples

julia> A = sparse(Float64[1.0 2.0; 0.0 3.0]);

julia> F = lu(A);

julia> B = sparse(Float64[1.0 1.0; 0.0 1.0]);

julia> lu!(F, B);

julia> F \ ones(2)
2-element Vector{Float64}:
 0.0
 1.0

rcond(F::UmfpackLU) -> Float64

Return UMFPACK's rough estimate of the reciprocal condition number of the factorized matrix, computed from the diagonal of the factor alone: the smallest entry of abs.(diag(F.U)) divided by the largest.

This is much cheaper than a norm-based estimate such as cond(A, 1), but also much cruder, and it describes the matrix UMFPACK actually factorized rather than A itself. UMFPACK scales the rows of A before factorizing by default (see F.Rs), so for instance every diagonal matrix reports 1. Unlike the Cholesky-based CHOLMOD.rcond, the value is neither an upper nor a lower bound on 1 / cond(A, 2). Use it to detect a singular or badly pivoted factorization, not to measure conditioning.

Returns 0 if the matrix is singular, and 1 if the matrix is 1-by-1.

Examples

julia> F = lu(sparse([1.0 3.0; 0.0 1.0]));

julia> SparseArrays.UMFPACK.rcond(F)
0.25

julia> minimum(abs, diag(F.U)) / maximum(abs, diag(F.U))
0.25

julia> SparseArrays.UMFPACK.rcond(lu(sparse([1.0 2.0; 0.0 0.0]); check=false))
0.0

UMFPACK.UmfpackWS(F::UmfpackLU)

Scratch space for ldiv!(x, F, b; workspace), which makes repeated solves allocation-free. Without it, ldiv! allocates its scratch space on each call. A workspace grows as needed, so it can be reused with any factorization, but not by two calls at once.

A solve leaves UMFPACK's statistics of it, such as the number of refinement steps taken, in workspace.info, which UMFPACK.show_umf_info(F, workspace) prints. They are those of the last system solved: the last column of a matrix right-hand side, and the imaginary part of a complex right-hand side with a real F.

copy(F::UmfpackLU)::UmfpackLU

Return an independent copy of F, with its own matrix, symbolic and numeric factors, control, info and lock; refactorizing either one with lu! does not affect the other. deepcopy(F) does the same.

SPQR.QRSparse{Tv,Ti} <: Factorization{Tv}

The QR factorization of a sparse matrix computed by SPQR, returned by qr. SPQR's output is copied into Julia arrays, so F holds no memory owned by the C library. Ti is Int32 or Int64 (only Int32 on 32-bit systems); see qr for the element types.

They satisfy F.Q * F.R == A[F.prow, F.pcol], where F.Q is square and only its leading columns enter the product.

F supports \ and ldiv! for least squares and minimum-norm solutions, rank, copy, and F', which is the LQ factorization AdjointQRSparse of A'. ldiv!, with F or F', takes an optional SPQR.SpqrWS to avoid allocating. F is never modified once qr has returned and a solve only reads it, so several tasks can solve with one F at the same time.

Examples

julia> A = sparse([1.0 0.0; 1.0 1.0; 0.0 1.0]);

julia> F = qr(A);

julia> propertynames(F)
(:R, :Q, :prow, :pcol)

julia> F.Q * F.R ≈ A[F.prow, F.pcol]
true

qr(A::SparseMatrixCSC; tol=_default_tol(A), ordering=ORDERING_DEFAULT) -> QRSparse

Compute the QR factorization of a sparse matrix A. Fill-reducing row and column permutations are used such that F.R = F.Q'*A[F.prow,F.pcol]. The main application of this type is to solve least squares or underdetermined problems with \. The function calls the C library SPQR^[ACM933].

• tol: a column whose norm falls to tol or below during the factorization is treated as zero, which is how rank deficiency is detected. The default is 20 * (m + n) * eps() * maximum(norm, eachcol(A)) for an m × n matrix.

• ordering: the fill-reducing column ordering, one of the constants SPQR.ORDERING_DEFAULT (SPQR's choice), ORDERING_FIXED and ORDERING_NATURAL (no fill-reducing ordering), ORDERING_COLAMD, ORDERING_AMD (on A'A), ORDERING_METIS, ORDERING_CHOLMOD, ORDERING_BEST (the best of COLAMD, AMD and METIS) and ORDERING_BESTAMD (the best of COLAMD and AMD). With ORDERING_FIXED, F.pcol is the identity unless A is rank deficient, in which case the columns that SPQR finds dependent are moved to the end.

Examples

julia> A = sparse([1,2,3,4], [1,1,2,2], [1.0,1.0,1.0,1.0])
4×2 SparseMatrixCSC{Float64, Int64} with 4 stored entries:
 1.0   ⋅
 1.0   ⋅
  ⋅   1.0
  ⋅   1.0

julia> qr(A)
SparseArrays.SPQR.QRSparse{Float64, Int64}
Q factor:
4×4 SparseArrays.SPQR.QRSparseQ{Float64, Int64}:
 -0.707107   0.0        0.0       -0.707107
  0.0       -0.707107  -0.707107   0.0
  0.0       -0.707107   0.707107   0.0
 -0.707107   0.0        0.0        0.707107
R factor:
2×2 SparseMatrixCSC{Float64, Int64} with 2 stored entries:
 -1.41421      ⋅
     ⋅     -1.41421
Row permutation:
4-element Vector{Int64}:
 1
 3
 4
 2
Column permutation:
2-element Vector{Int64}:
 1
 2

lq(A::SparseMatrixCSC; tol=_default_tol(A'), ordering=ORDERING_DEFAULT) -> AdjointQRSparse

Compute the LQ factorization of a sparse matrix A as the adjoint of the sparse QR factorization of A', that is qr(A')', using SPQR. See qr for the keyword arguments and the sparse Q.

The factorization F satisfies A[F.prow, F.pcol] == F.L * F.Q, where F.L is a lower triangular sparse matrix and F.Q the adjoint of the Q of the QR factorization. F \ b solves the underdetermined system A * x == b for a wide A and returns the minimum-norm solution, as for dense lq. F' is the QR factorization of A', and lq(A') reuses qr(A) without a copy.

Examples

julia> A = sparse([1.0 0 1 0; 0 1 0 1]);

julia> F = lq(A);

julia> F.L * F.Q ≈ A[F.prow, F.pcol]
true

julia> F \ [1.0, 2.0] ≈ Matrix(A) \ [1.0, 2.0]
true

rank(F::QRSparse{Tv,Ti}) -> Ti

Return the numerical rank of the matrix that F factorizes, as SPQR found it with the tol given to qr.

Examples

julia> A = sparse([1.0 1; 1 1; 1 1 + 1e-10]);

julia> rank(qr(A)), rank(qr(A; tol = 1e-8))
(2, 1)

rank(S::SparseMatrixCSC{Tv,Ti}; [tol::Real]) -> Ti

Return the numerical rank of S, computed as rank(qr(S; tol)). A column whose norm falls to tol or below during the factorization is treated as zero; see qr for the default.

SPQR.SpqrWS(F::QRSparse)

Scratch space for ldiv!(x, F, b; workspace), which makes repeated solves allocation-free. Without it, ldiv! allocates its scratch space on each call. A workspace grows as needed, so it can be reused with any factorization of the same element type, but not by two calls at once.

(\)(F::QRSparse, B::StridedVecOrMat)

Solve the least squares problem $\min\|Ax - b\|^2$ or the linear system of equations $Ax=b$ when F is the sparse QR factorization of $A$. A basic solution is returned when the problem is underdetermined; A \ b and factorize(A) \ b instead return the minimum-norm solution through lq, as for dense matrices.

Examples

julia> A = sparse([1,2,4], [1,1,1], [1.0,1.0,1.0], 4, 2)
4×2 SparseMatrixCSC{Float64, Int64} with 3 stored entries:
 1.0  ⋅
 1.0  ⋅
  ⋅   ⋅
 1.0  ⋅

julia> qr(A)\fill(1.0, 4)
2-element Vector{Float64}:
 1.0
 0.0

(\)(F::AdjointFactorization{<:Any,<:QRSparse}, B::StridedVecOrMat)

Solve the underdetermined system $A^*x=b$ when F is the sparse QR factorization of the tall matrix $A$, i.e. F = qr(A) with size(A, 1) >= size(A, 2). The minimum-norm solution is returned; when $A$ is rank deficient, the equations corresponding to the dependent columns of $A$ are dropped, mirroring the basic solution returned by F \ B. Overdetermined systems are not supported here as they would require a factorization of $A^*$ rather than of $A$.

Examples

julia> A = sparse([1,2,3,4,1,2,3,4], [1,1,1,1,2,2,2,2], [1.0,1.0,1.0,1.0,1.0,-1.0,1.0,-1.0])
4×2 SparseMatrixCSC{Float64, Int64} with 8 stored entries:
 1.0   1.0
 1.0  -1.0
 1.0   1.0
 1.0  -1.0

julia> x = qr(A)'\[4.0, 0.0]
4-element Vector{Float64}:
 1.0
 1.0
 1.0
 1.0

julia> A'x
2-element Vector{Float64}:
 4.0
 0.0

copy(F::QRSparse)

Return an independent copy of F, with its own factors and permutations. deepcopy(F) does the same. Solving in parallel does not need a copy: tasks can share one F.

LibSuiteSparse.set_libdir!(dir)
LibSuiteSparse.set_libdir!(nothing)

Load the SuiteSparse libraries from dir instead of the copies bundled with Julia, or restore the bundled copies with nothing. dir must contain the whole set of libraries (libsuitesparseconfig, libamd, libcamd, libcolamd, libccolamd, libcholmod, libspqr and libumfpack) under the same file names as the bundled ones, built from the same major SuiteSparse version.

The libraries are loaded on first use, so this must be called before the first solver call. It throws once any of them has been loaded. The directory can also be set with the JULIA_SUITESPARSE_LIBDIR environment variable; set_libdir! takes precedence.

LibSuiteSparse.libdir()

Return the directory the SuiteSparse libraries are loaded from, or will be loaded from if none has been loaded yet. See set_libdir!.

LibSuiteSparse.init_suitesparse

Internal function which is used to initialize the SuiteSparse libraries to the correct memory management functions. Any package which directly wraps one of the following SuiteSparse libraries must ensure that this function is called before the use of that library: AMD, CAMD, COLAMD, CCOLAMD, UMFPACK, CXSparse, CHOLMOD, KLU, BTF, LDL, RBio, SPQR, SPEX, and ParU

Notes:

• Currently this function only sets the memory management functions of SuiteSparse_config, however there are also override functions for printf, hypot, and divcomplex.
• SuiteSparse_config, and this initialization function, is not a dependency of CSparse, GraphBLAS, or LAGraph.