Symmetry sectors and fusion trees
Type hierarchy
TensorKitSectors.Sector
— Typeabstract type Sector end
Abstract type for representing the (isomorphism classes of) simple objects in (unitary and pivotal) (pre-)fusion categories, e.g. the irreducible representations of a finite or compact group. Subtypes I<:Sector
as the set of labels of a GradedSpace
.
Every new I<:Sector
should implement the following methods:
one(::Type{I})
: unit element ofI
conj(a::I)
: $a̅$, conjugate or dual label of $a$⊗(a::I, b::I)
: iterable with unique fusion outputs of $a ⊗ b$ (i.e. don't repeat in case of multiplicities)Nsymbol(a::I, b::I, c::I)
: number of timesc
appears ina ⊗ b
, i.e. the multiplicityFusionStyle(::Type{I})
:UniqueFusion()
,SimpleFusion()
orGenericFusion()
BraidingStyle(::Type{I})
:Bosonic()
,Fermionic()
,Anyonic()
, ...Fsymbol(a::I, b::I, c::I, d::I, e::I, f::I)
: F-symbol: scalar (in case ofUniqueFusion
/SimpleFusion
) or matrix (in case ofGenericFusion
)Rsymbol(a::I, b::I, c::I)
: R-symbol: scalar (in case ofUniqueFusion
/SimpleFusion
) or matrix (in case ofGenericFusion
)
and optionally
dim(a::I)
: quantum dimension of sectora
frobeniusschur(a::I)
: Frobenius-Schur indicator ofa
Bsymbol(a::I, b::I, c::I)
: B-symbol: scalar (in case ofUniqueFusion
/SimpleFusion
) or matrix (in case ofGenericFusion
)twist(a::I)
-> twist of sectora
Furthermore, iterate
and Base.IteratorSize
should be made to work for the singleton type SectorValues{I}
.
TensorKitSectors.SectorValues
— Typestruct SectorValues{I<:Sector}
Singleton type to represent an iterator over the possible values of type I
, whose instance is obtained as values(I)
. For a new I::Sector
, the following should be defined
Base.iterate(::SectorValues{I}[, state])
: iterate over the valuesBase.IteratorSize(::Type{SectorValues{I}})
:HasLenght()
,SizeUnkown()
orIsInfinite()
depending on whether the number of values of typeI
is finite (and sufficiently small) or infinite; for a large number of values,SizeUnknown()
is recommend because this will trigger the use ofGenericGradedSpace
.
If IteratorSize(I) == HasLength()
, also the following must be implemented:
Base.length(::SectorValues{I})
: the number of different valuesBase.getindex(::SectorValues{I}, i::Int)
: a mapping between an indexi
and an instance ofI
findindex(::SectorValues{I}, c::I)
: reverse mapping between a valuec::I
and an indexi::Integer ∈ 1:length(values(I))
TensorKitSectors.FusionStyle
— TypeFusionStyle(::Sector)
FusionStyle(I::Type{<:Sector})
Trait to describe the fusion behavior of sectors of type I
, which can be either
UniqueFusion()
: single fusion output when fusing two sectors;SimpleFusion()
: multiple outputs, but every output occurs at most one, also known as multiplicity-free (e.g. irreps of $SU(2)$);GenericFusion()
: multiple outputs that can occur more than once (e.g. irreps of $SU(3)$).
There is an abstract supertype MultipleFusion
of which both SimpleFusion
and GenericFusion
are subtypes. Furthermore, there is a type alias MultiplicityFreeFusion
for those fusion types which do not require muliplicity labels, i.e. MultiplicityFreeFusion = Union{UniqueFusion,SimpleFusion}
.
TensorKitSectors.BraidingStyle
— TypeBraidingStyle(::Sector) -> ::BraidingStyle
BraidingStyle(I::Type{<:Sector}) -> ::BraidingStyle
Return the type of braiding and twist behavior of sectors of type I
, which can be either
Bosonic()
: symmetric braiding with trivial twist (i.e. identity)Fermionic()
: symmetric braiding with non-trivial twist (squares to identity)Anyonic()
: general $R_(a,b)^c$ phase or matrix (depending onSimpleFusion
orGenericFusion
fusion) and arbitrary twists
Note that Bosonic
and Fermionic
are subtypes of SymmetricBraiding
, which means that braids are in fact equivalent to crossings (i.e. braiding twice is an identity: isone(Rsymbol(b,a,c)*Rsymbol(a,b,c)) == true
) and permutations are uniquely defined.
TensorKitSectors.AbstractIrrep
— Typeabstract type AbstractIrrep{G<:Group} <: Sector end
Abstract supertype for sectors which corresponds to irreps (irreducible representations) of a group G
. As we assume unitary representations, these would be finite groups or compact Lie groups. Note that this could also include projective rather than linear representations.
Actual concrete implementations of those irreps can be obtained as Irrep[G]
, or via their actual name, which generically takes the form (asciiG)Irrep
, i.e. the ASCII spelling of the group name followed by Irrep
.
All irreps have BraidingStyle
equal to Bosonic()
and thus trivial twists.
TensorKitSectors.Trivial
— TypeTrivial
Singleton type to represent the trivial sector, i.e. the trivial representation of the trivial group. This is equivalent to Rep[ℤ₁]
, or the unit object of the category Vect
of ordinary vector spaces.
TensorKitSectors.ZNIrrep
— Typestruct ZNIrrep{N} <: AbstractIrrep{ℤ{N}}
ZNIrrep{N}(n::Integer)
Irrep[ℤ{N}](n::Integer)
Represents irreps of the group $ℤ_N$ for some value of N<64
. (We need 2*(N-1) <= 127 in order for a ⊗ b to work correctly.) For N
equals 2
, 3
or 4
, ℤ{N}
can be replaced by ℤ₂
, ℤ₃
, ℤ₄
. An arbitrary Integer
n
can be provided to the constructor, but only the value mod(n, N)
is relevant.
Fields
n::Int8
: the integer label of the irrep, moduloN
.
TensorKitSectors.U1Irrep
— Typestruct U1Irrep <: AbstractIrrep{U₁}
U1Irrep(charge::Real)
Irrep[U₁](charge::Real)
Represents irreps of the group $U₁$. The irrep is labelled by a charge, which should be an integer for a linear representation. However, it is often useful to allow half integers to represent irreps of $U₁$ subgroups of $SU₂$, such as the Sz of spin-1/2 system. Hence, the charge is stored as a HalfInt
from the package HalfIntegers.jl, but can be entered as arbitrary Real
. The sequence of the charges is: 0, 1/2, -1/2, 1, -1, ...
Fields
charge::HalfInt
: the label of the irrep, which can be any half integer.
TensorKitSectors.SU2Irrep
— Typestruct SU2Irrep <: AbstractIrrep{SU₂}
SU2Irrep(j::Real)
Irrep[SU₂](j::Real)
Represents irreps of the group $SU₂$. The irrep is labelled by a half integer j
which can be entered as an abitrary Real
, but is stored as a HalfInt
from the HalfIntegers.jl package.
Fields
j::HalfInt
: the label of the irrep, which can be any non-negative half integer.
TensorKitSectors.CU1Irrep
— Typestruct CU1Irrep <: AbstractIrrep{CU₁}
CU1Irrep(j, s = ifelse(j>zero(j), 2, 0))
Irrep[CU₁](j, s = ifelse(j>zero(j), 2, 0))
Represents irreps of the group $U₁ ⋊ C$ ($U₁$ and charge conjugation or reflection), which is also known as just O₂
.
Fields
j::HalfInt
: the value of the $U₁$ charge.s::Int
: the representation of charge conjugation.
They can take values:
- if
j == 0
,s = 0
(trivial charge conjugation) ors = 1
(non-trivial charge conjugation) - if
j > 0
,s = 2
(two-dimensional representation)
TensorKitSectors.ProductSector
— TypeProductSector{T<:SectorTuple}
Represents the Deligne tensor product of sectors. The type parameter T
is a tuple of the component sectors. The recommended way to construct a ProductSector
is using the deligneproduct
(⊠
) operator on the components.
TensorKitSectors.FermionParity
— TypeFermionParity <: Sector
Represents sectors with fermion parity. The fermion parity is a ℤ₂ quantum number that yields an additional sign when two odd fermions are exchanged.
Fields
isodd::Bool
: indicates whether the fermion parity is odd (true
) or even (false
).
See also: FermionNumber
, FermionSpin
TensorKitSectors.FermionNumber
— Typeconst FermionNumber = U1Irrep ⊠ FermionParity
FermionNumber(a::Int)
Represents the fermion number as the direct product of a $U₁$ irrep a
and a fermion parity, with the restriction that the fermion parity is odd if and only if a
is odd.
See also: U1Irrep
, FermionParity
TensorKitSectors.FermionSpin
— Typeconst FermionSpin = SU2Irrep ⊠ FermionParity
FermionSpin(j::Real)
Represents the fermion spin as the direct product of a $SU₂$ irrep j
and a fermion parity, with the restriction that the fermion parity is odd if 2 * j
is odd.
See also: SU2Irrep
, FermionParity
TensorKitSectors.FibonacciAnyon
— Typestruct FibonacciAnyon <: Sector
FibonacciAnyon(s::Symbol)
Represents the anyons of the Fibonacci modular fusion category. It can take two values, corresponding to the trivial sector FibonacciAnyon(:I)
and the non-trivial sector FibonacciAnyon(:τ)
with fusion rules $τ ⊗ τ = 1 ⊕ τ$.
Fields
isone::Bool
: indicates whether the sector corresponds the to trivial anyon:I
(true
), or the non-trivial anyon:τ
(false
).
TensorKitSectors.IsingAnyon
— Typestruct IsingAnyon <: Sector
IsingAnyon(s::Symbol)
Represents the anyons of the Ising modular 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 ⊕ ψ$.
Fields
s::Symbol
: the label of the represented anyon, which can be:I
,:σ
, or:ψ
.
Useful constants
TensorKitSectors.Irrep
— Constantconst Irrep
A constant of a singleton type used as Irrep[G]
with G<:Group
a type of group, to construct or obtain a concrete subtype of AbstractIrrep{G}
that implements the data structure used to represent irreducible representations of the group G
.
Methods for defining and characterizing Sector
subtypes
Base.one
— Methodone(::Sector) -> Sector
one(::Type{<:Sector}) -> Sector
Return the unit element within this type of sector.
TensorKitSectors.dual
— Methoddual(a::Sector) -> Sector
Return the conjugate label conj(a)
.
TensorKitSectors.Nsymbol
— FunctionNsymbol(a::I, b::I, c::I) where {I<:Sector} -> Integer
Return an Integer
representing the number of times c
appears in the fusion product a ⊗ b
. Could be a Bool
if FusionStyle(I) == UniqueFusion()
or SimpleFusion()
.
TensorKitSectors.:⊗
— Function⊗(a::I, b::I...) where {I<:Sector}
otimes(a::I, b::I...) where {I<:Sector}
Return an iterable of elements of c::I
that appear in the fusion product a ⊗ b
.
Note that every element c
should appear at most once, fusion degeneracies (if FusionStyle(I) == GenericFusion()
) should be accessed via Nsymbol(a, b, c)
.
TensorKitSectors.Fsymbol
— FunctionFsymbol(a::I, b::I, c::I, d::I, e::I, f::I) where {I<:Sector}
Return the F-symbol $F^{abc}_d$ that associates the two different fusion orders of sectors a
, b
and c
into an ouput sector d
, using either an intermediate sector $a ⊗ b → e$ or $b ⊗ c → f$:
a-<-μ-<-e-<-ν-<-d a-<-λ-<-d
∨ ∨ -> Fsymbol(a,b,c,d,e,f)[μ,ν,κ,λ] ∨
b c f
v
b-<-κ
∨
c
If FusionStyle(I)
is UniqueFusion
or SimpleFusion
, the F-symbol is a number. Otherwise it is a rank 4 array of size (Nsymbol(a, b, e), Nsymbol(e, c, d), Nsymbol(b, c, f), Nsymbol(a, f, d))
.
TensorKitSectors.Rsymbol
— FunctionRsymbol(a::I, b::I, c::I) where {I<:Sector}
Returns the R-symbol $R^{ab}_c$ that maps between $c → a ⊗ b$ and $c → b ⊗ a$ as in
a -<-μ-<- c b -<-ν-<- c
∨ -> Rsymbol(a,b,c)[μ,ν] v
b a
If FusionStyle(I)
is UniqueFusion()
or SimpleFusion()
, the R-symbol is a number. Otherwise it is a square matrix with row and column size Nsymbol(a,b,c) == Nsymbol(b,a,c)
.
TensorKitSectors.Bsymbol
— FunctionBsymbol(a::I, b::I, c::I) where {I<:Sector}
Return the value of $B^{ab}_c$ which appears in transforming a splitting vertex into a fusion vertex using the transformation
a -<-μ-<- c a -<-ν-<- c
∨ -> √(dim(c)/dim(a)) * Bsymbol(a,b,c)[μ,ν] ∧
b dual(b)
If FusionStyle(I)
is UniqueFusion()
or SimpleFusion()
, the B-symbol is a number. Otherwise it is a square matrix with row and column size Nsymbol(a, b, c) == Nsymbol(c, dual(b), a)
.
TensorKitSectors.dim
— Methoddim(a::Sector)
Return the (quantum) dimension of the sector a
.
TensorKitSectors.frobeniusschur
— Functionfrobeniusschur(a::Sector)
Return the Frobenius-Schur indicator of a sector a
.
TensorKitSectors.twist
— Methodtwist(a::Sector)
Return the twist of a sector a
Base.isreal
— Methodisreal(::Type{<:Sector}) -> Bool
Return whether the topological data (Fsymbol, Rsymbol) of the sector is real or not (in which case it is complex).
TensorKitSectors.sectorscalartype
— Functionsectorscalartype(I::Type{<:Sector}) -> Type
Return the scalar type of the topological data (Fsymbol, Rsymbol) of the sector I
.
TensorKitSectors.deligneproduct
— Method⊠(s₁::Sector, s₂::Sector)
deligneproduct(s₁::Sector, s₂::Sector)
Given two sectors s₁
and s₂
, which label an isomorphism class of simple objects in a fusion category $C₁$ and $C₂$, s1 ⊠ s2
(obtained as \boxtimes+TAB
) labels the isomorphism class of simple objects in the Deligne tensor product category $C₁ ⊠ C₂$.
The Deligne tensor product also works in the type domain and for spaces and tensors. For group representations, we have Irrep[G₁] ⊠ Irrep[G₂] == Irrep[G₁ × G₂]
.
Compile all revelant methods for a sector:
TensorKitSectors.precompile_sector
— Functionprecompile_sector(I::Type{<:Sector})
Precompile common methods for the given sector type.
Types and methods for groups
Types and constants:
# TODO: add documentation for the following types
Group
TensorKitSectors.AbelianGroup
U₁
ℤ{N} where N
SU{N} where N
const SU₂ = SU{2}
ProductGroup
Specific methods:
TensorKitSectors.:×
— Function×(G::Vararg{Type{<:Group}}) -> ProductGroup{Tuple{G...}}
times(G::Vararg{Type{<:Group}}) -> ProductGroup{Tuple{G...}}
Construct the direct product of a (list of) groups.
Methods for defining and generating fusion trees
TensorKit.FusionTree
— Typestruct FusionTree{I, N, M, L}
Represents a fusion tree of sectors of type I<:Sector
, fusing (or splitting) N
uncoupled sectors to a coupled sector. It actually represents a splitting tree, but fusion tree is a more common term.
Fields
uncoupled::NTuple{N,I}
: the uncoupled sectors coming out of the splitting tree, before the possible 𝑍 isomorphism (seeisdual
).coupled::I
: the coupled sector.isdual::NTuple{N,Bool}
: indicates whether a 𝑍 isomorphism is present (true
) or not (false
) for each uncoupled sector.innerlines::NTuple{M,I}
: the labels of theM=max(0, N-2)
inner lines of the splitting tree.vertices::NTuple{L,Int}
: the integer values of theL=max(0, N-1)
vertices of the splitting tree. IfFusionStyle(I) isa MultiplicityFreeFusion
, thenvertices
is simply equal to the constant valuentuple(n->1, L)
.
TensorKit.fusiontrees
— Methodfusiontrees(uncoupled::NTuple{N,I}[,
coupled::I=one(I)[, isdual::NTuple{N,Bool}=ntuple(n -> false, length(uncoupled))]])
where {N,I<:Sector} -> FusionTreeIterator{I,N,I}
Return an iterator over all fusion trees with a given coupled sector label coupled
and uncoupled sector labels and isomorphisms uncoupled
and isdual
respectively.
Methods for manipulating fusion trees
For manipulating single fusion trees, the following internal methods are defined:
TensorKit.insertat
— Functioninsertat(f::FusionTree{I, N₁}, i::Int, f₂::FusionTree{I, N₂})
-> <:AbstractDict{<:FusionTree{I, N₁+N₂-1}, <:Number}
Attach a fusion tree f₂
to the uncoupled leg i
of the fusion tree f₁
and bring it into a linear combination of fusion trees in standard form. This requires that f₂.coupled == f₁.uncoupled[i]
and f₁.isdual[i] == false
.
TensorKit.split
— Functionsplit(f::FusionTree{I, N}, M::Int)
-> (::FusionTree{I, M}, ::FusionTree{I, N-M+1})
Split a fusion tree into two. The first tree has as uncoupled sectors the first M
uncoupled sectors of the input tree f
, whereas its coupled sector corresponds to the internal sector between uncoupled sectors M
and M+1
of the original tree f
. The second tree has as first uncoupled sector that same internal sector of f
, followed by remaining N-M
uncoupled sectors of f
. It couples to the same sector as f
. This operation is the inverse of insertat
in the sense that if f₁, f₂ = split(t, M) ⇒ f == insertat(f₂, 1, f₁)
.
TensorKit.merge
— Functionmerge(f₁::FusionTree{I, N₁}, f₂::FusionTree{I, N₂}, c::I, μ = 1)
-> <:AbstractDict{<:FusionTree{I, N₁+N₂}, <:Number}
Merge two fusion trees together to a linear combination of fusion trees whose uncoupled sectors are those of f₁
followed by those of f₂
, and where the two coupled sectors of f₁
and f₂
are further fused to c
. In case of FusionStyle(I) == GenericFusion()
, also a degeneracy label μ
for the fusion of the coupled sectors of f₁
and f₂
to c
needs to be specified.
TensorKit.elementary_trace
— Functionelementary_trace(f::FusionTree{I,N}, i) where {I<:Sector,N} -> <:AbstractDict{FusionTree{I,N-2}, <:Number}
Perform an elementary trace of neighbouring uncoupled indices i
and i+1
on a fusion tree f
, and returns the result as a dictionary of output trees and corresponding coefficients.
TensorKit.planar_trace
— Methodplanar_trace(f::FusionTree{I,N}, q1::IndexTuple{N₃}, q2::IndexTuple{N₃}) where {I<:Sector,N,N₃}
-> <:AbstractDict{FusionTree{I,N-2*N₃}, <:Number}
Perform a planar trace of the uncoupled indices of the fusion tree f
at q1
with those at q2
, where q1[i]
is connected to q2[i]
for all i
. The result is returned as a dictionary of output trees and corresponding coefficients.
TensorKit.artin_braid
— Functionartin_braid(f::FusionTree, i; inv::Bool = false) -> <:AbstractDict{typeof(f), <:Number}
Perform an elementary braid (Artin generator) of neighbouring uncoupled indices i
and i+1
on a fusion tree f
, and returns the result as a dictionary of output trees and corresponding coefficients.
The keyword inv
determines whether index i
will braid above or below index i+1
, i.e. applying artin_braid(f′, i; inv = true)
to all the outputs f′
of artin_braid(f, i; inv = false)
and collecting the results should yield a single fusion tree with non-zero coefficient, namely f
with coefficient 1
. This keyword has no effect if BraidingStyle(sectortype(f)) isa SymmetricBraiding
.
TensorKit.braid
— Methodbraid(f::FusionTree{<:Sector, N}, levels::NTuple{N, Int}, p::NTuple{N, Int})
-> <:AbstractDict{typeof(t), <:Number}
Perform a braiding of the uncoupled indices of the fusion tree f
and return the result as a <:AbstractDict
of output trees and corresponding coefficients. The braiding is determined by specifying that the new sector at position k
corresponds to the sector that was originally at the position i = p[k]
, and assigning to every index i
of the original fusion tree a distinct level or depth levels[i]
. This permutation is then decomposed into elementary swaps between neighbouring indices, where the swaps are applied as braids such that if i
and j
cross, $τ_{i,j}$ is applied if levels[i] < levels[j]
and $τ_{j,i}^{-1}$ if levels[i] > levels[j]
. This does not allow to encode the most general braid, but a general braid can be obtained by combining such operations.
TensorKit.permute
— Methodpermute(f::FusionTree, p::NTuple{N, Int}) -> <:AbstractDict{typeof(t), <:Number}
Perform a permutation of the uncoupled indices of the fusion tree f
and returns the result as a <:AbstractDict
of output trees and corresponding coefficients; this requires that BraidingStyle(sectortype(f)) isa SymmetricBraiding
.
These can be composed to implement elementary manipulations of fusion-splitting tree pairs, according to the following methods
# TODO: add documentation for the following methods
TensorKit.bendright
TensorKit.bendleft
TensorKit.foldright
TensorKit.foldleft
TensorKit.cycleclockwise
TensorKit.cycleanticlockwise
Finally, these are used to define large manipulations of fusion-splitting tree pairs, which are then used in the index manipulation of AbstractTensorMap
objects. The following methods defined on fusion splitting tree pairs have an associated definition for tensors.
TensorKit.repartition
— Functionrepartition(f₁::FusionTree{I, N₁}, f₂::FusionTree{I, N₂}, N::Int) where {I, N₁, N₂}
-> <:AbstractDict{Tuple{FusionTree{I, N}, FusionTree{I, N₁+N₂-N}}, <:Number}
Input is a double fusion tree that describes the fusion of a set of incoming uncoupled sectors to a set of outgoing uncoupled sectors, represented using the individual trees of outgoing (f₁
) and incoming sectors (f₂
) respectively (with identical coupled sector f₁.coupled == f₂.coupled
). Computes new trees and corresponding coefficients obtained from repartitioning the tree by bending incoming to outgoing sectors (or vice versa) in order to have N
outgoing sectors.
repartition(tsrc::AbstractTensorMap{S}, N₁::Int, N₂::Int; copy::Bool=false) where {S}
-> tdst::AbstractTensorMap{S,N₁,N₂}
Return tensor tdst
obtained by repartitioning the indices of t
. The codomain and domain of tdst
correspond to the first N₁
and last N₂
spaces of t
, respectively.
If copy=false
, tdst
might share data with tsrc
whenever possible. Otherwise, a copy is always made.
To repartition into an existing destination, see repartition!.
Base.transpose
— Methodtranspose(f₁::FusionTree{I}, f₂::FusionTree{I},
p1::NTuple{N₁, Int}, p2::NTuple{N₂, Int}) where {I, N₁, N₂}
-> <:AbstractDict{Tuple{FusionTree{I, N₁}, FusionTree{I, N₂}}, <:Number}
Input is a double fusion tree that describes the fusion of a set of incoming uncoupled sectors to a set of outgoing uncoupled sectors, represented using the individual trees of outgoing (t1
) and incoming sectors (t2
) respectively (with identical coupled sector t1.coupled == t2.coupled
). Computes new trees and corresponding coefficients obtained from repartitioning and permuting the tree such that sectors p1
become outgoing and sectors p2
become incoming.
TensorKit.braid
— Methodbraid(f₁::FusionTree{I}, f₂::FusionTree{I},
levels1::IndexTuple, levels2::IndexTuple,
p1::IndexTuple{N₁}, p2::IndexTuple{N₂}) where {I<:Sector, N₁, N₂}
-> <:AbstractDict{Tuple{FusionTree{I, N₁}, FusionTree{I, N₂}}, <:Number}
Input is a fusion-splitting tree pair that describes the fusion of a set of incoming uncoupled sectors to a set of outgoing uncoupled sectors, represented using the splitting tree f₁
and fusion tree f₂
, such that the incoming sectors f₂.uncoupled
are fused to f₁.coupled == f₂.coupled
and then to the outgoing sectors f₁.uncoupled
. Compute new trees and corresponding coefficients obtained from repartitioning and braiding the tree such that sectors p1
become outgoing and sectors p2
become incoming. The uncoupled indices in splitting tree f₁
and fusion tree f₂
have levels (or depths) levels1
and levels2
respectively, which determines how indices braid. In particular, if i
and j
cross, $τ_{i,j}$ is applied if levels[i] < levels[j]
and $τ_{j,i}^{-1}$ if levels[i] > levels[j]
. This does not allow to encode the most general braid, but a general braid can be obtained by combining such operations.
TensorKit.permute
— Methodpermute(f₁::FusionTree{I}, f₂::FusionTree{I},
p1::NTuple{N₁, Int}, p2::NTuple{N₂, Int}) where {I, N₁, N₂}
-> <:AbstractDict{Tuple{FusionTree{I, N₁}, FusionTree{I, N₂}}, <:Number}
Input is a double fusion tree that describes the fusion of a set of incoming uncoupled sectors to a set of outgoing uncoupled sectors, represented using the individual trees of outgoing (t1
) and incoming sectors (t2
) respectively (with identical coupled sector t1.coupled == t2.coupled
). Computes new trees and corresponding coefficients obtained from repartitioning and permuting the tree such that sectors p1
become outgoing and sectors p2
become incoming.