Reference

HSL MA57 and MA97 can be used for the solution of symmetric, possibly indefinite, linear systems. They are often used for the solution of saddle-point systems, i.e., systems of the form

\[\begin{bmatrix} H & A^T \\ A & \end{bmatrix} \begin{bmatrix} x \\ y \end{bmatrix} = \begin{bmatrix} b \\ c \end{bmatrix}.\]

Special cases occur when $H$ is the identity matrix and either $b = 0$ or $c = 0$, which correspond to the least-squares problem

\[\min_y \ \|A^T y - b\|\]

and to the least-norm problem

\[\min_x \ \|x\| \ \text{subject to } Ax = c\]

respectively.

HSL.Ma57 — Type.

Instantiate an object of type Ma57 and perform the symbolic analysis on a sparse Julia matrix.

M = Ma57(A; kwargs...)

Input arguments

A::SparseMatrixCSC{T<:Ma57Data,Int}: input matrix. The lower triangle will be extracted.

Keyword arguments

All keyword arguments are passed directly to ma57_coord().

Example:


julia> T = Float64;

julia> rows = Int32[1, 1, 2, 2, 3, 3, 5]; cols = Int32[1, 2, 3, 5, 3, 4, 5];

julia> vals = T[2, 3, 4, 6, 1, 5, 1];

julia> A = sparse(rows, cols, vals); A = A + triu(A, 1)';

julia> M = Ma57(A)

HSL.Ma57{Float64}(5, 7, Int32[1, 1, 2, 2, 3, 3, 5], Int32[1, 2, 3, 5, 3, 4, 5], [2.0, 3.0, 4.0, 6.0, 1.0, 5.0, 1.0], HSL.Ma57_Control{Float64}(Int32[6, 6, 6, -1, 0, 5, 1, 0, 10, 1, 16, 16, 10, 100, 1, 0, 0, 0, 0, 0], [0.01, 1.0e-20, 0.5, 0.0, 0.0]), HSL.Ma57_Info{Float64}(Int32[0, 0, 0, 0, 12, 13, 4, 2, 48, 53  …  0, 0, 0, 0, 0, 2, 0, 0, 0, 0], [10.0, 34.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], 0, 0, 0, 0, 0, 0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0), 1.1, 81, Int32[5, 4, 3, 2, 1, 2, 9, 0, 0, 0  …  4, 3, 3, 2, 2, 1, 1, 0, 0, 0], 48, Float64[], 53, Int32[])

HSL.Ma57Exception — Type.

Exception type raised in case of error.

HSL.Ma57_Control — Type.

Main control type for MA57.

Ma57_Control(; kwargs...)

Keyword arguments:

print_level::Int: integer controling the verbosit level. Accepted values are:
- <0: no printing
- 0: errors and warnings only (default)
- 1: errors, warnings and basic diagnostics
- 2: errors, warning and full diagnostics
unit_diagnostics::Int: Fortran file unit for diagnostics (default: 6)
unit_error::Int: Fortran file unit for errors (default: 6)
unit_warning::Int: Fortran file unit for warnings (default: 6)

Example:

julia> using HSL
julia> Ma57_Control{Float64}(print_level=1)
HSL.Ma57_Control{Float64}(Int32[6, 6, 6, -1, 1, 5, 1, 0, 10, 1, 16, 16, 10, 100, 1, 0, 0, 0, 0, 0], [0.01, 1.0e-20, 0.5, 0.0, 0.0])

HSL.Ma57_Info — Type.

Main info type for MA57

info = Ma57_Info{T <: Ma97Real}()

An info variable is used to collect statistics on the analysis, factorization, and solve.

Example:

julia> using HSL
julia> T = Float64;
julia> rows = Int32[1, 1, 2, 2, 3, 3, 5]; cols = Int32[1, 2, 3, 5, 3, 4, 5];
julia> vals = T[2, 3, 4, 6, 1, 5, 1];
julia> A = sparse(rows, cols, vals); A = A + triu(A, 1)';
julia> M = Ma57(A);
julia> M.info
HSL.Ma57_Info{Float64}(Int32[0, 0, 0, 0, 12, 13, 4, 2, 48, 53  …  0, 0, 0, 0, 0, 2, 0, 0, 0, 0], [10.0, 34.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], 0, 0, 0, 0, 0, 0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0)

HSL.Ma97 — Method.

Instantiate and perform symbolic analysis on a sparse Julia matrix

M = Ma97(A; kwargs...)

Instantiate an object of type Ma97 and perform the symbolic analysis on a sparse Julia matrix.

Input arguments

A::SparseMatrixCSC{T<:Ma97Data,Int}: input matrix. The lower triangle will be extracted.

Keyword arguments

All keyword arguments are passed directly to ma97_csc().

HSL.Ma97Exception — Type.

Exception type raised in case of error.

HSL.Ma97_Control — Type.

Main control type for MA97.

Ma97_Control(; kwargs...)

Keyword arguments:

print_level::Int: integer controling the verbosit level. Accepted values are:
- <0: no printing
- 0: errors and warnings only (default)
- 1: errors, warnings and basic diagnostics
- 2: errors, warning and full diagnostics
unit_diagnostics::Int: Fortran file unit for diagnostics (default: 6)
unit_error::Int: Fortran file unit for errors (default: 6)
unit_warning::Int: Fortran file unit for warnings (default: 6)

HSL.Ma97_Info — Type.

Main info type for MA97

info = Ma97_Info{T <: Ma97Real}()

An info variable is used to collect statistics on the analysis, factorization and solve.

HSL.ma57_alter_d — Function.

Alter block diagonal matrix D

Input arguments:

M::Ma57: Ma57 object
F::Array{Float64,2}: diagonal adjustment

Example:

julia> using HSL

julia> T = Float64;

julia> rows = Int32[1, 1, 2, 2, 3, 3, 5]; cols = Int32[1, 2, 3, 5, 3, 4, 5];

julia> vals = T[2, 3, 4, 6, 1, 5, 1];

julia> A = sparse(rows, cols, vals); A = A + triu(A, 1)';

julia> b = T[8, 45, 31, 15, 17]

julia> ϵ = sqrt(eps(eltype(A)))

julia> xexact = T[1, 2, 3, 4, 5]

julia> M = Ma57(A)

julia> ma57_factorize(M)

julia> (L, D, s, p) = ma57_get_factors(M)

julia> d1 = abs.(diag(D))

julia> d2 = [diag(D, 1) ; 0]

julia> F = [full(d1)' ; full(d2)']

julia> ma57_alter_d(M, F)

HSL.ma57_factorize — Function.

Factorize Ma57 object.

ma57_factorize(M)

Input arguments:

M::Ma57: Ma57 object

Return values:

(none)

Stored information:

M.info.largest_front::Int: order of largest frontal matrix
M.info.num_2x2_pivots::Int: number of 2x2 pivots used in factorization
M.info.num_delayed_pivots::Int: total number of fully-summed variables that were passed to the father node because of pivoting considerations
M.info.num_negative_eigs::Int: number of negative eigenvalues in factorization of M
M.info.rank::Int: rank of factorization of M
M.info.num_pivot_sign_changes::Int: number of sign changes of pivot when icntl(7) = 3 (ie, no pivoting)

Factorize a sparse matrix.

Input arguments:

A::SparseMatrixCSC{T<:Ma57Data,Int}: sparse matrix

Return values:

M::Ma57: factorized Ma57 object

Stored information:

M.info.largest_front::Int: order of largest frontal matrix
M.info.num_2x2_pivots::Int: number of 2x2 pivots used in factorization
M.info.num_delayed_pivots::Int: total number of fully-summed variables that were passed to the father node because of pivoting considerations
M.info.num_negative_eigs::Int: number of negative eigenvalues in factorization of A
M.info.rank::Int: rank of factorization of A
M.info.num_pivot_sign_changes::Int: number of sign changes of pivot when icntl(7) = 3 (ie, no pivoting)

Example:


julia> using HSL

julia> T = Float64;

julia> rows = Int32[1, 1, 2, 2, 3, 3, 5]; cols = Int32[1, 2, 3, 5, 3, 4, 5];

julia> vals = T[2, 3, 4, 6, 1, 5, 1];

julia> A = sparse(rows, cols, vals); A = A + triu(A, 1)';

julia> M = Ma57(A)

julia> ma57_factorize(M)      ## factorize `Ma57` object in place

julia> F = ma57_factorize(A)  ## factorize sparse matrix and return `Ma57` object

julia> M.info.largest_front

4

julia> A.info.largest_front   ## same result

4

HSL.ma57_get_factors — Function.

Retrieve factors, scaling, and permutation.

L, D, s, iperm = ma57getfactors(M)

Function will retrieve a unit triangular matrix L, a block-diagonal matrix D a scaling vector s and a permutation vector p such that

P * S * A * S * P' = L * D * L'

where S = spdiagm(s) and P = speye(n)[p,:].

Note that the numerical factorization must have been performed and have succeeded.

Input arguments

M::Ma57: factorized Ma57 object

Return values

L::SparseMatrixCSC{T<:Ma57Data,Int}: L factor
D::SparseMatrixCSC{T<:Ma57Data,Int}: D factor
s::Array{T,1}: diagonal components of scaling matrix S
iperm::Array{Int,1}: array representin permutation matrix P

Example:


julia> using HSL

julia> T = Float64;

julia> rows = Int32[1, 1, 2, 2, 3, 3, 5]; cols = Int32[1, 2, 3, 5, 3, 4, 5];

julia> vals = T[2, 3, 4, 6, 1, 5, 1];

julia> A = sparse(rows, cols, vals); A = A + triu(A, 1)';

julia> b = T[8, 45, 31, 15, 17]

julia> ϵ = sqrt(eps(eltype(A)))

julia> xexact = T[1, 2, 3, 4, 5]

julia> M = Ma57(A)

julia> ma57_factorize(M)

julia> (L, D, s, p) = ma57_get_factors(M)

julia> S = spdiagm(s)

julia> P = speye(M.n)[p, :]

julia> vecnorm(P * S * A * S * P' - L * D * L') ≤ ϵ * vecnorm(A)

true

HSL.ma57_least_squares — Method.

Solve the least-squares problem minimize ‖Ax - b‖, where A has shape m-by-n with m > n, by solving the saddle-point system [ I A ] [ r ] [ b ] [ A' ] [ x ] = [ 0 ].

HSL.ma57_min_norm — Method.

Solve the minimum-norm problem minimize ‖x‖ subject to Ax=b, where A has shape m-by-n with m < n, by solving the saddle-point system [ I A' ] [ x ] [ 0 ] [ A ] [ y ] = [ b ].

HSL.ma57_solve — Function.

System solve.

Solve after factorization without iterative refinement

ma57_solve(ma57, b; kwargs...)

Input arguments:

ma57::Ma57{T<:Ma57Data}: an Ma57 structure for which the analysis and factorization have been performed
b::Array{T}: vector of array of right-hand sides. Note that b will be overwritten. To solve a system with multiple right-hand sides, b should have size n by nrhs.

Keyword arguments:

job::Symbol=:A: task to perform. Accepted values are
- :A: solve Ax = b
- :LS: solve LPSx = PSb
- :DS: solve DPS⁻¹x = PSb
- :LPS: solve L'PS⁻¹x = PS⁻¹b

Return values:

x::Array{T}: an array of the same size as b containing the solutions.

Solve after factorization with iterative refinement

ma57_solve(A, b, nitref; kwargs...)

Input arguments:

A::SparseMatrixCSC{T<:Ma57Data,Int}: input matrix. A full matrix will be converted to sparse.
b::Array{T}: vector of array of right-hand sides. Note that b will be overwritten. To solve a system with multiple right-hand sides, b should have size n by nrhs.
nitref::Int: number of iterative refinement steps

Return values:

x::Array{T}: an array of the same size as b containing the solutions.

Stored information:

Accessible through the Ma57 matrix object's `info` attribute

ma57.info.backward_error1::T: max{i} |b - Ax|i / (|b| + |A| |x|)_i
ma57.info.backward_error2::T: max{i} |b - Ax|i / ((|A| |x|)i + ||Ai||{∞} ||x||{∞})
ma57.info.matrix_inf_norm::T: ||A||_{∞}
ma57.info.solution_inf_norm::T: ||x||_{∞}
ma57.info.scaled_residuals::T: norm of scaled residuals = max{i} |sumj a{ij} xj - bi| / ||A||{∞} ||x||_{∞})
ma57.info.cond1::T: condition number as defined in Arioli, M. Demmel, J. W., and Duff, I. S. (1989). Solving sparse linear systems with sparse backward error. SIAM J.Matrix Anal. and Applics. 10, 165-190.
ma57.info.cond2::T: condition number as defined in Arioli, M. Demmel, J. W., and Duff, I. S. (1989). Solving sparse linear systems with sparse backward error. SIAM J.Matrix Anal. and Applics. 10, 165-190.
ma57.info.error_inf_norm::T: upper bound for the infinity norm of the error in the solution

Example:


julia> using HSL

julia> T = Float64;

julia> rows = Int32[1, 1, 2, 2, 3, 3, 5]; cols = Int32[1, 2, 3, 5, 3, 4, 5];

julia> vals = T[2, 3, 4, 6, 1, 5, 1];

julia> A = sparse(rows, cols, vals); A = A + triu(A, 1)';

julia> b = T[8, 45, 31, 15, 17]

julia> ϵ = sqrt(eps(eltype(A)))

julia> xexact = T[1, 2, 3, 4, 5]

julia> M = Ma57(A)

julia> ma57_factorize(M)

julia> x = ma57_solve(M, b)      ## solve without iterative refinement

julia> norm(x - xexact) ≤ ϵ * norm(xexact)

true

julia> xx = ma57_solve(M, b, 2)  ## solve with iterative refinement

julia> norm(xx - xexact) ≤ ϵ * norm(xexact)

true

HSL.ma57_solve! — Method.

Solve a symmetric linear system. The symbolic analysis and numerical factorization must have been performed and must have succeeded.

HSL.ma57_solve! — Method.

Solve a symmetric linear system. The symbolic analysis and numerical factorization must have been performed and must have succeeded.

HSL.ma97_coord — Function.

Instantiate and perform symbolic analysis using coordinate arrays

M = ma97_coord(n, cols, rows, nzval; kwargs...)

Instantiate an object of type Ma97 and perform the symbolic analysis on a matrix described in sparse coordinate format.

Input arguments

n::Int: the matrix size
cols::Vector{T<:Integer}: array of column indices for the lower triangle
rows::Vector{T<:Integer}: array of row indices for the lower triangle
nzval::Vector{T<:Ma97Data}: array of values for the lower triangle

Keyword arguments

All keyword arguments are passed directly to the Ma97_Control constructor.

HSL.ma97_csc — Function.

Instantiate and perform symbolic analysis using CSC arrays

M = ma97_csc(n, colptr, rowval, nzval; kwargs...)

Instantiate an object of type Ma97 and perform the symbolic analysis on a matrix described in sparse CSC format.

Input arguments

n::Int: the matrix size
colptr::Vector{T<:Integer}: CSC colptr array for the lower triangle
rowval::Vector{T<:Integer}: CSC rowval array for the lower triangle
nzval::Vector{T<:Ma97Data}: CSC nzval array for the lower triangle

Keyword arguments

All keyword arguments are passed directly to the Ma97_Control constructor.

HSL.ma97_enquire — Function.

ma97enquire: see the documentation for `ma97inquire`.

HSL.ma97_factorise — Function.

ma97factorise: see the documentation for `ma97factorize`.

HSL.ma97_factorise! — Function.

ma97factorise!: see the documentation for `ma97factorize!`.

HSL.ma97_factorize — Function.

Combined Analysis and factorization

M = ma97_factorize(A; kwargs...)

Convenience method that combines the symbolic analysis and numerical factorization phases. An MA97 instance is returned, that can subsequently be passed to other functions, e.g., ma97_solve().

Input Arguments

A::SparseMatrixCSC{T<:Ma97Data,Int}: Julia sparse matrix

Keyword Arguments

matrix_type::Symbol=:real_indef: indicates the matrix type. Accepted values are
- :real_spd for a real symmetric and positive definite matrix
- :real_indef for a real symmetric and indefinite matrix.

HSL.ma97_factorize! — Function.

Perform numerical factorization.

ma97_factorize!(ma97; kwargs...)

The symbolic analysis must have been performed and must have succeeded.

Input Arguments

ma97::Ma97{T<:Ma97Data}:: an Ma97 structure for which the analysis has been performed

Keyword Arguments

matrix_type::Symbol=:real_indef: indicates the matrix type. Accepted values are
- :real_spd for a real symmetric and positive definite matrix
- :real_indef for a real symmetric and indefinite matrix.

HSL.ma97_inquire — Function.

Inquire about a factorization or solve

ma97_inquire(ma97; kwargs...)

Obtain information on the pivots after a successful factorization or solve.

Input Arguments

ma97::Ma97{T<:Ma97Data}: an Ma97 structure for which the analysis and factorization have been performed

Keyword arguments

matrix_type::Symbol=:real_indef: indicates the matrix type. Accepted values are
- :real_spd for a real symmetric and positive definite matrix
- :real_indef for a real symmetric and indefinite matrix.

Return values

An inquiry on a real or complex indefinite matrix returns two vectors:

piv_order: contains the pivot sequence; a negative value indicates that the corresponding variable is part of a 2x2 pivot,
d: a 2 by n array whose first row contains the diagonal of D⁻¹ in the factorization, and whose nonzeros in the second row contain the off-diagonals.

An inquiry on a positive definite matrix returns one vector with the pivot values.

HSL.ma97_least_squares — Method.

Solve least-squares problem

ma97_least_squares(A, b)

Solve the least-squares problem

minimize ‖Ax - b‖₂

where A has shape m-by-n with m > n, by solving the saddle-point system

[ I   A ] [ r ]   [ b ]
[ A'    ] [ x ] = [ 0 ].

Input arguments

A::SparseMatrixCSC{T<:Ma97Data,Int}: input matrix of shape m-by-n with m > n. A full matrix will be converted to sparse.
b::Vector{T}: right-hand side vector

Return value

x::Vector{T}: solution vector.

HSL.ma97_min_norm — Method.

Solve a minimum-norm problem

ma97_min_norm(A, b)

solves

minimize ‖x‖₂  subject to Ax=b,

where A has shape m-by-n with m < n, by solving the saddle-point system

[ I  A' ] [ x ]   [ 0 ]
[ A     ] [ y ] = [ b ].

Input arguments

A::SparseMatrixCSC{T<:Ma97Data,Int}: input matrix of shape m-by-n with m < n. A full matrix will be converted to sparse.
b::Vector{T}: right-hand side vector

Return value

x::Vector{T}: solution vector.

HSL.ma97_solve — Function.

System solve

Solve after factorization

ma97_solve(ma97, b; kwargs...)

Input arguments

ma97::Ma97{T<:Ma97Data}: an Ma97 structure for which the analysis and factorization have been performed
b::Array{T}: vector of array of right-hand sides. Note that b will be overwritten. To solve a system with multiple right-hand sides, b should have size n by nrhs.

Keyword arguments

job::Symbol=:A: task to perform. Accepted values are
- :A: solve Ax = b
- :PL: solve PLx = Sb
- :D: solve Dx = b
- :LPS: solve L'P'S⁻¹x = b
- :DLPS: solve DL'P'S⁻¹x = b.

Return values

x::Array{T}: an array of the same size as b containing the solutions.

Combined analysis, factorization and solve

ma97_solve(A, b; kwargs...)

Input arguments

A::SparseMatrixCSC{T<:Ma97Data,Int}: input matrix. A full matrix will be converted to sparse.
b::Array{T}: vector of array of right-hand sides. Note that b will be overwritten. To solve a system with multiple right-hand sides, b should have size n by nrhs.

Keyword arguments

matrix_type::Symbol=:real_indef: indicates the matrix type. Accepted values are
- :real_spd for a real symmetric and positive definite matrix
- :real_indef for a real symmetric and indefinite matrix.

Return values

x::Array{T}: an array of the same size as b containing the solutions.

HSL.ma97_solve! — Function.

In-place system solve

See the documentation for ma97_solve(). The only difference is that the right-hand side b is overwritten with the solution.