Prevent construction of multifusion tensors with incompatible coloring

src/spaces/gradedspace.jl

@@ -122,6 +122,8 @@ function flip(V::GradedSpace{I}) where {I <: Sector}
 end
 
 function unitspace(S::Type{<:GradedSpace{I}}) where {I <: Sector}
+    UnitStyle(I) isa GenericUnit &&
+        throw(ArgumentError("Cannot construct unit space for sector types with semisimple unit structure."))
     return S(unit => 1 for unit in allunits(I))
 end
 zerospace(S::Type{<:GradedSpace}) = S()

src/spaces/homspace.jl

@@ -9,6 +9,43 @@ to denote categories and their objects, and keep `HomSpace` distinct.
 struct HomSpace{S <: ElementarySpace, P1 <: CompositeSpace{S}, P2 <: CompositeSpace{S}}
     codomain::P1
     domain::P2
+    function HomSpace{S, P1, P2}(codomain::P1, domain::P2) where {S <: ElementarySpace, P1 <: CompositeSpace{S}, P2 <: CompositeSpace{S}}
+        _check_unit_compatibility(codomain, domain)
+        return new{S, P1, P2}(codomain, domain)
+    end
+end
+function HomSpace(codomain::P1, domain::P2) where {S, P1 <: CompositeSpace{S}, P2 <: CompositeSpace{S}}
+    return HomSpace{S, P1, P2}(codomain, domain)
+end
+
+function _check_unit_compatibility(codomain::ProductSpace{S}, domain::ProductSpace{S}) where {S <: ElementarySpace}
+    UnitStyle(sectortype(S)) isa GenericUnit || return nothing
+    N₁, N₂ = length(codomain), length(domain)
+    N₁ == N₂ == 0 && return nothing # one() ← one()
+
+    # product spaces themselves already check that their factors are compatible
+    if N₁ == 0 # codomain is empty, domain is non-empty
+        VdomL, VdomR = domain[1], domain[N₂]
+        _leftunit(VdomL) == _rightunit(VdomR) ||
+            throw(ArgumentError("cannot construct HomSpace: domain has incompatible left and right units"))
+        return nothing
+    elseif N₂ == 0 # domain is empty, codomain is non-empty
+        VcodL, VcodR = codomain[1], codomain[N₁]
+        _leftunit(VcodL) == _rightunit(VcodR) ||
+            throw(ArgumentError("cannot construct HomSpace: codomain has incompatible left and right units"))
+        return nothing
+    end
+
+    # codomain and domain are non-empty
+    # just need to check coupled charge compatibility
+    VcodL, VdomL = codomain[1], domain[1]
+    _leftunit(VcodL) == _leftunit(VdomL) ||
+        throw(ArgumentError("cannot construct HomSpace: codomain and domain have incompatible left units"))
+
+    VcodR, VdomR = codomain[N₁], domain[N₂]
+    _rightunit(VcodR) == _rightunit(VdomR) ||
+        throw(ArgumentError("cannot construct HomSpace: codomain and domain have incompatible right units"))
+    return nothing
 end
 
 function HomSpace(codomain::S, domain::CompositeSpace{S}) where {S <: ElementarySpace}

src/spaces/productspace.jl

@@ -7,7 +7,22 @@ Only tensor products between [`ElementarySpace`](@ref) objects of the same type
 """
 struct ProductSpace{S <: ElementarySpace, N} <: CompositeSpace{S}
     spaces::NTuple{N, S}
-    ProductSpace{S, N}(spaces::NTuple{N, S}) where {S <: ElementarySpace, N} = new{S, N}(spaces)
+    function ProductSpace{S, N}(spaces::NTuple{N, S}) where {S <: ElementarySpace, N}
+        _check_unit_compatibility(S, spaces)
+        return new{S, N}(spaces)
+    end
+end
+
+# check that every product of elementary spaces is compatible color-wise
+function _check_unit_compatibility(::Type{S}, spaces::NTuple{N, S}) where {N, S <: ElementarySpace}
+    UnitStyle(sectortype(S)) isa GenericUnit || return nothing
+    N == 0 && return nothing
+
+    @inbounds for i in 1:(N - 1)
+        _rightunit(spaces[i]) == _leftunit(spaces[i + 1]) ||
+            throw(ArgumentError(lazy"spaces $(i) and $(i + 1) have incompatible coloring"))
+    end
+    return nothing
 end
 
 function ProductSpace{S, N}(spaces::Vararg{S, N}) where {S <: ElementarySpace, N}
@@ -64,7 +79,7 @@ Base.axes(P::ProductSpace) = map(axes, P)
 Base.axes(P::ProductSpace, n::Int) = axes(P[n])
 
 dual(P::ProductSpace{<:ElementarySpace, 0}) = P
-dual(P::ProductSpace) = ProductSpace(map(dual, reverse(P)))
+dual(P::ProductSpace) = ProductSpace(map(dual, reverse(P.spaces)))
 Base.conj(P::ProductSpace{<:ElementarySpace, 0}) = P
 Base.conj(P::ProductSpace) = ProductSpace(map(conj, P))
 

src/spaces/vectorspaces.jl

@@ -134,9 +134,8 @@ Always returns `false` for spaces where `V == conj(V)`, i.e. vector spaces over
 
 Return the corresponding vector space of type `S` that represents the trivial
 one-dimensional space, i.e. the space that is isomorphic to the corresponding field.
-For vector spaces where `I = sectortype(S)` has a semi-simple unit structure
-(`UnitStyle(I) == GenericUnit()`), this returns a multi-dimensional space corresponding to all unit sectors:
-`dim(unitspace(V), s) == 1` for all `s in allunits(I)`. 
+For vector spaces where `I = sectortype(S)` has a semisimple unit structure
+(`UnitStyle(I) == GenericUnit()`), this errors. 
 
 !!! note
     `unitspace(V)`is different from `one(V)`. The latter returns the empty product space
@@ -170,15 +169,17 @@ function leftunitspace(V::ElementarySpace)
     if UnitStyle(I) isa SimpleUnit
         return unitspace(typeof(V))
     else
-        !isempty(sectors(V)) || throw(ArgumentError("Cannot determine the left unit of an empty space"))
-        _allequal(leftunit, sectors(V)) ||
-            throw(ArgumentError(lazy"sectors of $V do not have the same left unit"))
-
-        sector = leftunit(first(sectors(V)))
-        return spacetype(V)(sector => 1)
+        return spacetype(V)(_leftunit(V) => 1)
     end
 end
 
+function _leftunit(V::ElementarySpace)
+    !isempty(sectors(V)) || throw(ArgumentError("Cannot determine the left unit of an empty space"))
+    _allequal(leftunit, sectors(V)) ||
+        throw(ArgumentError(lazy"sectors of $V do not have the same left unit"))
+    return leftunit(first(sectors(V)))
+end
+
 """
     rightunitspace(V::S) where {S <: ElementarySpace} -> S
 
@@ -192,22 +193,24 @@ function rightunitspace(V::ElementarySpace)
     if UnitStyle(I) isa SimpleUnit
         return unitspace(typeof(V))
     else
-        !isempty(sectors(V)) || throw(ArgumentError("Cannot determine the right unit of an empty space"))
-        _allequal(rightunit, sectors(V)) ||
-            throw(ArgumentError(lazy"sectors of $V do not have the same right unit"))
-
-        sector = rightunit(first(sectors(V)))
-        return spacetype(V)(sector => 1)
+        return spacetype(V)(_rightunit(V) => 1)
     end
 end
 
+function _rightunit(V::ElementarySpace)
+    !isempty(sectors(V)) || throw(ArgumentError("Cannot determine the right unit of an empty space"))
+    _allequal(rightunit, sectors(V)) ||
+        throw(ArgumentError(lazy"sectors of $V do not have the same right unit"))
+    return rightunit(first(sectors(V)))
+end
+
 """
     isunitspace(V::S) where {S <: ElementarySpace} -> Bool
 
 Return whether the elementary space `V` is a unit space, i.e. is isomorphic to the
 trivial one-dimensional space. For vector spaces of type `GradedSpace{I}` where `Sector` `I` has a
-semi-simple unit structure, this returns `true` if `V` is isomorphic to either the left, right or
-semi-simple unit space.
+semisimple unit structure, this returns `true` if `V` is isomorphic to either the left, right or
+semisimple unit space.
 """
 function isunitspace(V::ElementarySpace)
     I = sectortype(V)

test/setup.jl

@@ -7,6 +7,7 @@ export random_fusion
 export sectorlist, fast_sectorlist
 # export dim_isapprox
 export default_spacelist, factorization_spacelist, ad_spacelist
+export VIBM, VIBMRepA4
 export test_ad_rrule
 export _isunitary, _isone
 

test/symmetries/spaces.jl

@@ -233,12 +233,13 @@ end
     @test_throws ArgumentError("Sector $(randunit) appears multiple times") GradedSpace(randunit => 1, randunit => 3)
 
     @test isunitspace(W)
-    @test @constinferred(unitspace(V)) == W == unitspace(typeof(V))
     if UnitStyle(I) isa SimpleUnit
+        @test @constinferred(unitspace(V)) == W == unitspace(typeof(V))
         @test @constinferred(leftunitspace(V)) == W == @constinferred(rightunitspace(V))
     else
         @test_throws ArgumentError leftunitspace(V)
         @test_throws ArgumentError rightunitspace(V)
+        @test_throws ArgumentError unitspace(V)
     end
     @test eval_show(W) == W
     @test isa(V, VectorSpace)
@@ -261,8 +262,10 @@ end
     @test @constinferred(⊕(V, zerospace(V))) == V
     @test @constinferred(⊕(V, V)) == Vect[I](c => 2dim(V, c) for c in sectors(V))
     @test @constinferred(⊕(V, V, V, V)) == Vect[I](c => 4dim(V, c) for c in sectors(V))
-    @test @constinferred(⊕(V, unitspace(V))) == Vect[I](c => isunit(c) + dim(V, c) for c in sectors(V))
-    @test @constinferred(fuse(V, unitspace(V))) == V
+    if UnitStyle(I) isa SimpleUnit
+        @test @constinferred(⊕(V, unitspace(V))) == Vect[I](c => isunit(c) + dim(V, c) for c in sectors(V))
+        @test @constinferred(fuse(V, unitspace(V))) == V
+    end
     d = Dict{I, Int}()
     for a in sectors(V), b in sectors(V)
         for c in a ⊗ b
@@ -470,15 +473,44 @@ end
             @test @constinferred(insertrightunit(one(V1) ← V1, 0)) == (unitspace(V1) ← V1)
             @test_throws BoundsError insertleftunit(one(V1) ← V1, 0)
         else
-            @test_throws ArgumentError insertrightunit(one(V1) ← V1, 0)
-            @test_throws ArgumentError insertleftunit(one(V1) ← V1, 0)
+            errmsg = "cannot insert a sensible unit space in the empty product space"
+            @test_throws ArgumentError(errmsg) insertrightunit(one(V1) ← V1, 0)
+            @test_throws ArgumentError(errmsg) insertleftunit(one(V1) ← V1, 0)
         end
         @test (V1 ⊗ V2 ← V1 ⊗ V2) == @constinferred TensorKit.compose(W, W')
         @test W == @constinferred permute(W, ((1, 2), (3, 4, 5)))
         @test permute(W, ((2, 5, 4), (1, 3))) == (V2 ⊗ V3 ⊗ V4 ← V1' ⊗ V5') # cyclic permutation
     end
 end
 
+@timedtestset "ProductSpace and HomSpace: GenericUnit coloring" for V in (VIBM, VIBMRepA4)
+    @test UnitStyle(sectortype(V[1])) isa GenericUnit
+    V1, V2, V3, V4, V5 = V
+
+    @test @constinferred(one(V1)) == ProductSpace{typeof(V1)}(())
+
+    @test rightunitspace(V1) == leftunitspace(V2)
+    P1 = @constinferred ProductSpace(V1, V2)
+    @test @constinferred(⊗(V1, V2)) == P1
+
+    @test rightunitspace(V2) != leftunitspace(V1)
+    @test_throws ArgumentError (⊗(V2, V1))
+
+    @test rightunitspace(V3) == leftunitspace(V4)
+    @test rightunitspace(V4) == leftunitspace(V5)
+    P2 = @constinferred ProductSpace(V3, V4, V5)
+
+    @test leftunitspace(P1[1]) == rightunitspace(dual(P2[length(P2)]))
+    @test HomSpace(P1, P2') isa HomSpace
+    @test_throws ArgumentError P1 ← P2
+
+    @test (V1 ← one(V1)) isa HomSpace
+    @test (one(V1) ← one(V1)) isa HomSpace
+
+    @test leftunitspace(V2) != rightunitspace(V2)
+    @test_throws ArgumentError V2 ← one(V2)
+end
+
 @timedtestset "show and friends" begin
     V = U1Space(i => 1 for i in 1:3)
     @test string(V) == "$(type_repr(typeof(V)))(1 => 1, 2 => 1, 3 => 1)"