Mercurial > repos > public > sbplib_julia
diff src/LazyTensors/lazy_tensor_operations.jl @ 1023:52f07c77299d refactor/sbpoperators/inflation
Merge refactor/lazy_tensors
author | Jonatan Werpers <jonatan@werpers.com> |
---|---|
date | Mon, 21 Mar 2022 09:51:07 +0100 |
parents | bbbc31953367 f7a718bcb4da |
children | f857057e61e6 |
line wrap: on
line diff
--- a/src/LazyTensors/lazy_tensor_operations.jl Fri Mar 18 16:57:00 2022 +0100 +++ b/src/LazyTensors/lazy_tensor_operations.jl Mon Mar 21 09:51:07 2022 +0100 @@ -1,204 +1,137 @@ """ - LazyTensorMappingApplication{T,R,D} <: LazyArray{T,R} + LazyTensorApplication{T,R,D} <: LazyArray{T,R} -Struct for lazy application of a TensorMapping. Created using `*`. +Struct for lazy application of a LazyTensor. Created using `*`. -Allows the result of a `TensorMapping` applied to a vector to be treated as an `AbstractArray`. -With a mapping `m` and a vector `v` the LazyTensorMappingApplication object can be created by `m*v`. +Allows the result of a `LazyTensor` applied to a vector to be treated as an `AbstractArray`. +With a mapping `m` and a vector `v` the LazyTensorApplication object can be created by `m*v`. The actual result will be calcualted when indexing into `m*v`. """ -struct LazyTensorMappingApplication{T,R,D, TM<:TensorMapping{<:Any,R,D}, AA<:AbstractArray{<:Any,D}} <: LazyArray{T,R} +struct LazyTensorApplication{T,R,D, TM<:LazyTensor{<:Any,R,D}, AA<:AbstractArray{<:Any,D}} <: LazyArray{T,R} t::TM o::AA - function LazyTensorMappingApplication(t::TensorMapping{<:Any,R,D}, o::AbstractArray{<:Any,D}) where {R,D} + function LazyTensorApplication(t::LazyTensor{<:Any,R,D}, o::AbstractArray{<:Any,D}) where {R,D} + @boundscheck check_domain_size(t, size(o)) I = ntuple(i->1, range_dim(t)) T = typeof(apply(t,o,I...)) return new{T,R,D,typeof(t), typeof(o)}(t,o) end end -# TODO: Do boundschecking on creation! - -Base.getindex(ta::LazyTensorMappingApplication{T,R}, I::Vararg{Any,R}) where {T,R} = apply(ta.t, ta.o, I...) -Base.getindex(ta::LazyTensorMappingApplication{T,1}, I::CartesianIndex{1}) where {T} = apply(ta.t, ta.o, I.I...) # Would otherwise be caught in the previous method. -Base.size(ta::LazyTensorMappingApplication) = range_size(ta.t) -# TODO: What else is needed to implement the AbstractArray interface? -Base.:*(a::TensorMapping, v::AbstractArray) = LazyTensorMappingApplication(a,v) -Base.:*(a::TensorMapping, b::TensorMapping) = throw(MethodError(Base.:*,(a,b))) -Base.:*(a::TensorMapping, args::Union{TensorMapping, AbstractArray}...) = foldr(*,(a,args...)) +function Base.getindex(ta::LazyTensorApplication{T,R}, I::Vararg{Any,R}) where {T,R} + @boundscheck checkbounds(ta, Int.(I)...) + return apply(ta.t, ta.o, I...) +end +Base.getindex(ta::LazyTensorApplication{T,1} where T, I::CartesianIndex{1}) = ta[Tuple(I)...] # Would otherwise be caught in the previous method. +Base.size(ta::LazyTensorApplication) = range_size(ta.t) -# # We need the associativity to be a→b→c = a→(b→c), which is the case for '→' -# # Should we overload some other infix binary opesrator? -# →(tm::TensorMapping{T,R,D}, o::AbstractArray{T,D}) where {T,R,D} = LazyTensorMappingApplication(tm,o) -# TODO: We need to be really careful about good error messages. -# For example what happens if you try to multiply LazyTensorMappingApplication with a TensorMapping(wrong order)? """ - LazyTensorMappingTranspose{T,R,D} <: TensorMapping{T,D,R} + LazyTensorTranspose{T,R,D} <: LazyTensor{T,D,R} -Struct for lazy transpose of a TensorMapping. +Struct for lazy transpose of a LazyTensor. If a mapping implements the the `apply_transpose` method this allows working with -the transpose of mapping `m` by using `m'`. `m'` will work as a regular TensorMapping lazily calling +the transpose of mapping `m` by using `m'`. `m'` will work as a regular LazyTensor lazily calling the appropriate methods of `m`. """ -struct LazyTensorMappingTranspose{T,R,D, TM<:TensorMapping{T,R,D}} <: TensorMapping{T,D,R} +struct LazyTensorTranspose{T,R,D, TM<:LazyTensor{T,R,D}} <: LazyTensor{T,D,R} tm::TM end # # TBD: Should this be implemented on a type by type basis or through a trait to provide earlier errors? -# Jonatan 2020-09-25: Is the problem that you can take the transpose of any TensorMapping even if it doesn't implement `apply_transpose`? -Base.adjoint(tm::TensorMapping) = LazyTensorMappingTranspose(tm) -Base.adjoint(tmt::LazyTensorMappingTranspose) = tmt.tm +# Jonatan 2020-09-25: Is the problem that you can take the transpose of any LazyTensor even if it doesn't implement `apply_transpose`? +Base.adjoint(tm::LazyTensor) = LazyTensorTranspose(tm) +Base.adjoint(tmt::LazyTensorTranspose) = tmt.tm -apply(tmt::LazyTensorMappingTranspose{T,R,D}, v::AbstractArray{<:Any,R}, I::Vararg{Any,D}) where {T,R,D} = apply_transpose(tmt.tm, v, I...) -apply_transpose(tmt::LazyTensorMappingTranspose{T,R,D}, v::AbstractArray{<:Any,D}, I::Vararg{Any,R}) where {T,R,D} = apply(tmt.tm, v, I...) +apply(tmt::LazyTensorTranspose{T,R,D}, v::AbstractArray{<:Any,R}, I::Vararg{Any,D}) where {T,R,D} = apply_transpose(tmt.tm, v, I...) +apply_transpose(tmt::LazyTensorTranspose{T,R,D}, v::AbstractArray{<:Any,D}, I::Vararg{Any,R}) where {T,R,D} = apply(tmt.tm, v, I...) -range_size(tmt::LazyTensorMappingTranspose) = domain_size(tmt.tm) -domain_size(tmt::LazyTensorMappingTranspose) = range_size(tmt.tm) +range_size(tmt::LazyTensorTranspose) = domain_size(tmt.tm) +domain_size(tmt::LazyTensorTranspose) = range_size(tmt.tm) -struct LazyTensorMappingBinaryOperation{Op,T,R,D,T1<:TensorMapping{T,R,D},T2<:TensorMapping{T,R,D}} <: TensorMapping{T,D,R} +struct LazyTensorBinaryOperation{Op,T,R,D,T1<:LazyTensor{T,R,D},T2<:LazyTensor{T,R,D}} <: LazyTensor{T,D,R} tm1::T1 tm2::T2 - @inline function LazyTensorMappingBinaryOperation{Op,T,R,D}(tm1::T1,tm2::T2) where {Op,T,R,D, T1<:TensorMapping{T,R,D},T2<:TensorMapping{T,R,D}} + function LazyTensorBinaryOperation{Op,T,R,D}(tm1::T1,tm2::T2) where {Op,T,R,D, T1<:LazyTensor{T,R,D},T2<:LazyTensor{T,R,D}} + @boundscheck check_domain_size(tm2, domain_size(tm1)) + @boundscheck check_range_size(tm2, range_size(tm1)) return new{Op,T,R,D,T1,T2}(tm1,tm2) end end -# TODO: Boundschecking in constructor. -apply(tmBinOp::LazyTensorMappingBinaryOperation{:+,T,R,D}, v::AbstractArray{<:Any,D}, I::Vararg{Any,R}) where {T,R,D} = apply(tmBinOp.tm1, v, I...) + apply(tmBinOp.tm2, v, I...) -apply(tmBinOp::LazyTensorMappingBinaryOperation{:-,T,R,D}, v::AbstractArray{<:Any,D}, I::Vararg{Any,R}) where {T,R,D} = apply(tmBinOp.tm1, v, I...) - apply(tmBinOp.tm2, v, I...) +LazyTensorBinaryOperation{Op}(s,t) where Op = LazyTensorBinaryOperation{Op,eltype(s), range_dim(s), domain_dim(s)}(s,t) -range_size(tmBinOp::LazyTensorMappingBinaryOperation) = range_size(tmBinOp.tm1) -domain_size(tmBinOp::LazyTensorMappingBinaryOperation) = domain_size(tmBinOp.tm1) +apply(tmBinOp::LazyTensorBinaryOperation{:+,T,R,D}, v::AbstractArray{<:Any,D}, I::Vararg{Any,R}) where {T,R,D} = apply(tmBinOp.tm1, v, I...) + apply(tmBinOp.tm2, v, I...) +apply(tmBinOp::LazyTensorBinaryOperation{:-,T,R,D}, v::AbstractArray{<:Any,D}, I::Vararg{Any,R}) where {T,R,D} = apply(tmBinOp.tm1, v, I...) - apply(tmBinOp.tm2, v, I...) -Base.:+(tm1::TensorMapping{T,R,D}, tm2::TensorMapping{T,R,D}) where {T,R,D} = LazyTensorMappingBinaryOperation{:+,T,R,D}(tm1,tm2) -Base.:-(tm1::TensorMapping{T,R,D}, tm2::TensorMapping{T,R,D}) where {T,R,D} = LazyTensorMappingBinaryOperation{:-,T,R,D}(tm1,tm2) +range_size(tmBinOp::LazyTensorBinaryOperation) = range_size(tmBinOp.tm1) +domain_size(tmBinOp::LazyTensorBinaryOperation) = domain_size(tmBinOp.tm1) + """ - TensorMappingComposition{T,R,K,D} + LazyTensorComposition{T,R,K,D} -Lazily compose two `TensorMapping`s, so that they can be handled as a single `TensorMapping`. +Lazily compose two `LazyTensor`s, so that they can be handled as a single `LazyTensor`. """ -struct TensorMappingComposition{T,R,K,D, TM1<:TensorMapping{T,R,K}, TM2<:TensorMapping{T,K,D}} <: TensorMapping{T,R,D} +struct LazyTensorComposition{T,R,K,D, TM1<:LazyTensor{T,R,K}, TM2<:LazyTensor{T,K,D}} <: LazyTensor{T,R,D} t1::TM1 t2::TM2 - @inline function TensorMappingComposition(t1::TensorMapping{T,R,K}, t2::TensorMapping{T,K,D}) where {T,R,K,D} + function LazyTensorComposition(t1::LazyTensor{T,R,K}, t2::LazyTensor{T,K,D}) where {T,R,K,D} @boundscheck check_domain_size(t1, range_size(t2)) return new{T,R,K,D, typeof(t1), typeof(t2)}(t1,t2) end end -range_size(tm::TensorMappingComposition) = range_size(tm.t1) -domain_size(tm::TensorMappingComposition) = domain_size(tm.t2) +range_size(tm::LazyTensorComposition) = range_size(tm.t1) +domain_size(tm::LazyTensorComposition) = domain_size(tm.t2) -function apply(c::TensorMappingComposition{T,R,K,D}, v::AbstractArray{<:Any,D}, I::Vararg{Any,R}) where {T,R,K,D} +function apply(c::LazyTensorComposition{T,R,K,D}, v::AbstractArray{<:Any,D}, I::Vararg{Any,R}) where {T,R,K,D} apply(c.t1, c.t2*v, I...) end -function apply_transpose(c::TensorMappingComposition{T,R,K,D}, v::AbstractArray{<:Any,R}, I::Vararg{Any,D}) where {T,R,K,D} +function apply_transpose(c::LazyTensorComposition{T,R,K,D}, v::AbstractArray{<:Any,R}, I::Vararg{Any,D}) where {T,R,K,D} apply_transpose(c.t2, c.t1'*v, I...) end -Base.@propagate_inbounds Base.:∘(s::TensorMapping, t::TensorMapping) = TensorMappingComposition(s,t) - -""" - LazyLinearMap{T,R,D,...}(A, range_indicies, domain_indicies) - -TensorMapping defined by the AbstractArray A. `range_indicies` and `domain_indicies` define which indicies of A should -be considerd the range and domain of the TensorMapping. Each set of indices must be ordered in ascending order. - -For instance, if A is a m x n matrix, and range_size = (1,), domain_size = (2,), then the LazyLinearMap performs the -standard matrix-vector product on vectors of size n. -""" -struct LazyLinearMap{T,R,D, RD, AA<:AbstractArray{T,RD}} <: TensorMapping{T,R,D} - A::AA - range_indicies::NTuple{R,Int} - domain_indicies::NTuple{D,Int} - - function LazyLinearMap(A::AA, range_indicies::NTuple{R,Int}, domain_indicies::NTuple{D,Int}) where {T,R,D, RD, AA<:AbstractArray{T,RD}} - if !issorted(range_indicies) || !issorted(domain_indicies) - throw(DomainError("range_indicies and domain_indicies must be sorted in ascending order")) - end - - return new{T,R,D,RD,AA}(A,range_indicies,domain_indicies) - end -end - -range_size(llm::LazyLinearMap) = size(llm.A)[[llm.range_indicies...]] -domain_size(llm::LazyLinearMap) = size(llm.A)[[llm.domain_indicies...]] - -function apply(llm::LazyLinearMap{T,R,D}, v::AbstractArray{<:Any,D}, I::Vararg{Any,R}) where {T,R,D} - view_index = ntuple(i->:,ndims(llm.A)) - for i ∈ 1:R - view_index = Base.setindex(view_index, Int(I[i]), llm.range_indicies[i]) - end - A_view = @view llm.A[view_index...] - return sum(A_view.*v) -end - -function apply_transpose(llm::LazyLinearMap{T,R,D}, v::AbstractArray{<:Any,R}, I::Vararg{Any,D}) where {T,R,D} - apply(LazyLinearMap(llm.A, llm.domain_indicies, llm.range_indicies), v, I...) -end - """ - IdentityMapping{T,D} <: TensorMapping{T,D,D} - -The lazy identity TensorMapping for a given size. Usefull for building up higher dimensional tensor mappings from lower -dimensional ones through outer products. Also used in the Implementation for InflatedTensorMapping. -""" -struct IdentityMapping{T,D} <: TensorMapping{T,D,D} - size::NTuple{D,Int} -end - -IdentityMapping{T}(size::NTuple{D,Int}) where {T,D} = IdentityMapping{T,D}(size) -IdentityMapping{T}(size::Vararg{Int,D}) where {T,D} = IdentityMapping{T,D}(size) -IdentityMapping(size::Vararg{Int,D}) where D = IdentityMapping{Float64,D}(size) + LazyTensorComposition(tm, tmi::IdentityTensor) + LazyTensorComposition(tmi::IdentityTensor, tm) -range_size(tmi::IdentityMapping) = tmi.size -domain_size(tmi::IdentityMapping) = tmi.size - -apply(tmi::IdentityMapping{T,D}, v::AbstractArray{<:Any,D}, I::Vararg{Any,D}) where {T,D} = v[I...] -apply_transpose(tmi::IdentityMapping{T,D}, v::AbstractArray{<:Any,D}, I::Vararg{Any,D}) where {T,D} = v[I...] - +Composes a `Tensormapping` `tm` with an `IdentityTensor` `tmi`, by returning `tm` """ - Base.:∘(tm, tmi) - Base.:∘(tmi, tm) - -Composes a `Tensormapping` `tm` with an `IdentityMapping` `tmi`, by returning `tm` -""" -@inline function Base.:∘(tm::TensorMapping{T,R,D}, tmi::IdentityMapping{T,D}) where {T,R,D} +function LazyTensorComposition(tm::LazyTensor{T,R,D}, tmi::IdentityTensor{T,D}) where {T,R,D} @boundscheck check_domain_size(tm, range_size(tmi)) return tm end -@inline function Base.:∘(tmi::IdentityMapping{T,R}, tm::TensorMapping{T,R,D}) where {T,R,D} +function LazyTensorComposition(tmi::IdentityTensor{T,R}, tm::LazyTensor{T,R,D}) where {T,R,D} @boundscheck check_domain_size(tmi, range_size(tm)) return tm end -# Specialization for the case where tm is an IdentityMapping. Required to resolve ambiguity. -@inline function Base.:∘(tm::IdentityMapping{T,D}, tmi::IdentityMapping{T,D}) where {T,D} +# Specialization for the case where tm is an IdentityTensor. Required to resolve ambiguity. +function LazyTensorComposition(tm::IdentityTensor{T,D}, tmi::IdentityTensor{T,D}) where {T,D} @boundscheck check_domain_size(tm, range_size(tmi)) return tmi end """ - InflatedTensorMapping{T,R,D} <: TensorMapping{T,R,D} + InflatedLazyTensor{T,R,D} <: LazyTensor{T,R,D} -An inflated `TensorMapping` with dimensions added before and afer its actual dimensions. +An inflated `LazyTensor` with dimensions added before and afer its actual dimensions. """ -struct InflatedTensorMapping{T,R,D,D_before,R_middle,D_middle,D_after, TM<:TensorMapping{T,R_middle,D_middle}} <: TensorMapping{T,R,D} - before::IdentityMapping{T,D_before} +struct InflatedLazyTensor{T,R,D,D_before,R_middle,D_middle,D_after, TM<:LazyTensor{T,R_middle,D_middle}} <: LazyTensor{T,R,D} + before::IdentityTensor{T,D_before} tm::TM - after::IdentityMapping{T,D_after} + after::IdentityTensor{T,D_after} - function InflatedTensorMapping(before, tm::TensorMapping{T}, after) where T + function InflatedLazyTensor(before, tm::LazyTensor{T}, after) where T R_before = range_dim(before) R_middle = range_dim(tm) R_after = range_dim(after) @@ -211,36 +144,37 @@ return new{T,R,D,D_before,R_middle,D_middle,D_after, typeof(tm)}(before, tm, after) end end + """ - InflatedTensorMapping(before, tm, after) - InflatedTensorMapping(before,tm) - InflatedTensorMapping(tm,after) + InflatedLazyTensor(before, tm, after) + InflatedLazyTensor(before,tm) + InflatedLazyTensor(tm,after) -The outer product of `before`, `tm` and `after`, where `before` and `after` are `IdentityMapping`s. +The outer product of `before`, `tm` and `after`, where `before` and `after` are `IdentityTensor`s. -If one of `before` or `after` is left out, a 0-dimensional `IdentityMapping` is used as the default value. +If one of `before` or `after` is left out, a 0-dimensional `IdentityTensor` is used as the default value. -If `tm` already is an `InflatedTensorMapping`, `before` and `after` will be extended instead of -creating a nested `InflatedTensorMapping`. +If `tm` already is an `InflatedLazyTensor`, `before` and `after` will be extended instead of +creating a nested `InflatedLazyTensor`. """ -InflatedTensorMapping(::IdentityMapping, ::TensorMapping, ::IdentityMapping) +InflatedLazyTensor(::IdentityTensor, ::LazyTensor, ::IdentityTensor) -function InflatedTensorMapping(before, itm::InflatedTensorMapping, after) - return InflatedTensorMapping( - IdentityMapping(before.size..., itm.before.size...), +function InflatedLazyTensor(before, itm::InflatedLazyTensor, after) + return InflatedLazyTensor( + IdentityTensor(before.size..., itm.before.size...), itm.tm, - IdentityMapping(itm.after.size..., after.size...), + IdentityTensor(itm.after.size..., after.size...), ) end -InflatedTensorMapping(before::IdentityMapping, tm::TensorMapping{T}) where T = InflatedTensorMapping(before,tm,IdentityMapping{T}()) -InflatedTensorMapping(tm::TensorMapping{T}, after::IdentityMapping) where T = InflatedTensorMapping(IdentityMapping{T}(),tm,after) +InflatedLazyTensor(before::IdentityTensor, tm::LazyTensor{T}) where T = InflatedLazyTensor(before,tm,IdentityTensor{T}()) +InflatedLazyTensor(tm::LazyTensor{T}, after::IdentityTensor) where T = InflatedLazyTensor(IdentityTensor{T}(),tm,after) # Resolve ambiguity between the two previous methods -InflatedTensorMapping(I1::IdentityMapping{T}, I2::IdentityMapping{T}) where T = InflatedTensorMapping(I1,I2,IdentityMapping{T}()) +InflatedLazyTensor(I1::IdentityTensor{T}, I2::IdentityTensor{T}) where T = InflatedLazyTensor(I1,I2,IdentityTensor{T}()) -# TODO: Implement some pretty printing in terms of ⊗. E.g InflatedTensorMapping(I(3),B,I(2)) -> I(3)⊗B⊗I(2) +# TODO: Implement some pretty printing in terms of ⊗. E.g InflatedLazyTensor(I(3),B,I(2)) -> I(3)⊗B⊗I(2) -function range_size(itm::InflatedTensorMapping) +function range_size(itm::InflatedLazyTensor) return flatten_tuple( range_size(itm.before), range_size(itm.tm), @@ -248,7 +182,7 @@ ) end -function domain_size(itm::InflatedTensorMapping) +function domain_size(itm::InflatedLazyTensor) return flatten_tuple( domain_size(itm.before), domain_size(itm.tm), @@ -256,7 +190,7 @@ ) end -function apply(itm::InflatedTensorMapping{T,R,D}, v::AbstractArray{<:Any,D}, I::Vararg{Any,R}) where {T,R,D} +function apply(itm::InflatedLazyTensor{T,R,D}, v::AbstractArray{<:Any,D}, I::Vararg{Any,R}) where {T,R,D} dim_before = range_dim(itm.before) dim_domain = domain_dim(itm.tm) dim_range = range_dim(itm.tm) @@ -268,7 +202,7 @@ return apply(itm.tm, v_inner, inner_index...) end -function apply_transpose(itm::InflatedTensorMapping{T,R,D}, v::AbstractArray{<:Any,R}, I::Vararg{Any,D}) where {T,R,D} +function apply_transpose(itm::InflatedLazyTensor{T,R,D}, v::AbstractArray{<:Any,R}, I::Vararg{Any,D}) where {T,R,D} dim_before = range_dim(itm.before) dim_domain = domain_dim(itm.tm) dim_range = range_dim(itm.tm) @@ -281,87 +215,10 @@ end -""" - split_index(::Val{dim_before}, ::Val{dim_view}, ::Val{dim_index}, ::Val{dim_after}, I...) - -Splits the multi-index `I` into two parts. One part which is expected to be -used as a view, and one which is expected to be used as an index. -Eg. -``` -split_index(Val(1),Val(3),Val(2),Val(1),(1,2,3,4)) -> (1,:,:,:,4), (2,3) -``` - -`dim_view` controls how many colons are in the view, and `dim_index` controls -how many elements are extracted from the middle. -`dim_before` and `dim_after` decides the length of the index parts before and after the colons in the view index. - -Arguments should satisfy `length(I) == dim_before+B_domain+dim_after`. - -The returned values satisfy - * `length(view_index) == dim_before + dim_view + dim_after` - * `length(I_middle) == dim_index` -""" -function split_index(::Val{dim_before}, ::Val{dim_view}, ::Val{dim_index}, ::Val{dim_after}, I...) where {dim_before,dim_view, dim_index,dim_after} - I_before, I_middle, I_after = split_tuple(I, Val(dim_before), Val(dim_index)) - - view_index = (I_before..., ntuple((i)->:, dim_view)..., I_after...) - - return view_index, I_middle -end - -# TODO: Can this be replaced by something more elegant while still being type stable? 2020-10-21 -# See: -# https://github.com/JuliaLang/julia/issues/34884 -# https://github.com/JuliaLang/julia/issues/30386 -""" - slice_tuple(t, Val(l), Val(u)) - -Get a slice of a tuple in a type stable way. -Equivalent to `t[l:u]` but type stable. -""" -function slice_tuple(t,::Val{L},::Val{U}) where {L,U} - return ntuple(i->t[i+L-1], U-L+1) -end - -""" - split_tuple(t::Tuple{...}, ::Val{M}) where {N,M} - -Split the tuple `t` into two parts. the first part is `M` long. -E.g -```julia -split_tuple((1,2,3,4),Val(3)) -> (1,2,3), (4,) -``` -""" -function split_tuple(t::NTuple{N,Any},::Val{M}) where {N,M} - return slice_tuple(t,Val(1), Val(M)), slice_tuple(t,Val(M+1), Val(N)) -end - -""" - split_tuple(t::Tuple{...},::Val{M},::Val{K}) where {N,M,K} - -Same as `split_tuple(t::NTuple{N},::Val{M})` but splits the tuple in three parts. With the first -two parts having lenght `M` and `K`. -""" -function split_tuple(t::NTuple{N,Any},::Val{M},::Val{K}) where {N,M,K} - p1, tail = split_tuple(t, Val(M)) - p2, p3 = split_tuple(tail, Val(K)) - return p1,p2,p3 -end - - -""" - flatten_tuple(t) - -Takes a nested tuple and flattens the whole structure -""" -flatten_tuple(t::NTuple{N, Number} where N) = t -flatten_tuple(t::Tuple) = ((flatten_tuple.(t)...)...,) # simplify? -flatten_tuple(ts::Vararg) = flatten_tuple(ts) - @doc raw""" LazyOuterProduct(tms...) -Creates a `TensorMappingComposition` for the outerproduct of `tms...`. +Creates a `LazyTensorComposition` for the outerproduct of `tms...`. This is done by separating the outer product into regular products of outer products involving only identity mappings and one non-identity mapping. First let @@ -397,20 +254,19 @@ """ function LazyOuterProduct end -function LazyOuterProduct(tm1::TensorMapping{T}, tm2::TensorMapping{T}) where T - itm1 = InflatedTensorMapping(tm1, IdentityMapping{T}(range_size(tm2))) - itm2 = InflatedTensorMapping(IdentityMapping{T}(domain_size(tm1)),tm2) +function LazyOuterProduct(tm1::LazyTensor{T}, tm2::LazyTensor{T}) where T + itm1 = InflatedLazyTensor(tm1, IdentityTensor{T}(range_size(tm2))) + itm2 = InflatedLazyTensor(IdentityTensor{T}(domain_size(tm1)),tm2) return itm1∘itm2 end -LazyOuterProduct(t1::IdentityMapping{T}, t2::IdentityMapping{T}) where T = IdentityMapping{T}(t1.size...,t2.size...) -LazyOuterProduct(t1::TensorMapping, t2::IdentityMapping) = InflatedTensorMapping(t1, t2) -LazyOuterProduct(t1::IdentityMapping, t2::TensorMapping) = InflatedTensorMapping(t1, t2) +LazyOuterProduct(t1::IdentityTensor{T}, t2::IdentityTensor{T}) where T = IdentityTensor{T}(t1.size...,t2.size...) +LazyOuterProduct(t1::LazyTensor, t2::IdentityTensor) = InflatedLazyTensor(t1, t2) +LazyOuterProduct(t1::IdentityTensor, t2::LazyTensor) = InflatedLazyTensor(t1, t2) -LazyOuterProduct(tms::Vararg{TensorMapping}) = foldl(LazyOuterProduct, tms) +LazyOuterProduct(tms::Vararg{LazyTensor}) = foldl(LazyOuterProduct, tms) -⊗(a::TensorMapping, b::TensorMapping) = LazyOuterProduct(a,b) """ @@ -420,24 +276,41 @@ # TODO: Describe when it is useful """ -function inflate(tm::TensorMapping, sz, dir) - Is = IdentityMapping{eltype(tm)}.(sz) +function inflate(tm::LazyTensor, sz, dir) + Is = IdentityTensor{eltype(tm)}.(sz) parts = Base.setindex(Is, tm, dir) return foldl(⊗, parts) end -function check_domain_size(tm::TensorMapping, sz) +function check_domain_size(tm::LazyTensor, sz) if domain_size(tm) != sz - throw(SizeMismatch(tm,sz)) + throw(DomainSizeMismatch(tm,sz)) + end +end + +function check_range_size(tm::LazyTensor, sz) + if range_size(tm) != sz + throw(RangeSizeMismatch(tm,sz)) end end -struct SizeMismatch <: Exception - tm::TensorMapping +struct DomainSizeMismatch <: Exception + tm::LazyTensor sz end -function Base.showerror(io::IO, err::SizeMismatch) - print(io, "SizeMismatch: ") - print(io, "domain size $(domain_size(err.tm)) of TensorMapping not matching size $(err.sz)") +function Base.showerror(io::IO, err::DomainSizeMismatch) + print(io, "DomainSizeMismatch: ") + print(io, "domain size $(domain_size(err.tm)) of LazyTensor not matching size $(err.sz)") end + + +struct RangeSizeMismatch <: Exception + tm::LazyTensor + sz +end + +function Base.showerror(io::IO, err::RangeSizeMismatch) + print(io, "RangeSizeMismatch: ") + print(io, "range size $(range_size(err.tm)) of LazyTensor not matching size $(err.sz)") +end