IsingAnyons and a bugfix

src/TensorKit.jl

@@ -14,6 +14,7 @@ export Abelian, NonAbelian, SimpleNonAbelian, DegenerateNonAbelian,
 export Z2Irrep, Z3Irrep, Z4Irrep, ZNIrrep, U1Irrep, SU2Irrep, CU1Irrep
         #FermionParity, FermionNumber, FermionSpin # specific sectors
 export FibonacciAnyon
+export IsingAnyon
 
 export VectorSpace, Field, ElementarySpace, InnerProductSpace, EuclideanSpace # abstract vector spaces
 export ComplexSpace, CartesianSpace, GeneralSpace, GradedSpace, ZNSpace, Rep # concrete spaces

src/sectors/anyons.jl

@@ -1,22 +1,20 @@
 # FibonacciAnyons
 """
     struct FibonacciAnyon <: Sector
-    FibonacciAnyon(s::Union{Symbol, Integer})
+    FibonacciAnyon(s::Symbol)
 
 Represents the anyons (isomorphism classes of simple objects) of the Fibonacci fusion
 category. It can take two values, corresponding to the trivial sector
-`FibonacciAnyon(:I) == FibonacciAnyon(0)` and the non-trivial sector
-`FibonacciAnyon(:τ) = FibonacciAnyon(1)` with fusion rules ``τ ⊗ τ = 1 ⊕ τ``.
+`FibonacciAnyon(:I)` and the non-trivial sector `FibonacciAnyon(:τ)` with fusion rules
+``τ ⊗ τ = 1 ⊕ τ``.
 """
 struct FibonacciAnyon <: Sector
     isone::Bool
     function FibonacciAnyon(s::Symbol)
-        s == :I || s == :τ || throw(ArgumentError("Unknown FibonacciAnyon $s."))
+        s in (:I, :τ, :tau) || throw(ArgumentError("Unknown FibonacciAnyon $s."))
         new(s === :I)
     end
 end
-Fibonacci(i::Integer) = iszero(i) ? Fibonacci(:I) :
-                            (isone(i) ? Fibonacci(:τ) : error("unkown Fibonacci anyon"))
 
 Base.IteratorSize(::Type{SectorValues{FibonacciAnyon}}) = HasLength()
 Base.length(::SectorValues{FibonacciAnyon}) = 2
@@ -37,7 +35,7 @@ Base.convert(::Type{FibonacciAnyon}, s::Symbol) = FibonacciAnyon(s)
 Base.one(::Type{FibonacciAnyon}) = FibonacciAnyon(:I)
 Base.conj(s::FibonacciAnyon) = s
 
-const _goldenratio = (1 + sqrt(5)) / 2
+const _goldenratio = Float64(MathConstants.golden)
 dim(a::FibonacciAnyon) = isone(a) ? one(_goldenratio) : _goldenratio
 
 FusionStyle(::Type{FibonacciAnyon}) = SimpleNonAbelian()
@@ -95,11 +93,11 @@ function Fsymbol(a::FibonacciAnyon, b::FibonacciAnyon, c::FibonacciAnyon,
 end
 
 function Rsymbol(a::FibonacciAnyon, b::FibonacciAnyon, c::FibonacciAnyon)
-    Nsymbol(a, b, c) || return 0*exp((0π/1)*im)
+    Nsymbol(a, b, c) || return 0*cis(0π/1)
     if isone(a) || isone(b)
-        return exp((0π/1)*im)
+        return cis(0π/1)
     else
-        return isone(c) ? exp(+(4π/5)*im) : exp(-(3π/5)*im)
+        return isone(c) ? cis(4π/5) : cis(-3π/5)
     end
 end
 
@@ -112,3 +110,149 @@ end
 
 Base.hash(a::FibonacciAnyon, h::UInt) = hash(a.isone, h)
 Base.isless(a::FibonacciAnyon, b::FibonacciAnyon) = isless(!a.isone, !b.isone)
+
+# IsingAnyons
+"""
+    struct IsingAnyon <: Sector
+    IsingAnyon(s::Symbol)
+
+Represents the anyons (isomorphism classes of simple objects) of the Ising fusion category.
+It can take three values, corresponding to the trivial sector `IsingAnyon(:I)` and the
+non-trivial sectors `IsingAnyon(:σ)` and `IsingAnyon(:ψ)`, with fusion rules
+``ψ ⊗ ψ = 1``, ``σ ⊗ ψ = σ``, and ``σ ⊗ σ = 1 ⊕ ψ``.
+"""
+struct IsingAnyon <: Sector
+    s::Symbol
+    function IsingAnyon(s::Symbol)
+        s == :sigma && (s = :σ)
+        s == :psi && (s = :ψ)
+        if !(s in (:I, :σ, :ψ))
+            throw(ValueError("Unknown IsingAnyon $s."))
+        end
+        new(s)
+    end
+end
+
+const all_isinganyons = (IsingAnyon(:I), IsingAnyon(:σ), IsingAnyon(:ψ))
+
+Base.IteratorSize(::Type{SectorValues{IsingAnyon}}) = HasLength()
+Base.length(::SectorValues{IsingAnyon}) = length(all_isinganyons)
+Base.iterate(::SectorValues{IsingAnyon}, i = 1) = iterate(all_isinganyons, i)
+Base.getindex(S::SectorValues{IsingAnyon}, i) = getindex(all_isinganyons, i)
+
+function findindex(::SectorValues{IsingAnyon}, a::IsingAnyon)
+    a == all_isinganyons[1] && return 1
+    a == all_isinganyons[2] && return 2
+    return 3
+end
+
+Base.convert(::Type{IsingAnyon}, s::Symbol) = IsingAnyon(s)
+Base.one(::Type{IsingAnyon}) = IsingAnyon(:I)
+Base.conj(s::IsingAnyon) = s
+
+dim(a::IsingAnyon) = a.s == :σ ? sqrt(2) : 1.0
+
+FusionStyle(::Type{IsingAnyon}) = SimpleNonAbelian()
+BraidingStyle(::Type{IsingAnyon}) = Anyonic()
+Base.isreal(::Type{IsingAnyon}) = false
+
+⊗(a::IsingAnyon, b::IsingAnyon) = IsingIterator(a, b)
+
+struct IsingIterator
+    a::IsingAnyon
+    b::IsingAnyon
+end
+
+Base.IteratorSize(::Type{IsingIterator}) = Base.HasLength()
+Base.IteratorEltype(::Type{IsingIterator}) = Base.HasEltype()
+Base.eltype(::Type{IsingIterator}) = IsingAnyon
+
+function Base.length(iter::IsingIterator)
+    σ = IsingAnyon(:σ)
+    return (iter.a == σ && iter.b == σ) ? 2 : 1
+end
+
+function Base.iterate(iter::IsingIterator, state = 1)
+    I, σ, ψ = all_isinganyons
+    if state == 1 # first iteration
+        iter.a == I && return (iter.b, 2)
+        iter.b == I && return (iter.a, 2)
+        iter.a == σ && iter.b == ψ && return (σ, 2)
+        iter.a == ψ && iter.b == σ && return (σ, 2)
+        return (I, 2)
+    elseif state == 2
+        (iter.a == iter.b == σ) && return (ψ, 3)
+        return nothing
+    else
+        return nothing
+    end
+end
+
+function Nsymbol(a::IsingAnyon, b::IsingAnyon, c::IsingAnyon)
+    I, σ, ψ = all_isinganyons
+    return ((a == I && b == c)
+            || (b == I && a == c)
+            || (c == I && a == b)
+            || (a == σ && b == σ && c == ψ)
+            || (a == σ && b == ψ && c == σ)
+            || (a == ψ && b == σ && c == σ)
+           )
+end
+
+function Fsymbol(a::IsingAnyon, b::IsingAnyon, c::IsingAnyon,
+                 d::IsingAnyon, e::IsingAnyon, f::IsingAnyon)
+    Nsymbol(a, b, e) || return 0.0
+    Nsymbol(e, c, d) || return 0.0
+    Nsymbol(b, c, f) || return 0.0
+    Nsymbol(a, f, d) || return 0.0
+    I, σ, ψ = all_isinganyons
+    if a == b == c == d == σ
+        if e == f == ψ
+            return -1.0/sqrt(2.0)
+        else
+            return 1.0/sqrt(2.0)
+        end
+    end
+    if e == f == σ
+        if a == c == σ && b == d == ψ
+            return -1.0
+        elseif a == c == ψ && b == d == σ
+            return -1.0
+        end
+    end
+    return 1.0
+end
+
+function Rsymbol(a::IsingAnyon, b::IsingAnyon, c::IsingAnyon)
+    Nsymbol(a, b, c) || return complex(0.0)
+    I, σ, ψ = all_isinganyons
+    if c == I
+        if b == a == σ
+            return cis(-π/8)
+        elseif b == a == ψ
+            return complex(-1.0)
+        end
+    elseif c == σ && (a == σ && b == ψ || a == ψ && b == σ)
+        return -1.0im
+    elseif c == ψ && a == b == σ
+        return cis(3π/8)
+    end
+    return complex(1.0)
+end
+
+Base.show(io::IO, ::Type{IsingAnyon}) = print(io, "IsingAnyon")
+
+function Base.show(io::IO, a::IsingAnyon)
+    if get(io, :typeinfo, nothing) === IsingAnyon
+        return print(io, ":$(a.s)")
+    else
+        return print(io, "IsingAnyon(:$(a.s))")
+    end
+end
+
+Base.hash(s::IsingAnyon, h::UInt) = hash(s.s, h)
+
+function Base.isless(a::IsingAnyon, b::IsingAnyon)
+    vals = SectorValues{IsingAnyon}()
+    return isless(findindex(vals, a), findindex(vals, b))
+end

src/sectors/sectors.jl

@@ -68,7 +68,7 @@ Base.show(io::IO, ::Trivial) = print(io, "Trivial()")
 
 Base.IteratorSize(::Type{SectorValues{Trivial}}) = HasLength()
 Base.length(::SectorValues{Trivial}) = 1
-Base.iterate(::SectorValues{Trivial}, done = false) = return i ? nothing : (Trivial(), true)
+Base.iterate(::SectorValues{Trivial}, i = false) = return i ? nothing : (Trivial(), true)
 Base.getindex(::SectorValues{Trivial}, i::Int) = i == 1 ? Trivial() :
     throw(BoundsError(values(Trivial), i))
 findindex(::SectorValues{Trivial}, c::Trivial) = 1

test/runtests.jl

@@ -57,9 +57,9 @@ function hasfusiontensor(I::Type{<:Sector})
 end
 
 sectorlist = (Z2Irrep, Z3Irrep, Z4Irrep, U1Irrep, CU1Irrep, SU2Irrep, NewSU2Irrep,
-            FibonacciAnyon, Z3Irrep ⊠ Z4Irrep, U1Irrep ⊠ SU2Irrep, SU2Irrep ⊠ SU2Irrep,
-            NewSU2Irrep ⊠ NewSU2Irrep, NewSU2Irrep ⊠ SU2Irrep, SU2Irrep ⊠ NewSU2Irrep,
-            Z2Irrep ⊠ FibonacciAnyon ⊠ FibonacciAnyon)
+              FibonacciAnyon, IsingAnyon, Z3Irrep ⊠ Z4Irrep, U1Irrep ⊠ SU2Irrep, SU2Irrep ⊠
+              SU2Irrep, NewSU2Irrep ⊠ NewSU2Irrep, NewSU2Irrep ⊠ SU2Irrep, SU2Irrep ⊠
+              NewSU2Irrep, Z2Irrep ⊠ FibonacciAnyon ⊠ FibonacciAnyon)
 
 Ti = time()
 include("sectors.jl")