Sectors refactor

Project.toml

@@ -1,7 +1,7 @@
 name = "GradedArrays"
 uuid = "bc96ca6e-b7c8-4bb6-888e-c93f838762c2"
 authors = ["ITensor developers <[email protected]> and contributors"]
-version = "0.4.26"
+version = "0.5.0"
 
 [deps]
 ArrayLayouts = "4c555306-a7a7-4459-81d9-ec55ddd5c99a"
@@ -15,9 +15,16 @@ MatrixAlgebraKit = "6c742aac-3347-4629-af66-fc926824e5e4"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 SplitApplyCombine = "03a91e81-4c3e-53e1-a0a4-9c0c8f19dd66"
 TensorAlgebra = "68bd88dc-f39d-4e12-b2ca-f046b68fcc6a"
+TensorKitSectors = "13a9c161-d5da-41f0-bcbd-e1a08ae0647f"
 TensorProducts = "decf83d6-1968-43f4-96dc-fdb3fe15fc6d"
 TypeParameterAccessors = "7e5a90cf-f82e-492e-a09b-e3e26432c138"
 
+[weakdeps]
+SUNRepresentations = "1a50b95c-7aac-476d-a9ce-2bfc675fc617"
+
+[extensions]
+GradedArraysSUNRepresentationsExt = "SUNRepresentations"
+
 [compat]
 ArrayLayouts = "1"
 BlockArrays = "1.6"
@@ -28,8 +35,10 @@ HalfIntegers = "1.6"
 LinearAlgebra = "1.10"
 MatrixAlgebraKit = "0.2, 0.3, 0.4, 0.5"
 Random = "1.10"
+SUNRepresentations = "0.3"
 SplitApplyCombine = "1.2.3"
 TensorAlgebra = "0.3.2, 0.4"
+TensorKitSectors = "0.1, 0.2"
 TensorProducts = "0.1.3"
 TypeParameterAccessors = "0.4"
 julia = "1.10"

docs/Project.toml

@@ -5,5 +5,5 @@ Literate = "98b081ad-f1c9-55d3-8b20-4c87d4299306"
 
 [compat]
 Documenter = "1"
-GradedArrays = "0.4"
+GradedArrays = "0.5"
 Literate = "2"

examples/Project.toml

@@ -2,4 +2,4 @@
 GradedArrays = "bc96ca6e-b7c8-4bb6-888e-c93f838762c2"
 
 [compat]
-GradedArrays = "0.4"
+GradedArrays = "0.5"

ext/GradedArraysSUNRepresentationsExt.jl

@@ -0,0 +1,12 @@
+module GradedArraysSUNRepresentationsExt
+
+using SUNRepresentations: SUNIrrep
+using GradedArrays: GradedArrays, SectorRange
+
+function GradedArrays.SectorRange{SUNIrrep{N}}(λ::NTuple{M,Int}) where {N,M}
+  M + 1 == N || throw(ArgumentError("Length of λ must be N-1 for SU(N) irreps"))
+  return SectorRange(SUNIrrep((λ..., 0)))
+end
+GradedArrays.sector_label(c::SUNIrrep) = Base.front(c.I)
+
+end

src/GradedArrays.jl

@@ -1,20 +1,11 @@
 module GradedArrays
 
 include("gradedunitrange_interface.jl")
-include("symmetry_style.jl")
 
+include("abstractsector.jl")
 include("sectorunitrange.jl")
 include("gradedunitrange.jl")
 
-include("abstractsector.jl")
-include("sector_definitions/fib.jl")
-include("sector_definitions/ising.jl")
-include("sector_definitions/o2.jl")
-include("sector_definitions/trivial.jl")
-include("sector_definitions/su.jl")
-include("sector_definitions/su2k.jl")
-include("sector_definitions/u1.jl")
-include("sector_definitions/zn.jl")
 include("namedtuple_operations.jl")
 include("sector_product.jl")
 
@@ -23,10 +14,8 @@ include("gradedarray.jl")
 include("tensoralgebra.jl")
 include("factorizations.jl")
 
-export SU2,
-  U1,
-  Z,
-  dag,
+export TrivialSector, Z, Z2, U1, O2, SU2, Fib, Ising
+export dag,
   dual,
   flip,
   gradedrange,

src/abstractsector.jl

@@ -1,82 +1,211 @@
-# This file defines the abstract type AbstractSector
-# all fusion categories (Z{2}, SU2, Ising...) are subtypes of AbstractSector
-
+# This file defines the interface for type Sector
+# all fusion categories (Z{2}, SU2, Ising...) are subtypes of Sector
 using TensorProducts: TensorProducts, ⊗
+import TensorKitSectors as TKS
 
-abstract type AbstractSector end
+"""
+    SectorRange(sector::TKS.Sector)
 
-# ===================================  Base interface  =====================================
-function Base.isless(c1::C, c2::C) where {C<:AbstractSector}
-  return isless(sector_label(c1), sector_label(c2))
+Unit range with elements of type `Int` that additionally stores a sector to denote the grading.
+Equivalent to `Base.OneTo(length(sector))`.
+"""
+struct SectorRange{I<:TKS.Sector} <: AbstractUnitRange{Int}
+  label::I
 end
 
-Base.length(s::AbstractSector) = quantum_dimension(s)
+label(r::SectorRange) = r.label
+sector_type(I::Type{<:SectorRange}) = I
+
+# ===================================  Base interface  =====================================
+
+Base.length(r::SectorRange) = quantum_dimension(r)
+
+Base.isless(r1::SectorRange, r2::SectorRange) = isless(label(r1), label(r2))
+Base.isless(r1::SectorRange, r2::TKS.Sector) = isless(label(r1), r2)
+Base.isless(r1::TKS.Sector, r2::SectorRange) = isless(r1, label(r2))
+
+Base.isequal(r1::SectorRange, r2::SectorRange) = isequal(label(r1), label(r2))
+Base.:(==)(r1::SectorRange, r2::SectorRange) = label(r1) == label(r2)
+Base.:(==)(r1::SectorRange, r2::TKS.Sector) = label(r1) == r2
+Base.:(==)(r1::TKS.Sector, r2::SectorRange) = r1 == label(r2)
+
+Base.hash(r::SectorRange, h::UInt) = hash(label(r), h)
+
+Base.OneTo(r::SectorRange) = Base.OneTo(length(r))
+Base.first(r::SectorRange) = first(Base.OneTo(r))
+Base.last(r::SectorRange) = last(Base.OneTo(r))
+
+function Base.show(io::IO, r::SectorRange{I}) where {I}
+  show(io, typeof(r))
+  print(io, '(')
+  l = sector_label(r)
+  isnothing(l) || show(io, l)
+  print(io, ')')
+  return nothing
+end
 
 # =================================  Sectors interface  ====================================
+
 trivial(x) = trivial(typeof(x))
 function trivial(axis_type::Type{<:AbstractUnitRange})
   return gradedrange([trivial(sector_type(axis_type)) => 1])  # always returns nondual
 end
 function trivial(type::Type)
   return error("`trivial` not defined for type $(type).")
 end
-
-istrivial(c::AbstractSector) = (c == trivial(c))
-
-function sector_label(c::AbstractSector)
-  return error("method `sector_label` not defined for type $(typeof(c))")
+trivial(::Type{SectorRange{I}}) where {I} = SectorRange{I}(one(I))
+trivial(::Type{I}) where {I<:TKS.Sector} = one(I)
+
+istrivial(r::SectorRange) = isone(label(r))
+istrivial(r) = (r == trivial(r))
+
+sector_label(r::SectorRange) = sector_label(label(r))
+function sector_label(c::TKS.Sector)
+  return map(fieldnames(typeof(c))) do f
+    return getfield(c, f)
+  end
+  return c
 end
 
 quantum_dimension(g::AbstractUnitRange) = length(g)
-quantum_dimension(s::AbstractSector) = quantum_dimension(SymmetryStyle(s), s)
-
-function quantum_dimension(::NotAbelianStyle, c::AbstractSector)
-  return error("method `quantum_dimension` not defined for type $(typeof(c))")
-end
+quantum_dimension(r::SectorRange) = quantum_dimension(label(r))
+quantum_dimension(s::TKS.Sector) = TKS.dim(s)
 
-quantum_dimension(::AbelianStyle, ::AbstractSector) = 1
+to_sector(x::TKS.Sector) = SectorRange(x)
 
 # convert to range
-to_gradedrange(c::AbstractSector) = gradedrange([c => 1])
-to_gradedrange(sr::SectorUnitRange) = mortar_axis([sr])
-to_gradedrange(g::AbstractGradedUnitRange) = g
+to_gradedrange(c::SectorRange) = gradedrange([c => 1])
+to_gradedrange(c::TKS.Sector) = to_gradedrange(SectorRange(c))
 
-function nsymbol(s1::AbstractSector, s2::AbstractSector, s3::AbstractSector)
-  full_space = to_gradedrange(s1 ⊗ s2)
-  i = findfirst(==(s3), sectors(full_space))
-  isnothing(i) && return 0
-  return sector_multiplicities(full_space)[i]
+function nsymbol(s1::SectorRange, s2::SectorRange, s3::SectorRange)
+  return TKS.Nsymbol(label(s1), label(s2), label(s3))
 end
 
+dual(c::TKS.Sector) = TKS.dual(c)
+dual(r1::SectorRange) = typeof(r1)(dual(label(r1)))
+
 # ===============================  Fusion rule interface  ==================================
-function fusion_rule(c1::AbstractSector, c2::AbstractSector)
-  return fusion_rule(combine_styles(SymmetryStyle(c1), SymmetryStyle(c2)), c1, c2)
-end
 
-function fusion_rule(::NotAbelianStyle, c1::C, c2::C) where {C<:AbstractSector}
-  sector_degen_pairs = label_fusion_rule(C, sector_label(c1), sector_label(c2))
-  return gradedrange(sector_degen_pairs)
-end
+TKS.FusionStyle(::Type{SectorRange{I}}) where {I} = TKS.FusionStyle(I)
+TKS.BraidingStyle(::Type{SectorRange{I}}) where {I} = TKS.BraidingStyle(I)
 
-# abelian case: return Sector
-function fusion_rule(::AbelianStyle, c1::C, c2::C) where {C<:AbstractSector}
-  return only(sectors(fusion_rule(NotAbelianStyle(), c1, c2)))
-end
+abstract type SymmetryStyle end
+
+struct AbelianStyle <: SymmetryStyle end
+struct NotAbelianStyle <: SymmetryStyle end
 
-function label_fusion_rule(sector_type::Type{<:AbstractSector}, l1, l2)
-  return [abelian_label_fusion_rule(sector_type, l1, l2) => 1]
+SymmetryStyle(x) = SymmetryStyle(typeof(x))
+
+# default SymmetryStyle to AbelianStyle
+# allows for abelian-like slicing style for GradedUnitRange: assume length(::label) = 1
+# and preserve labels in any slicing operation
+SymmetryStyle(T::Type) = AbelianStyle()
+function SymmetryStyle(::Type{T}) where {T<:SectorRange}
+  if TKS.FusionStyle(T) === TKS.UniqueFusion() && TKS.BraidingStyle(T) === TKS.Bosonic()
+    return AbelianStyle()
+  else
+    return NotAbelianStyle()
+  end
+end
+SymmetryStyle(G::Type{<:AbstractUnitRange}) = SymmetryStyle(sector_type(G))
+
+combine_styles(::AbelianStyle, ::AbelianStyle) = AbelianStyle()
+combine_styles(::SymmetryStyle, ::SymmetryStyle) = NotAbelianStyle()
+
+function fusion_rule(r1::SectorRange, r2::SectorRange)
+  a = label(r1)
+  b = label(r2)
+  fstyle = TKS.FusionStyle(typeof(r1)) & TKS.FusionStyle(typeof(r2))
+  fstyle === TKS.UniqueFusion() && return SectorRange(only(TKS.otimes(a, b)))
+  return gradedrange(
+    vec([SectorRange(c) => TKS.Nsymbol(a, b, c) for c in TKS.otimes(a, b)])
+  )
 end
 
 # =============================  TensorProducts interface  =====--==========================
-TensorProducts.tensor_product(s::AbstractSector) = s
 
-function TensorProducts.tensor_product(c1::AbstractSector, c2::AbstractSector)
-  return fusion_rule(c1, c2)
+TensorProducts.tensor_product(s::SectorRange) = s
+TensorProducts.tensor_product(c1::SectorRange, c2::SectorRange) = fusion_rule(c1, c2)
+function TensorProducts.tensor_product(c1::TKS.Sector, c2::TKS.Sector)
+  return tensor_product(to_sector(c1), to_sector(c2))
+end
+function TensorProducts.tensor_product(c1::SectorRange, c2::TKS.Sector)
+  return tensor_product(c1, to_sector(c2))
+end
+function TensorProducts.tensor_product(c1::TKS.Sector, c2::SectorRange)
+  return tensor_product(to_sector(c1), c2)
+end
+
+# =====================================  Sectors ===========================================
+
+const TrivialSector = SectorRange{TKS.Trivial}
+TrivialSector() = TrivialSector(TKS.Trivial())
+sector_label(::TKS.Trivial) = nothing
+function fusion_rule(::TrivialSector, r::SectorRange)
+  return TKS.FusionStyle(label(r)) === TKS.UniqueFusion() ? r : to_gradedrange(r)
+end
+function fusion_rule(r::SectorRange, ::TrivialSector)
+  return TKS.FusionStyle(label(r)) === TKS.UniqueFusion() ? r : to_gradedrange(r)
+end
+fusion_rule(r::TrivialSector, ::TrivialSector) = r
+
+Base.:(==)(::TrivialSector, ::TrivialSector) = true
+Base.:(==)(::TrivialSector, r::SectorRange) = isone(label(r))
+Base.:(==)(r::SectorRange, ::TrivialSector) = isone(label(r))
+Base.isless(::TrivialSector, ::TrivialSector) = false
+Base.isless(::TrivialSector, r::SectorRange) = trivial(r) < r
+Base.isless(r::SectorRange, ::TrivialSector) = r < trivial(r)
+
+# use promotion to handle trivial sectors in tensor products
+Base.promote_rule(::Type{TrivialSector}, ::Type{T}) where {T<:SectorRange} = T
+Base.convert(::Type{T}, ::TrivialSector) where {T<:SectorRange} = trivial(T)
+
+const Z{N} = SectorRange{TKS.ZNIrrep{N}}
+sector_label(c::TKS.ZNIrrep) = c.n
+modulus(::Z{N}) where {N} = N
+const Z2 = Z{2}
+
+const U1 = SectorRange{TKS.U1Irrep}
+sector_label(c::TKS.U1Irrep) = c.charge
+Base.isless(r1::U1, r2::U1) = isless(sector_label(r1), sector_label(r2))
+
+const O2 = SectorRange{TKS.CU1Irrep}
+function O2(l::Real)
+  j = max(l, zero(l))
+  s = if l == 0
+    0
+  elseif l == -1
+    1
+  else
+    2
+  end
+  return O2(TKS.CU1Irrep(j, s))
+end
+function sector_label(c::TKS.CU1Irrep)
+  return if c.s == 0
+    c.j
+  elseif c.s == 1
+    oftype(c.j, -1)
+  else
+    c.j
+  end
 end
 
-# ================================  GradedUnitRanges interface  ==================================
-sector_type(S::Type{<:AbstractSector}) = S
+const SU2 = SectorRange{TKS.SU2Irrep}
+sector_label(c::TKS.SU2Irrep) = c.j
+
+const Fib = SectorRange{TKS.FibonacciAnyon}
+function Fib(s::AbstractString)
+  s == "1" && return Fib(:I)
+  s == "τ" && return Fib(:τ)
+  throw(ArgumentError("Unrecognized input `$s`"))
+end
+sector_label(c::TKS.FibonacciAnyon) = c.isone ? "1" : "τ"
 
-function findfirstblock(g::AbstractGradedUnitRange, s::AbstractSector)
-  return findfirstblock_sector(g::AbstractGradedUnitRange, s)
+const Ising = SectorRange{TKS.IsingAnyon}
+function Ising(s::AbstractString)
+  s in ("1", "σ", "ψ") || throw(ArgumentError("Unrecognized input `$s`"))
+  sym = s == "1" ? :I : Symbol(s)
+  return Ising(sym)
 end
+sector_label(c::TKS.IsingAnyon) = Symbol(c.s) == :I ? "1" : String(c.s)

src/factorizations.jl

@@ -58,9 +58,9 @@ function fluxify(A, Aaxes, charge; side::Symbol=:domain)
   istrivial(charge) && return TensorAlgebra.unmatricize(A, (Aaxes[1],), (Aaxes[2],))
 
   if side === :domain
-    A′ = TensorAlgebra.unmatricize(A, (Aaxes[1],), (Aaxes[2], to_gradedrange(charge)))
+    A′ = TensorAlgebra.unmatricize(A, (Aaxes[1],), (Aaxes[2], to_gradedrange(dual(charge))))
   else
-    A′ = TensorAlgebra.unmatricize(A, (to_gradedrange(charge), Aaxes[1]), (Aaxes[2],))
+    A′ = TensorAlgebra.unmatricize(A, (to_gradedrange(dual(charge)), Aaxes[1]), (Aaxes[2],))
   end
 
   A″ = similar(A′, Aaxes)

src/fusion.jl

@@ -28,22 +28,22 @@ end
 
 # allow to fuse a Sector with a GradedUnitRange
 function TensorProducts.tensor_product(
-  s::Union{AbstractSector,SectorUnitRange}, g::AbstractGradedUnitRange
+  s::Union{SectorRange,SectorUnitRange}, g::AbstractGradedUnitRange
 )
   return to_gradedrange(s) ⊗ g
 end
 
 function TensorProducts.tensor_product(
-  g::AbstractGradedUnitRange, s::Union{AbstractSector,SectorUnitRange}
+  g::AbstractGradedUnitRange, s::Union{SectorRange,SectorUnitRange}
 )
   return g ⊗ to_gradedrange(s)
 end
 
-function TensorProducts.tensor_product(sr::SectorUnitRange, s::AbstractSector)
+function TensorProducts.tensor_product(sr::SectorUnitRange, s::SectorRange)
   return sr ⊗ sectorrange(s, 1)
 end
 
-function TensorProducts.tensor_product(s::AbstractSector, sr::SectorUnitRange)
+function TensorProducts.tensor_product(s::SectorRange, sr::SectorUnitRange)
   return sectorrange(s, 1) ⊗ sr
 end