Linear operators

A LinearOperator is a structure representing a linear operator from a VectorSpace to an other. More precisely, a LinearOperator is comprised of the three fields domain::VectorSpace, codomain::VectorSpace and coefficients::AbsractMatrix with matching dimensions and size.

julia> A = LinearOperator(Taylor(1), Taylor(1), [1 2 ; 3 4])LinearOperator : Taylor(1) → Taylor(1) with coefficients Matrix{Int64}: 1  2 3  4

The three fields domain, codomain and coefficients are accessible via the respective functions of the same name.

julia> domain(A)Taylor(1)julia> codomain(A)Taylor(1)julia> coefficients(A)2×2 Matrix{Int64}: 1  2 3  4

For convenience, the methods zeros, ones, fill and fill! are available:

julia> dom, codom = Taylor(1), Taylor(2)(Taylor(1), Taylor(2))julia> zeros(dom, codom)LinearOperator : Taylor(1) → Taylor(2) with coefficients Matrix{Float64}: 0.0  0.0 0.0  0.0 0.0  0.0julia> ones(dom, codom)LinearOperator : Taylor(1) → Taylor(2) with coefficients Matrix{Float64}: 1.0  1.0 1.0  1.0 1.0  1.0julia> fill(2, dom, codom)LinearOperator : Taylor(1) → Taylor(2) with coefficients Matrix{Int64}: 2  2 2  2 2  2julia> fill!(zeros(dom, codom), 2)LinearOperator : Taylor(1) → Taylor(2) with coefficients Matrix{Float64}: 2.0  2.0 2.0  2.0 2.0  2.0

The coefficients of a LinearOperator are indexed according to the indices of the domain and codomain (as given by indices).

julia> A[0:1,0:1] # indices(domain(A)), indices(codomain(A))2×2 Matrix{Int64}: 1  2 3  4

When the domain and/or the codomain of a LinearOperator is a CartesianSpace, its coefficients can be thought of as a block matrix. The function component extracts a LinearOperator composing the cartesian space; its coefficients are a view into the parent operator, so the two stay glued. The function unpack unglues the operator into the Matrix of its components.

julia> B = LinearOperator(ScalarSpace() × Taylor(1)^2, ScalarSpace() × Taylor(1)^2, reshape(1:25, 5, 5))LinearOperator : 𝕂 × Taylor(1)² → 𝕂 × Taylor(1)² with coefficients Base.ReshapedArray{Int64, 2, UnitRange{Int64}, Tuple{}}: 1   6  11  16  21 2   7  12  17  22 3   8  13  18  23 4   9  14  19  24 5  10  15  20  25julia> B[1:5,1:5] # indices(domain(B)), indices(codomain(B))5×5 Matrix{Int64}: 1   6  11  16  21 2   7  12  17  22 3   8  13  18  23 4   9  14  19  24 5  10  15  20  25julia> component(B, 1, 1) # extract the linear operator associated with the domain ScalarSpace() and codomain ScalarSpace()LinearOperator : 𝕂 → 𝕂 with coefficients SubArray{Int64, 2, Base.ReshapedArray{Int64, 2, UnitRange{Int64}, Tuple{}}, Tuple{UnitRange{Int64}, UnitRange{Int64}}, false}: 1julia> component(B, 2, 2) # extract the linear operator associated with the domain Taylor(1)^2 and codomain Taylor(1)^2LinearOperator : Taylor(1)² → Taylor(1)² with coefficients SubArray{Int64, 2, Base.ReshapedArray{Int64, 2, UnitRange{Int64}, Tuple{}}, Tuple{UnitRange{Int64}, UnitRange{Int64}}, false}:  7  12  17  22  8  13  18  23  9  14  19  24 10  15  20  25julia> component(component(B, 2, 2), 1, 1)LinearOperator : Taylor(1) → Taylor(1) with coefficients SubArray{Int64, 2, Base.ReshapedArray{Int64, 2, UnitRange{Int64}, Tuple{}}, Tuple{UnitRange{Int64}, UnitRange{Int64}}, false}: 7  12 8  13julia> component(component(B, 2, 2), 2, 2)LinearOperator : Taylor(1) → Taylor(1) with coefficients SubArray{Int64, 2, Base.ReshapedArray{Int64, 2, UnitRange{Int64}, Tuple{}}, Tuple{UnitRange{Int64}, UnitRange{Int64}}, false}: 19  24 20  25julia> unpack(B)2×2 Matrix{LinearOperator{T, S, SubArray{Int64, 2, Base.ReshapedArray{Int64, 2, UnitRange{Int64}, Tuple{}}, Tuple{UnitRange{Int64}, UnitRange{Int64}}, false}} where {T<:VectorSpace, S<:VectorSpace}}: [1;;]           [6 11 16 21] [2; 3; 4; 5;;]  [7 12 17 22; 8 13 18 23; 9 14 19 24; 10 15 20 25]

Similarly, the function eachcomponent returns a Generator whose iterates yield each LinearOperator composing the cartesian space.

Arithmetic

The addition and subtraction operations are implemented as the + and - functions respectively.

julia> C = LinearOperator(Taylor(1), Taylor(1), [1 2 ; 3 4])LinearOperator : Taylor(1) → Taylor(1) with coefficients Matrix{Int64}: 1  2 3  4julia> D = LinearOperator(Taylor(1), Taylor(2), [1 2 ; 3 4 ; 5 6])LinearOperator : Taylor(1) → Taylor(2) with coefficients Matrix{Int64}: 1  2 3  4 5  6julia> C + DLinearOperator : Taylor(1) → Taylor(2) with coefficients Matrix{Int64}: 2  4 6  8 5  6julia> C - DLinearOperator : Taylor(1) → Taylor(2) with coefficients Matrix{Int64}:  0   0  0   0 -5  -6julia> C + ILinearOperator(Taylor(1), Taylor(1), [1 2; 3 4]) + UniformScalingOperator{Bool}(true)julia> C - ILinearOperator(Taylor(1), Taylor(1), [1 2; 3 4]) + UniformScalingOperator{Int64}(-1)

The product between LinearOperator is implemented as the * and ^ functions. The division between LinearOperator is implemented as the \ method.

julia> C * DLinearOperator : Taylor(1) → Taylor(1) with coefficients Matrix{Int64}:  7  10 15  22julia> C ^ 3LinearOperator : Taylor(1) → Taylor(1) with coefficients Matrix{Int64}: 37   54 81  118julia> C \ CLinearOperator : Taylor(1) → Taylor(1) with coefficients Matrix{Float64}: 1.0  0.0 0.0  1.0

The action of a LinearOperator is performed by the right product * of a LinearOperator with a Sequence; alternatively, LinearOperator defines a method on a Sequence representing *. Naturally, the resulting sequence is an element of the codomain of the LinearOperator.

Conversely, the operator \ between a LinearOperator and a Sequence corresponds to the action of the inverse of the LinearOperator; the output sequence is an element of the domain of the LinearOperator.

julia> x = Sequence(Taylor(2), [1, 1, 1])Sequence in Taylor(2) with coefficients Vector{Int64}: 1 1 1julia> C * x # C(x)Sequence in Taylor(1) with coefficients Vector{Int64}: 3 7julia> D \ xSequence in Taylor(1) with coefficients Vector{Float64}: -1.0000000000000018  1.0000000000000013

API

AbstractLinearOperator

LinearOperator{T<:VectorSpace,S<:VectorSpace,R<:AbstractMatrix} <: AbstractLinearOperator

julia> LinearOperator(Taylor(1), Taylor(1), [1 2 ; 3 4])LinearOperator : Taylor(1) → Taylor(1) with coefficients Matrix{Int64}: 1  2 3  4julia> LinearOperator(Taylor(2), ScalarSpace(), [1.0 0.5 0.25])LinearOperator : Taylor(2) → 𝕂 with coefficients Matrix{Float64}: 1.0  0.5  0.25julia> LinearOperator([1 2 3 ; 4 5 6])LinearOperator : 𝕂³ → 𝕂² with coefficients Matrix{Int64}: 1  2  3 4  5  6

component(A::LinearOperator, i, j)

unpack(A::LinearOperator)