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.0
julia> ones(dom, codom)LinearOperator : Taylor(1) → Taylor(2) with coefficients Matrix{Float64}:
 1.0  1.0
 1.0  1.0
 1.0  1.0
julia> fill(2, dom, codom)LinearOperator : Taylor(1) → Taylor(2) with coefficients Matrix{Int64}:
 2  2
 2  2
 2  2
julia> 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.

julia> B = LinearOperator(ParameterSpace() × Taylor(1)^2, ParameterSpace() × 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  25
julia> 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  25
julia> component(B, 1, 1) # extract the linear operator associated with the domain ParameterSpace() and codomain ParameterSpace()LinearOperator : 𝕂 → 𝕂 with coefficients SubArray{Int64, 2, Base.ReshapedArray{Int64, 2, UnitRange{Int64}, Tuple{}}, Tuple{UnitRange{Int64}, UnitRange{Int64}}, false}:
 1
julia> 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  25
julia> 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  13
julia> 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  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. Their bar counterparts add_bar (unicode alias +\bar<tab>) and sub_bar (unicode alias -\bar<tab>) give the result projected in the smallest compatible domain and codomain between the operands.

julia> C = LinearOperator(Taylor(1), Taylor(1), [1 2 ; 3 4])LinearOperator : Taylor(1) → Taylor(1) with coefficients Matrix{Int64}:
 1  2
 3  4
julia> 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  6
julia> C + DLinearOperator : Taylor(1) → Taylor(2) with coefficients Matrix{Int64}:
 2  4
 6  8
 5  6
julia> C - DLinearOperator : Taylor(1) → Taylor(2) with coefficients Matrix{Int64}:
  0   0
  0   0
 -5  -6
julia> add_bar(C, D) # project(C + D, Taylor(1), Taylor(1))LinearOperator : Taylor(1) → Taylor(1) with coefficients Matrix{Int64}:
 2  4
 6  8
julia> sub_bar(C, D) # project(C - D, Taylor(1), Taylor(1))LinearOperator : Taylor(1) → Taylor(1) with coefficients Matrix{Int64}:
 0  0
 0  0
julia> C + ILinearOperator : Taylor(1) → Taylor(1) with coefficients Matrix{Int64}:
 2  2
 3  5
julia> C - ILinearOperator : Taylor(1) → Taylor(1) with coefficients Matrix{Int64}:
 0  2
 3  3

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  22
julia> C ^ 3LinearOperator : Taylor(1) → Taylor(1) with coefficients Matrix{Int64}:
 37   54
 81  118
julia> 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
 1
julia> C * x # C(x)Sequence in Taylor(1) with coefficients Vector{Int64}:
 3
 7
julia> D \ xSequence in Taylor(1) with coefficients Vector{Float64}:
 -1.0000000000000018
  1.0000000000000013

API

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

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

julia> LinearOperator(Taylor(2), ParameterSpace(), [1.0 0.5 0.25])
LinearOperator : Taylor(2) → 𝕂 with coefficients Matrix{Float64}:
 1.0  0.5  0.25

julia> LinearOperator([1 2 3 ; 4 5 6])
LinearOperator : 𝕂³ → 𝕂² with coefficients Matrix{Int64}:
 1  2  3
 4  5  6