Set last dimension to 0 for fusiontensor of forbidden fusions

.gitignore

@@ -2,5 +2,6 @@
 *.jl.cov
 *.jl.mem
 Manifest.toml
+.vscode/
 .DS_Store
 build

src/fermions.jl

@@ -55,7 +55,8 @@ twist(a::FermionParity) = a.isodd ? -1 : +1
 
 function fusiontensor(a::I, b::I, c::I) where {I <: FermionParity}
     @warn "FermionParity Arrays do not preserve categorical properties." maxlog = 1
-    return fill(Int(Nsymbol(a, b, c)), (1, 1, 1, 1))
+    Nabc = Nsymbol(a, b, c)
+    return fill(Int(Nabc), (1, 1, 1, Nabc))
 end
 
 function Base.show(io::IO, a::FermionParity)

src/irreps/cu1irrep.jl

@@ -218,8 +218,9 @@ function Rsymbol(a::CU1Irrep, b::CU1Irrep, c::CU1Irrep)
 end
 
 function fusiontensor(a::CU1Irrep, b::CU1Irrep, c::CU1Irrep)
-    C = fill(zero(sectorscalartype(CU1Irrep)), dim(a), dim(b), dim(c), 1)
-    !Nsymbol(a, b, c) && return C
+    Nabc = Nsymbol(a, b, c)
+    C = fill(zero(sectorscalartype(CU1Irrep)), dim(a), dim(b), dim(c), Nabc)
+    Nabc || return C
     if c.j == 0
         if a.j == b.j == 0
             C[1, 1, 1, 1] = 1.0

src/irreps/dnirrep.jl

@@ -198,8 +198,9 @@ end
 
 function fusiontensor(a::I, b::I, c::I) where {N, I <: DNIrrep{N}}
     T = sectorscalartype(I)
-    C = zeros(T, dim(a), dim(b), dim(c), 1)
-    Nsymbol(a, b, c) || return C
+    Nabc = Nsymbol(a, b, c)
+    C = zeros(T, dim(a), dim(b), dim(c), Nabc)
+    Nabc || return C
 
     if c.j == 0
         if a.j == b.j == 0 || (2 * a.j == 2 * b.j == N)

src/irreps/irreps.jl

@@ -66,7 +66,8 @@ Rsymbol(a::I, b::I, c::I) where {I <: AbelianIrrep} = braidingscalartype(I)(Nsym
 
 function fusiontensor(a::I, b::I, c::I) where {I <: AbelianIrrep}
     T = sectorscalartype(I)
-    return fill(T(Nsymbol(a, b, c)), (1, 1, 1, 1))
+    Nabc = Nsymbol(a, b, c)
+    return fill(T(Nabc), (1, 1, 1, Nabc))
 end
 
 include("znirrep.jl")

src/irreps/su2irrep.jl

@@ -72,7 +72,9 @@ end
 
 function fusiontensor(a::SU2Irrep, b::SU2Irrep, c::SU2Irrep)
     T = sectorscalartype(SU2Irrep)
-    C = Array{T}(undef, dim(a), dim(b), dim(c), 1)
+    Nabc = Nsymbol(a, b, c)
+    C = Array{T}(undef, dim(a), dim(b), dim(c), Nabc)
+    Nabc || return C
     ja, jb, jc = a.j, b.j, c.j
 
     for kc in 1:dim(c), kb in 1:dim(b), ka in 1:dim(a)

test/newsectors.jl

@@ -65,7 +65,9 @@ end
 dim(s::NewSU2Irrep) = twice(s.j) + 1
 
 function fusiontensor(a::NewSU2Irrep, b::NewSU2Irrep, c::NewSU2Irrep)
-    C = Array{Float64}(undef, dim(a), dim(b), dim(c), 1)
+    Nabc = Nsymbol(a, b, c)
+    C = Array{Float64}(undef, dim(a), dim(b), dim(c), Nabc)
+    Nabc == 0 && return C
     ja, jb, jc = a.j, b.j, c.j
 
     for kc in 1:dim(c), kb in 1:dim(b), ka in 1:dim(a)