Fermion Refactor

docs/src/lib/sectors.md

@@ -15,6 +15,7 @@ ZNIrrep
 U1Irrep
 SU2Irrep
 CU1Irrep
+FermionParity
 FibonacciAnyon
 FusionTree
 ```

docs/src/man/sectors.md

@@ -1244,10 +1244,17 @@ the corresponding tensors.
 
 ## Fermions
 
-TODO
-
-(Support for fermionic sectors and corresponding super vector spaces is on its way. This
-section will be completed when the implementation is finished.)
+TODO: Update the documentation for this section.
+
+Fermionic sectors are represented by the type [`FermionParity`](@ref), which effectively
+behaves like a ℤ₂ sector, but with two modifications. Firstly, the exchange of two sectors
+with odd fermion parity should yield a minus sign, which is taken care of by virtue of the
+R-symbol. This ensures that permuting tensors behave as expected. Secondly, diagrams with
+self-crossing lines (aka twists) give rise to a minus sign for odd fermion parity. This is
+in essence equivalent to having supertraces, which is what ensures that `@tensor` has a
+result that is invariant under permutation of its input tensors. This does however lead to
+unwanted minus signs for certain types of diagrams. To avoid this, the `@planar` macro does
+not include a supertrace, but requires a manual resolution of all crossings in the diagram.
 
 ## Anyons
 

src/TensorKit.jl

@@ -12,10 +12,9 @@ export Sector, AbstractIrrep, Irrep
 export FusionStyle, UniqueFusion, MultipleFusion, MultiplicityFreeFusion,
        SimpleFusion, GenericFusion
 export BraidingStyle, SymmetricBraiding, Bosonic, Fermionic, Anyonic
-export Z2Irrep, Z3Irrep, Z4Irrep, ZNIrrep, U1Irrep, SU2Irrep, CU1Irrep
-export Fermion, FermionParity, FermionNumber, FermionSpin
-export FibonacciAnyon
-export IsingAnyon
+export Trivial, Z2Irrep, Z3Irrep, Z4Irrep, ZNIrrep, U1Irrep, SU2Irrep, CU1Irrep
+export FermionParity, FermionNumber, FermionSpin
+export FibonacciAnyon, IsingAnyon
 
 export VectorSpace, Field, ElementarySpace # abstract vector spaces
 export InnerProductStyle, NoInnerProduct, HasInnerProduct, EuclideanProduct
@@ -38,7 +37,7 @@ export infimum, supremum, isisomorphic, ismonomorphic, isepimorphic
 export sectortype, sectors, hassector, Nsymbol, Fsymbol, Rsymbol, Bsymbol,
        frobeniusschur, twist
 export fusiontrees, braid, permute, transpose
-export Trivial, ZNSpace, SU2Irrep, U1Irrep, CU1Irrep # Fermion
+export ZNSpace, SU2Irrep, U1Irrep, CU1Irrep
 # other fusion tree manipulations, should not be exported:
 # export insertat, split, merge, repartition, artin_braid,
 #        bendleft, bendright, foldleft, foldright, cycleclockwise, cycleanticlockwise

src/sectors/fermions.jl

@@ -1,96 +1,82 @@
-struct Fermion{P,I<:Sector} <: Sector
-    sector::I
-    function Fermion{P,I}(sector::I) where {P,I<:Sector}
-        @assert BraidingStyle(I) isa Bosonic
-        return new{P,I}(sector)
-    end
-end
-Fermion{P}(sector::I) where {P,I<:Sector} = Fermion{P,I}(sector)
-Fermion{P,I}(sector) where {P,I<:Sector} = Fermion{P,I}(convert(I, sector))
-Base.convert(::Type{Fermion{P,I}}, a::Fermion{P,I}) where {P,I<:Sector} = a
-Base.convert(::Type{Fermion{P,I}}, a) where {P,I<:Sector} = Fermion{P,I}(convert(I, a))
+"""
+    FermionParity <: Sector
 
-fermionparity(f::Fermion{P}) where {P} = P(f.sector)
+Represents sectors with fermion parity. The fermion parity is a ℤ₂ quantum number that
+yields an additional sign when two odd fermions are exchanged.
 
-function Base.IteratorSize(::Type{SectorValues{Fermion{P,I}}}) where {P,I<:Sector}
-    return Base.IteratorSize(SectorValues{I})
+See also: `FermionNumber`, `FermionSpin`
+"""
+struct FermionParity <: Sector
+    isodd::Bool
 end
-Base.length(::SectorValues{Fermion{P,I}}) where {P,I<:Sector} = length(values(I))
-function Base.iterate(::SectorValues{Fermion{P,I}}) where {P,I<:Sector}
-    next = iterate(values(I))
-    @assert next !== nothing
-    value, state = next
-    return Fermion{P}(value), state
-end
-function Base.iterate(::SectorValues{Fermion{P,I}}, state) where {P,I<:Sector}
-    next = iterate(values(I), state)
-    if next === nothing
-        return nothing
-    else
-        value, state = next
-        return Fermion{P}(value), state
-    end
-end
-function Base.getindex(::SectorValues{Fermion{P,I}}, i) where {P,I<:Sector}
-    return Fermion{P}(values(I)[i])
+const fℤ₂ = FermionParity
+fermionparity(f::FermionParity) = f.isodd
+
+Base.convert(::Type{FermionParity}, a::FermionParity) = a
+Base.convert(::Type{FermionParity}, a) = FermionParity(a)
+
+Base.IteratorSize(::Type{SectorValues{FermionParity}}) = HasLength()
+Base.length(::SectorValues{FermionParity}) = 2
+function Base.iterate(::SectorValues{FermionParity}, i=0)
+    return i == 2 ? nothing : (FermionParity(i), i + 1)
 end
-function findindex(::SectorValues{Fermion{P,I}}, f::Fermion{P,I}) where {P,I<:Sector}
-    return findindex(values(I), f.sector)
+function Base.getindex(::SectorValues{FermionParity}, i)
+    return 1 <= i <= 2 ? FermionParity(i - 1) : throw(BoundsError(values(FermionParity), i))
 end
+findindex(::SectorValues{FermionParity}, f::FermionParity) = f.isodd ? 2 : 1
 
-Base.one(::Type{Fermion{P,I}}) where {P,I<:Sector} = Fermion{P}(one(I))
-Base.conj(f::Fermion{P}) where {P} = Fermion{P}(conj(f.sector))
+Base.one(::Type{FermionParity}) = FermionParity(false)
+Base.conj(f::FermionParity) = f
+dim(f::FermionParity) = 1
 
-dim(f::Fermion) = dim(f.sector)
+FusionStyle(::Type{FermionParity}) = UniqueFusion()
+BraidingStyle(::Type{FermionParity}) = Fermionic()
+Base.isreal(::Type{FermionParity}) = true
 
-FusionStyle(::Type{<:Fermion{<:Any,I}}) where {I<:Sector} = FusionStyle(I)
-BraidingStyle(::Type{<:Fermion}) = Fermionic()
-Base.isreal(::Type{Fermion{<:Any,I}}) where {I<:Sector} = isreal(I)
+⊗(a::FermionParity, b::FermionParity) = (FermionParity(a.isodd ⊻ b.isodd),)
 
-⊗(a::F, b::F) where {F<:Fermion} = SectorSet{F}(a.sector ⊗ b.sector)
-
-Nsymbol(a::F, b::F, c::F) where {F<:Fermion} = Nsymbol(a.sector, b.sector, c.sector)
+function Nsymbol(a::FermionParity, b::FermionParity, c::FermionParity)
+    return (a.isodd ⊻ b.isodd) == c.isodd
+end
+function Fsymbol(a::I, b::I, c::I, d::I, e::I, f::I) where {I<:FermionParity}
+    return Int(Nsymbol(a, b, e) * Nsymbol(e, c, d) * Nsymbol(b, c, f) * Nsymbol(a, f, d))
+end
 
-function Fsymbol(a::F, b::F, c::F, d::F, e::F, f::F) where {F<:Fermion}
-    return Fsymbol(a.sector, b.sector, c.sector, d.sector, e.sector, f.sector)
+function Rsymbol(a::F, b::F, c::F) where {F<:FermionParity}
+    return a.isodd && b.isodd ? -Int(Nsymbol(a, b, c)) : Int(Nsymbol(a, b, c))
 end
+twist(a::FermionParity) = a.isodd ? -1 : +1
 
-function Rsymbol(a::F, b::F, c::F) where {F<:Fermion}
-    if fermionparity(a) && fermionparity(b)
-        return -Rsymbol(a.sector, b.sector, c.sector)
+function Base.show(io::IO, a::FermionParity)
+    if get(io, :typeinfo, nothing) === FermionParity
+        print(io, Int(a.isodd))
     else
-        return +Rsymbol(a.sector, b.sector, c.sector)
+        print(io, "FermionParity(", Int(a.isodd), ")")
     end
 end
+type_repr(::Type{FermionParity}) = "FermionParity"
 
-twist(a::Fermion) = ifelse(fermionparity(a), -1, +1) * twist(a.sector)
+Base.hash(f::FermionParity, h::UInt) = hash(f.isodd, h)
+Base.isless(a::FermionParity, b::FermionParity) = isless(a.isodd, b.isodd)
 
-type_repr(::Type{Fermion{P,I}}) where {P,I<:Sector} = "Fermion{$P, " * type_repr(I) * "}"
+# Common fermionic combinations
+# -----------------------------
 
-function Base.show(io::IO, a::Fermion{P,I}) where {P,I<:Sector}
-    if get(io, :typeinfo, nothing) !== Fermion{P,I}
-        print(io, type_repr(typeof(a)), "(")
-    end
-    print(IOContext(io, :typeinfo => I), a.sector)
-    if get(io, :typeinfo, nothing) !== Fermion{P,I}
-        print(io, ")")
-    end
-end
-
-Base.hash(f::Fermion, h::UInt) = hash(f.sector, h)
-Base.isless(a::F, b::F) where {F<:Fermion} = isless(a.sector, b.sector)
+const FermionNumber = U1Irrep ⊠ FermionParity
+const fU₁ = FermionNumber
+type_repr(::Type{FermionNumber}) = "FermionNumber"
 
-_fermionparity(a::Z2Irrep) = isodd(a.n)
-_fermionnumber(a::U1Irrep) = isodd(convert(Int, a.charge))
-_fermionspin(a::SU2Irrep) = isodd(twice(a.j))
+# convenience default converter -> allows Vect[FermionNumber](1 => 1)
+function Base.convert(::Type{FermionNumber}, a::Int)
+    return U1Irrep(a) ⊠ FermionParity(isodd(a))
+end
 
-const FermionParity = Fermion{_fermionparity,Z2Irrep}
-const FermionNumber = Fermion{_fermionnumber,U1Irrep}
-const FermionSpin = Fermion{_fermionspin,SU2Irrep}
-const fℤ₂ = FermionParity
-const fU₁ = FermionNumber
+const FermionSpin = SU2Irrep ⊠ FermionParity
 const fSU₂ = FermionSpin
-
-type_repr(::Type{FermionParity}) = "FermionParity"
-type_repr(::Type{FermionNumber}) = "FermionNumber"
 type_repr(::Type{FermionSpin}) = "FermionSpin"
+
+# convenience default converter -> allows Vect[FermionSpin](1 => 1)
+function Base.convert(::Type{FermionSpin}, a::Real)
+    s = SU2Irrep(a)
+    return s ⊠ FermionParity(isodd(twice(s.j)))
+end

src/sectors/sectors.jl

@@ -452,6 +452,6 @@ end
 # possible sectors
 include("groups.jl")
 include("irreps.jl") # irreps of symmetry groups, with bosonic braiding
+include("product.jl") # direct product of different sectors
 include("fermions.jl") # irreps with defined fermionparity and fermionic braiding
 include("anyons.jl") # non-group sectors
-include("product.jl") # direct product of different sectors

test/runtests.jl

@@ -18,7 +18,6 @@ const TK = TensorKit
 Random.seed!(1234)
 
 smallset(::Type{I}) where {I<:Sector} = take(values(I), 5)
-smallset(::Type{FermionNumber}) = FermionNumber.((0, +1, -1, +2, -2))
 function smallset(::Type{ProductSector{Tuple{I1,I2}}}) where {I1,I2}
     iter = product(smallset(I1), smallset(I2))
     s = collect(i ⊠ j for (i, j) in iter if dim(i) * dim(j) <= 6)
@@ -51,10 +50,11 @@ function hasfusiontensor(I::Type{<:Sector})
 end
 
 sectorlist = (Z2Irrep, Z3Irrep, Z4Irrep, U1Irrep, CU1Irrep, SU2Irrep, NewSU2Irrep, # SU3Irrep,
-              FibonacciAnyon, IsingAnyon, FermionParity, FermionNumber, FermionSpin,
-              FermionParity ⊠ FermionParity, Z3Irrep ⊠ Z4Irrep, FermionNumber ⊠ SU2Irrep,
-              FermionSpin ⊠ SU2Irrep, NewSU2Irrep ⊠ NewSU2Irrep, NewSU2Irrep ⊠ SU2Irrep,
-              FermionSpin ⊠ NewSU2Irrep, Z2Irrep ⊠ FibonacciAnyon ⊠ FibonacciAnyon)
+              FibonacciAnyon, IsingAnyon, FermionParity, FermionParity ⊠ FermionParity,
+              Z3Irrep ⊠ Z4Irrep, FermionParity ⊠ U1Irrep ⊠ SU2Irrep,
+              FermionParity ⊠ SU2Irrep ⊠ SU2Irrep, NewSU2Irrep ⊠ NewSU2Irrep,
+              NewSU2Irrep ⊠ SU2Irrep, FermionParity ⊠ SU2Irrep ⊠ NewSU2Irrep,
+              Z2Irrep ⊠ FibonacciAnyon ⊠ FibonacciAnyon)
 
 Ti = time()
 include("sectors.jl")

test/tensors.jl

@@ -44,10 +44,10 @@ VfSU₂ = (ℂ[FermionSpin](0 => 3, 1 // 2 => 1),
          ℂ[FermionSpin](0 => 2, 1 // 2 => 2),
          ℂ[FermionSpin](0 => 1, 1 // 2 => 1, 3 // 2 => 1)')
 # VSU₃ = (ℂ[SU3Irrep]((0, 0, 0) => 3, (1, 0, 0) => 1),
-#         ℂ[SU3Irrep]((0, 0, 0) => 3, (2, 0, 0) => 1)',
-#         ℂ[SU3Irrep]((1, 1, 0) => 1, (2, 1, 0) => 1),
-#         ℂ[SU3Irrep]((1, 0, 0) => 1, (2, 0, 0) => 1),
-#         ℂ[SU3Irrep]((0, 0, 0) => 1, (1, 0, 0) => 1, (1, 1, 0) => 1)')
+#     ℂ[SU3Irrep]((0, 0, 0) => 3, (2, 0, 0) => 1)',
+#     ℂ[SU3Irrep]((1, 1, 0) => 1, (2, 1, 0) => 1),
+#     ℂ[SU3Irrep]((1, 0, 0) => 1, (2, 0, 0) => 1),
+#     ℂ[SU3Irrep]((0, 0, 0) => 1, (1, 0, 0) => 1, (1, 1, 0) => 1)')
 
 for V in (Vtr, Vℤ₂, Vfℤ₂, Vℤ₃, VU₁, VfU₁, VCU₁, VSU₂, VfSU₂)#, VSU₃)
     V1, V2, V3, V4, V5 = V