Propagate Tuples in InterlaceOperator (#259)

* propagate tuples in InterlateOperator

* Add vector parameter to ArraySpace

* Bugfix to Fun MatrixSppace constructor

* relax MatrixSpace sig in ProductSpaceOperators

Author: jishnub

src/Operators/Operator.jl

@@ -14,6 +14,8 @@ Abstract type to represent linear operators between spaces.
 """
 abstract type Operator{T} end #T is the entry type, Float64 or Complex{Float64}
 
+const VectorOrTupleOfOp{O<:Operator} = Union{AbstractVector{O}, Tuple{O, Vararg{O}}}
+
 eltype(::Operator{T}) where {T} = T
 eltype(::Type{<:Operator{T}}) where {T} = T
 eltype(::Type{OT}) where {OT<:Operator} = eltype(supertype(OT))

src/Operators/general/InterlaceOperator.jl

@@ -31,7 +31,7 @@ function spacescompatible(A::AbstractMatrix{T}) where T<:Operator
     true
 end
 
-spacescompatible(A::AbstractVector{T}) where {T<:Operator} = spacescompatible(map(domainspace,A))
+spacescompatible(A::VectorOrTupleOfOp) = spacescompatible(map(domainspace,A))
 
 function domainspace(A::AbstractMatrix{T}) where T<:Operator
     if !spacescompatible(A)
@@ -42,13 +42,12 @@ function domainspace(A::AbstractMatrix{T}) where T<:Operator
     Space(spl)
 end
 
-function rangespace(A::AbstractVector{T}) where T<:Operator
+function rangespace(A::VectorOrTupleOfOp)
     if !spacescompatible(A)
         error("Cannot construct rangespace for $A as domain spaces are not compatible")
     end
-
-    spl=map(rangespace,A)
-    Space(spl)
+    spl=map(rangespace, A)
+    ArraySpace(_convert_vector_or_svector(spl), first(spl))
 end
 
 promotespaces(A::AbstractMatrix{<:Operator}) = promotespaces(Matrix(A))
@@ -110,15 +109,16 @@ function InterlaceOperator(ops::AbstractMatrix{<:Operator},ds::Space,rs::Space)
         l,u = (1-dimension(rs),dimension(ds)-1)  # not banded
     end
 
-    opsm = convert(Matrix{Operator{mapreduce(eltype, promote_type, ops)}}, ops)
+    MT = Matrix{Operator{mapreduce(eltype, promote_type, ops)}}
+    opsm = strictconvert(MT, ops)
     InterlaceOperator(opsm,ds,rs,
                         cache(dsi),
                         cache(rsi),
                         (l,u))
 end
 
 
-function InterlaceOperator(ops::Vector{<:Operator},ds::Space,rs::Space)
+function InterlaceOperator(ops::VectorOrTupleOfOp, ds::Space, rs::Space)
     # calculate bandwidths
     p=size(ops,1)
     if all(isbanded,ops)
@@ -134,7 +134,8 @@ function InterlaceOperator(ops::Vector{<:Operator},ds::Space,rs::Space)
         l,u = (1-dimension(rs),dimension(ds)-1)  # not banded
     end
 
-    opsv = convert(Vector{Operator{mapreduce(eltype, promote_type, ops)}}, ops)
+    VT = Vector{Operator{mapreduce(eltype, promote_type, ops)}}
+    opsv = strictconvert(VT, _convert_vector(ops))
     InterlaceOperator(opsv,ds,rs,
                         cache(BlockInterlacer(tuple(blocklengths(ds)))),
                         cache(interlacer(rs)),
@@ -176,18 +177,18 @@ InterlaceOperator(opsin::AbstractMatrix{<:Operator},::Type{DS}) where {DS<:Space
 InterlaceOperator(opsin::AbstractMatrix,S...) =
     InterlaceOperator(Matrix{Operator{mapreduce(eltype,promote_type,opsin)}}(promotespaces(opsin)),S...)
 
-_convertVector(v::AbstractVector) = convert(Vector, v)
-_convertVector(t::Tuple) = [t...]
+_convert_vector(v::AbstractVector) = convert(Vector, v)
+_convert_vector_or_svector(v::AbstractVector) = _convert_vector(v)
+_convert_vector(t::Tuple) = [t...]
+_convert_vector_or_svector(t::Tuple) = SVector{length(t), mapreduce(typeof, typejoin, t)}(t)
 
 function InterlaceOperator(opsin::AbstractVector{<:Operator})
-    ops = promotedomainspace(opsin)
-    opsv = _convertVector(ops)
-    InterlaceOperator(opsv, domainspace(first(ops)), rangespace(opsv))
+    ops = _convert_vector(promotedomainspace(opsin))
+    InterlaceOperator(ops, domainspace(first(ops)), rangespace(ops))
 end
 function InterlaceOperator(opsin::Tuple{Operator, Vararg{Operator}})
     ops = promotedomainspace(opsin)
-    opsv = _convertVector(ops)
-    InterlaceOperator(opsv, domainspace(first(ops)), rangespace(opsv))
+    InterlaceOperator(ops, domainspace(first(ops)), rangespace(ops))
 end
 
 InterlaceOperator(ops::AbstractArray) =

src/Operators/spacepromotion.jl

@@ -86,8 +86,6 @@ promotedomainspace(P::Operator,sp::Space,cursp::Space) =
 
 
 
-const VectorOrTupleOfOp{O<:Operator} = Union{AbstractVector{O}, Tuple{O, Vararg{O}}}
-
 __maybetypedmap(f, k, ops) = map(op->f(op,k), ops)
 _maybetypedmap(f, k, O, ops::AbstractVector) =
     strictconvert(Vector{O}, __maybetypedmap(f, k, ops))

src/Spaces/ArraySpace.jl

@@ -10,28 +10,26 @@ f = Fun(x->[exp(x),sin(x)],-1..1)
 space(f) == ArraySpace(Chebyshev(),2)
 ```
 """
-struct ArraySpace{S,n,DD,RR} <: DirectSumSpace{NTuple{n,S},DD,Array{RR,n}}
-     spaces::Array{S,n}
+struct ArraySpace{S,n,DD,RR,A<:AbstractArray{S,n}} <: DirectSumSpace{NTuple{n,S},DD,Array{RR,n}}
+     spaces::A
 end
 
-const VectorSpace{S,DD,RR} = ArraySpace{S,1,DD,RR}
-const MatrixSpace{S,DD,RR} = ArraySpace{S,2,DD,RR}
+const VectorSpace{S,DD,RR,A<:AbstractVector{S}} = ArraySpace{S,1,DD,RR,A}
+const MatrixSpace{S,DD,RR,A<:AbstractMatrix{S}} = ArraySpace{S,2,DD,RR,A}
 
 #TODO: Think through domain/domaindominsion
-ArraySpace(sp::AbstractArray{SS,N}) where {SS<:Space,N} =
-    ArraySpace{SS,N,domaintype(first(sp)),mapreduce(rangetype,promote_type,sp)}(sp)
+ArraySpace(sp::AbstractArray{SS,N}, f = first(sp)) where {SS<:Space,N} =
+    ArraySpace{SS,N,domaintype(f),mapreduce(rangetype,promote_type,sp),typeof(sp)}(sp)
 ArraySpace(S::Space,n::NTuple{N,Int}) where {N} = ArraySpace(fill(S,n))
 ArraySpace(S::Space,n::Integer) = ArraySpace(S,(n,))
 ArraySpace(S::Space,n,m) = ArraySpace(fill(S,(n,m)))
 ArraySpace(d::Domain,n...) = ArraySpace(Space(d),n...)
 
 Space(sp::AbstractArray{<:Space}) = ArraySpace(sp)
-convert(::Type{Array}, sp::ArraySpace) = sp.spaces
-convert(::Type{Vector}, sp::VectorSpace) = sp.spaces
-convert(::Type{Matrix}, sp::MatrixSpace) = sp.spaces
-Array(sp::ArraySpace) = sp.spaces
-Vector(sp::VectorSpace) = sp.spaces
-Matrix(sp::MatrixSpace) = sp.spaces
+convert(::Type{A}, sp::ArraySpace) where {A<:Array} = A(sp.spaces)::A
+Array(sp::ArraySpace) = convert(Array, sp.spaces)
+Vector(sp::VectorSpace) = convert(Vector, sp.spaces)
+Matrix(sp::MatrixSpace) = convert(Matrix, sp.spaces)
 
 
 BlockInterlacer(sp::ArraySpace) = BlockInterlacer(blocklengths.(tuple(sp.spaces...)))
@@ -107,8 +105,7 @@ end
 
 
 Base.vec(AS::ArraySpace) = ArraySpace(vec(AS.spaces))
-Base.vec(f::Fun{ArraySpace{S,n,DD,RR}}) where {S,n,DD,RR} =
-    [f[j] for j=1:length(f.space)]
+Base.vec(f::Fun{<:ArraySpace}) = [f[j] for j=1:length(f.space)]
 
 repeat(A::ArraySpace,n,m) = ArraySpace(repeat(A.spaces,n,m))
 
@@ -117,14 +114,14 @@ component(A::MatrixSpace,k::Integer,j::Integer) = A.spaces[k,j]
 Base.getindex(f::Fun{DSS},k::Integer) where {DSS<:ArraySpace} = component(f,k)
 
 
-Base.getindex(f::Fun{MatrixSpace{S,DD,RR}},k::Integer,j::Integer) where {S,DD,RR} =
+Base.getindex(f::Fun{<:MatrixSpace},k::Integer,j::Integer) =
     f[k+stride(f,2)*(j-1)]
 
 Base.getindex(f::Fun{DSS},kj::CartesianIndex{1}) where {DSS<:ArraySpace} = f[kj[1]]
 Base.getindex(f::Fun{DSS},kj::CartesianIndex{2}) where {DSS<:ArraySpace} = f[kj[1],kj[2]]
 
 
-function Fun(A::AbstractArray{Fun{VectorSpace{S,DD,RR},V,VV},2}) where {S,V,VV,DD,RR}
+function Fun(A::AbstractMatrix{<:Fun{<:VectorSpace{S},V,VV}}) where {S,V,VV}
     @assert size(A,1)==1
 
     M = Matrix{Fun{S,V,VV}}(undef, length(space(A[1])),size(A,2))

src/Spaces/ProductSpaceOperators.jl

@@ -65,7 +65,7 @@ Evaluation(S::SumSpace,x,order) =
         InterlaceOperator(RowVector(vnocat(map(s->Evaluation(s,x,order),components(S))...)),SumSpace))
 
 
-ToeplitzOperator(G::Fun{MatrixSpace{S,DD,RR},V}) where {S,RR,V,DD} = interlace(map(ToeplitzOperator,Array(G)))
+ToeplitzOperator(G::Fun{<:MatrixSpace}) = interlace(map(ToeplitzOperator,Array(G)))
 
 ## Sum Space
 
@@ -253,10 +253,14 @@ choosedomainspace(M::CalculusOperator{UnsetSpace},sp::SumSpace)=mapreduce(s->cho
 
 ## Multiplcation for Array*Vector
 
-function Multiplication(f::Fun{MatrixSpace{S,DD,RR}},sp::VectorSpace{S2,DD2,RR2}) where {S,DD,RR,S2,DD2,RR2}
+function Multiplication(f::Fun{<:MatrixSpace}, sp::VectorSpace)
     @assert size(space(f),2)==length(sp)
     m=Array(f)
-    MultiplicationWrapper(f,interlace(Operator{promote_type(cfstype(f),prectype(sp))}[Multiplication(m[k,j],sp[j]) for k=1:size(m,1),j=1:size(m,2)]))
+    MultiplicationWrapper(f, interlace(
+            Operator{promote_type(cfstype(f),prectype(sp))}[
+                Multiplication(m[k,j],sp[j]) for k=1:size(m,1),j=1:size(m,2)]
+        )
+    )
 end