Reference

Contents

Index

Base.MatrixMethod
A = Matrix(op)

Materialize an operator as a dense array using op.ncol products.

source
LinearOperators.DiagonalQNType

Implementation of the diagonal quasi-Newton approximations described in

M. Zhu, J. L. Nazareth and H. Wolkowicz The Quasi-Cauchy Relation and Diagonal Updating. SIAM Journal on Optimization, vol. 9, number 4, pp. 1192-1204, 1999. https://doi.org/10.1137/S1052623498331793.

and

Andrei, N. A diagonal quasi-Newton updating method for unconstrained optimization. https://doi.org/10.1007/s11075-018-0562-7

source
LinearOperators.DiagonalQNMethod
DiagonalQN(d)

Construct a linear operator that represents a diagonal quasi-Newton approximation. The approximation satisfies the weak secant equation and is not guaranteed to be positive definite.

Arguments

  • d::AbstractVector: initial diagonal approximation;
  • psb::Bool: whether to use the diagonal PSB update or the Andrei update.
source
LinearOperators.LBFGSOperatorMethod
LBFGSOperator(T, n; [mem=5, scaling=true])
LBFGSOperator(n; [mem=5, scaling=true])

Construct a limited-memory BFGS approximation in forward form. If the type T is omitted, then Float64 is used.

source
LinearOperators.LSR1OperatorMethod
LSR1Operator(T, n; [mem=5, scaling=false)
LSR1Operator(n; [mem=5, scaling=false)

Construct a limited-memory SR1 approximation in forward form. If the type T is omitted, then Float64 is used.

source
LinearOperators.LinearOperatorType

Base type to represent a linear operator. The usual arithmetic operations may be applied to operators to combine or otherwise alter them. They can be combined with other operators, with matrices and with scalars. Operators may be transposed and conjugate-transposed using the usual Julia syntax.

source
LinearOperators.LinearOperatorMethod
LinearOperator(M::AbstractMatrix{T}; symmetric=false, hermitian=false, S = Vector{T}) where {T}

Construct a linear operator from a dense or sparse matrix. Use the optional keyword arguments to indicate whether the operator is symmetric and/or hermitian. Change S to use LinearOperators on GPU.

source
LinearOperators.LinearOperatorMethod
LinearOperator(type::Type{T}, nrow, ncol, symmetric, hermitian, prod!,
                [tprod!=nothing, ctprod!=nothing],
                S = Vector{T}) where {T}

Construct a linear operator from functions where the type is specified as the first argument. Change S to use LinearOperators on GPU. Notice that the linear operator does not enforce the type, so using a wrong type can result in errors. For instance,

A = [im 1.0; 0.0 1.0] # Complex matrix
function mulOp!(res, M, v, α, β)
  mul!(res, M, v, α, β)
end
op = LinearOperator(Float64, 2, 2, false, false, 
                    (res, v, α, β) -> mulOp!(res, A, v, α, β), 
                    (res, u, α, β) -> mulOp!(res, transpose(A), u, α, β), 
                    (res, w, α, β) -> mulOp!(res, A', w, α, β))
Matrix(op) # InexactError

The error is caused because Matrix(op) tries to create a Float64 matrix with the contents of the complex matrix A.

Using * may generate a vector that contains NaN values. This can also happen if you use the 3-args mul! function with a preallocated vector such as Vector{Float64}(undef, n). To fix this issue you will have to deal with the cases β == 0 and β != 0 separately:

d1 = [2.0; 3.0]
function mulSquareOpDiagonal!(res, d, v, α, β::T) where T
  if β == zero(T)
    res .= α .* d .* v
  else 
    res .= α .* d .* v .+ β .* res
  end
end
op = LinearOperator(Float64, 2, 2, true, true, 
                    (res, v, α, β) -> mulSquareOpDiagonal!(res, d, v, α, β))

It is possible to create an operator with the 3-args mul!. In this case, using the 5-args mul! will generate storage vectors.

A = rand(2, 2)
op = LinearOperator(Float64, 2, 2, false, false, 
                    (res, v) -> mul!(res, A, v),
                    (res, w) -> mul!(res, A', w))
source
LinearOperators.LinearOperatorMethod
LinearOperator(M::Hermitian{T}, S = Vector{T}) where {T}

Constructs a linear operator from a Hermitian matrix. If its elements are real, it is also symmetric. Change S to use LinearOperators on GPU.

source
LinearOperators.LinearOperatorMethod
LinearOperator(M::SymTridiagonal{T}, S = Vector{T}) where {T}

Constructs a linear operator from a symmetric tridiagonal matrix. If its elements are real, it is also Hermitian, otherwise complex symmetric. Change S to use LinearOperators on GPU.

source
LinearOperators.LinearOperatorMethod

LinearOperator(M::Symmetric{T}, S = Vector{T}) where {T <: Real} =

Construct a linear operator from a symmetric real square matrix M. Change S to use LinearOperators on GPU.

source
LinearOperators.SpectralGradientType

Implementation of a spectral gradient quasi-Newton approximation described in

Birgin, E. G., Martínez, J. M., & Raydan, M. Spectral Projected Gradient Methods: Review and Perspectives. https://doi.org/10.18637/jss.v060.i03

source
LinearOperators.SpectralGradientMethod
    SpectralGradient(σ, n)

Construct a spectral gradient Hessian approximation. The approximation is defined as σI.

Arguments

  • σ::Real: initial positive multiple of the identity;
  • n::Int: operator size.
source
LinearOperators.opEyeMethod
opEye(T, n; S = Vector{T})
opEye(n)

Identity operator of order n and of data type T (defaults to Float64). Change S to use LinearOperators on GPU.

source
LinearOperators.opEyeMethod
opEye(T, nrow, ncol; S = Vector{T})
opEye(nrow, ncol)

Rectangular identity operator of size nrowxncol and of data type T (defaults to Float64). Change S to use LinearOperators on GPU.

source
Base.kronMethod

kron(A, B)

Kronecker tensor product of A and B in linear operator form, if either or both are linear operators. If both A and B are matrices, then Base.kron is used.

source
Base.push!Method
push!(op, s, y)

Push a new {s,y} pair into a L-SR1 operator.

source
Base.push!Method
push!(op, s, y)
push!(op, s, y, Bs)
push!(op, s, y, α, g)
push!(op, s, y, α, g, Bs)

Push a new {s,y} pair into a L-BFGS operator. The second calling sequence is used for forward updating damping, using the preallocated vector Bs. If the operator is damped, the first call will create Bs and call the second call. The third and fourth calling sequences are used in inverse LBFGS updating in conjunction with damping, where α is the most recent steplength and g the gradient used when solving d=-Hg.

source
Base.showMethod
show(io, op)

Display basic information about a linear operator.

source
Base.sizeMethod
m = size(op, d)

Return the size of a linear operator along dimension d.

source
Base.sizeMethod
m, n = size(op)

Return the size of a linear operator as a tuple.

source
LinearOperators.BlockDiagonalOperatorMethod
BlockDiagonalOperator(M1, M2, ..., Mn; S = promote_type(storage_type.(M1, M2, ..., Mn)))

Creates a block-diagonal linear operator:

[ M1           ]
[    M2        ]
[       ...    ]
[           Mn ]

Change S to use LinearOperators on GPU.

source
LinearOperators.InverseLBFGSOperatorMethod
InverseLBFGSOperator(T, n, [mem=5; scaling=true])
InverseLBFGSOperator(n, [mem=5; scaling=true])

Construct a limited-memory BFGS approximation in inverse form. If the type T is omitted, then Float64 is used.

source
LinearOperators.has_args5Method
has_args5(op)

Determine whether the operator can work with the 5-args mul!. If false, storage vectors will be generated at the first call of the 5-args mul!. No additional vectors are generated when using the 3-args mul!.

source
LinearOperators.normestFunction

normest(S) estimates the matrix 2-norm of S. This function is an adaptation of Matlab's built-in NORMEST. This method allocates.


Inputs: S –- Matrix or LinearOperator type, tol –- relative error tol, default(or -1) Machine eps maxiter –- maximum iteration, default 100

Returns: e –- the estimated norm cnt –- the number of iterations used

source
LinearOperators.opCholeskyMethod
opCholesky(M, [check=false])

Inverse of a Hermitian and positive definite matrix as a linear operator using its Cholesky factorization. The factorization is computed only once. The optional check argument will perform cheap hermicity and definiteness checks. This Operator is not in-place when using mul!.

source
LinearOperators.opExtensionMethod
Z = opExtension(I, ncol)
Z = opExtension(:, ncol)

Creates a LinearOperator extending a vector of size length(I) to size ncol, where the position of the elements on the new vector are given by the indices I. The operation w = Z * v is equivalent to w = zeros(ncol); w[I] = v.

Z = opExtension(k, ncol)

Alias for opExtension([k], ncol).

source
LinearOperators.opInverseMethod
opInverse(M; symm=false, herm=false)

Inverse of a matrix as a linear operator using \. Useful for triangular matrices. Note that each application of this operator applies \. This Operator is not in-place when using mul!.

source
LinearOperators.opLDLMethod
opLDL(M, [check=false])

Inverse of a symmetric matrix as a linear operator using its LDLᵀ factorization if it exists. The factorization is computed only once. The optional check argument will perform a cheap hermicity check.

If M is sparse and real, then only the upper triangle should be stored in order to use LDLFactorizations.jl:

triu!(M)
opLDL(Symmetric(M, :U))
source
LinearOperators.opOnesMethod
opOnes(T, nrow, ncol; S = Vector{T})
opOnes(nrow, ncol)

Operator of all ones of size nrow-by-ncol of data type T (defaults to Float64). Change S to use LinearOperators on GPU.

source
LinearOperators.opRestrictionMethod
Z = opRestriction(I, ncol)
Z = opRestriction(:, ncol)

Creates a LinearOperator restricting a ncol-sized vector to indices I. The operation Z * v is equivalent to v[I]. I can be :.

Z = opRestriction(k, ncol)

Alias for opRestriction([k], ncol).

source
LinearOperators.opZerosMethod
opZeros(T, nrow, ncol; S = Vector{T})
opZeros(nrow, ncol)

Zero operator of size nrow-by-ncol, of data type T (defaults to Float64). Change S to use LinearOperators on GPU.

source