diff --git a/docs/src/lib/tensors.md b/docs/src/lib/tensors.md index 4d1c5c9cc..676a323b7 100644 --- a/docs/src/lib/tensors.md +++ b/docs/src/lib/tensors.md @@ -212,6 +212,17 @@ contract! ⊗(::AbstractTensorMap, ::AbstractTensorMap) ``` +### Planar tensor contractions + +For fermionic and anyonic tensors, planar contractions are needed to correctly handle +non-trivial braiding statistics. See the manual section on +[Fermionic tensor contractions](@ref) and [Anyonic tensor contractions](@ref) for +detailed information and examples. + +The `@planar` macro performs planar tensor contractions that respect non-trivial braiding, +while `@plansor` automatically dispatches between standard and planar contractions based +on the `BraidingStyle` of the tensors. + ## `TensorMap` factorizations The factorisation methods are powered by [MatrixAlgebraKit.jl](https://github.com/QuantumKitHub/MatrixAlgebraKit.jl) diff --git a/docs/src/man/tensors.md b/docs/src/man/tensors.md index 5a8aaf45a..d358b120d 100644 --- a/docs/src/man/tensors.md +++ b/docs/src/man/tensors.md @@ -1265,8 +1265,158 @@ Some examples: ## Fermionic tensor contractions -TODO +When working with fermionic tensors, i.e., tensors with `BraidingStyle(sectortype(tensor)) isa Fermionic`, +the non-trivial braiding statistics of fermionic particles must be accounted for in tensor +contractions. Unlike bosonic tensors where the order of operations doesn't matter (up to +trivial permutations), fermionic tensors require careful tracking of how tensor factors are +reordered and how indices are permuted. + +### The `@planar` macro + +For fermionic tensors, use the `@planar` macro instead of `@tensor`. The `@planar` macro +correctly handles the signs that arise from reordering fermionic indices: + +```julia +# Example with fermionic tensors +V = Vect[FermionParity](0 => 2, 1 => 2) # Fermionic vector space +A = randn(V ← V) +B = randn(V ← V') +C = randn(V' ← V) + +# Use @planar for fermionic contractions +@planar D[i; j] := A[i; j] + B[i; k] * C[k; j] +``` + +The key difference from `@tensor` is that `@planar` uses planar operations (`planaradd!`, +`planarcontract!`, `planartrace!`) that respect the fermionic braiding when permuting indices. +These operations automatically insert the appropriate signs (``-1`` raised to the appropriate +power) when indices with odd fermion parity are exchanged. + +### Braiding tensors + +The braiding of fermionic indices can be made explicit using braiding tensors, denoted by `τ`. +A braiding tensor `τ[a b; c d]` represents the operation that permutes indices, inserting +the appropriate fermionic signs. The braiding tensor must always have exactly 2 input and 2 output indices. + +The braiding tensor `τ` is special: TensorKit.jl automatically constructs it based on the +spaces of the indices it connects. You don't need to define it yourself; simply use the +symbol `τ` in your `@planar` expressions. For examples of using braiding tensors, see the +test suite in `test/planar.jl`. + +### Key differences from bosonic contractions + +1. **Order matters**: In `@tensor`, the order of factors doesn't affect the result (modulo + permutations). In `@planar`, reordering factors can introduce minus signs. + +2. **Index permutations**: When indices are permuted, `@planar` automatically accounts for + the fermionic statistics, inserting signs when odd-parity indices are exchanged. + +3. **Performance**: `@planar` has some overhead compared to `@tensor` due to the additional + bookkeeping required for braiding. For bosonic tensors, stick to `@tensor`. + +For more details on the mathematical framework underlying fermionic tensor contractions, +see the paper [arXiv:2404.14611](https://arxiv.org/abs/2404.14611), which describes the +implementation of fermionic tensor networks in TensorKit.jl. ## Anyonic tensor contractions -TODO +Anyonic tensors arise in the study of quantum systems with anyonic excitations, such as +fractional quantum Hall systems or topological quantum computing. In these systems, the +braiding of anyonic particles is described by a non-trivial braiding (or R-matrix), and +the tensor contractions must respect this braiding structure. + +### Using `@planar` for anyonic tensors + +Just like fermionic tensors, anyonic tensors require the use of `@planar` instead of `@tensor`: + +```julia +# Example with Fibonacci anyons +V = Vect[FibonacciAnyon](:I => 2, :τ => 2) +A = randn(V ← V) +B = randn(V ← V') +C = randn(V' ← V) + +# Use @planar for anyonic contractions +@planar D[i; j] := A[i; j] + B[i; k] * C[k; j] +``` + +For anyonic systems, the braiding is even more general than for fermions: instead of just +signs, the braiding tensors can have non-trivial matrix elements that depend on the +fusion-splitting structure of the anyons. The `@planar` macro handles all of this +automatically through the braiding tensor `τ`. + +### Anyonic braiding tensors + +The braiding tensors for anyons encode the non-trivial R-matrices (or F-matrices for the +fusion structure) of the anyonic theory. Similar to the fermionic case, these can be used +explicitly: + +```julia +# Explicit braiding in anyonic contractions +V = Vect[FibonacciAnyon](:I => 2, :τ => 2) +A = randn(V ← V ⊗ V) +B = randn(V ⊗ V ← V) +# Using braiding tensor to permute indices +@planar C[i; j] := A[i; k l] * τ[k l; m n] * B[m n; j] +``` + +The matrix elements of `τ` are determined by the specific anyonic theory (e.g., Fibonacci, +Ising, etc.) and are computed automatically by TensorKit.jl based on the sector type. + +### Differences from bosonic and fermionic contractions + +1. **General braiding**: Anyonic braiding is the most general case, with fermionic being + a special case (ℤ₂-graded braiding) and bosonic being the trivial case. + +2. **Fusion rules**: Anyonic tensors have more complex fusion rules than fermions or bosons, + which affects how tensor contractions are computed. + +3. **Computational cost**: Anyonic tensor contractions can be more expensive than fermionic + or bosonic contractions due to the more complex braiding structure. + +## Using `@plansor` for generic code + +When writing generic code that should work for bosonic, fermionic, and anyonic tensors, +use the `@plansor` macro. This macro automatically dispatches to the appropriate +implementation based on the `BraidingStyle` of the tensors: + +```julia +function generic_contraction(A::AbstractTensorMap, B::AbstractTensorMap) + @plansor C[i; j] := A[i; k] * B[k; j] + return C +end + +# Works correctly for all three cases: +# 1. Bosonic tensors: dispatches to @tensor (fast) +V_bosonic = ℂ^2 +A_bosonic = randn(V_bosonic ← V_bosonic) +B_bosonic = randn(V_bosonic ← V_bosonic) +C_bosonic = generic_contraction(A_bosonic, B_bosonic) + +# 2. Fermionic tensors: dispatches to @planar (correct signs) +V_fermionic = Vect[FermionParity](0 => 2, 1 => 2) +A_fermionic = randn(V_fermionic ← V_fermionic) +B_fermionic = randn(V_fermionic ← V_fermionic) +C_fermionic = generic_contraction(A_fermionic, B_fermionic) + +# 3. Anyonic tensors: dispatches to @planar (correct braiding) +V_anyonic = Vect[FibonacciAnyon](:I => 2, :τ => 2) +A_anyonic = randn(V_anyonic ← V_anyonic) +B_anyonic = randn(V_anyonic ← V_anyonic) +C_anyonic = generic_contraction(A_anyonic, B_anyonic) +``` + +The `@plansor` macro incurs a small runtime overhead for bosonic tensors due to the dispatch, +but ensures correctness for all braiding types. + +## Summary: `@tensor` vs `@planar` vs `@plansor` + +| Macro | Use case | Braiding handling | Performance | +|-------|----------|-------------------|-------------| +| `@tensor` | Bosonic tensors (no symmetries or abelian symmetries) | Ignores braiding (assumes trivial) | Fastest | +| `@planar` | Fermionic or anyonic tensors | Correctly handles non-trivial braiding | Slower due to braiding | +| `@plansor` | Generic code (unknown braiding type) | Runtime dispatch based on `BraidingStyle` | Small dispatch overhead | + +**Recommendation**: Use `@tensor` for bosonic tensors in performance-critical code, +`@planar` for fermionic/anyonic tensors, and `@plansor` for generic library code or when +the braiding type is not known at compile time. diff --git a/src/planar/macros.jl b/src/planar/macros.jl index c2afcb596..1bebfda2c 100644 --- a/src/planar/macros.jl +++ b/src/planar/macros.jl @@ -1,4 +1,71 @@ # NEW MACROS: @planar and @plansor +""" + @planar [kwargs...] expr + +Perform a planar tensor contraction, i.e., a tensor contraction that respects the +non-trivial braiding of fermionic or anyonic tensor indices. This macro is the +planar counterpart of the `@tensor` macro from TensorOperations.jl. + +The `@planar` macro should be used when working with tensors that have fermionic or +anyonic symmetries (i.e., `BraidingStyle(sectortype(tensor)) isa Fermionic` or +`BraidingStyle(sectortype(tensor)) isa Anyonic`). For such tensors, the order of +tensor contractions and index permutations matters due to non-trivial braiding, and +the `@planar` macro correctly handles the braiding tensors `τ` that encode these +braidings. + +# Syntax + +The syntax is identical to `@tensor`, with the addition that braiding tensors can be +explicitly specified using `τ[i j; k l]` to denote a braiding between indices. + +```julia +@planar C[i; j] := A[i; k] * B[k; j] +@planar C[i; j] := A[i; k] * τ[k l; m n] * B[m; n] * τ[n o; j l] +``` + +# Keyword arguments + +The same keyword arguments as `@tensor` are supported: +- `opt=true` or `opt=optexpr`: optimize the contraction order +- `backend=backend`: specify the backend for the contraction +- `allocator=allocator`: specify the allocator for temporary tensors +- `order=(...)`: manually specify the contraction order +- `contractcheck=true`: add runtime checks for contraction validity + +# Key differences from `@tensor` + +While `@tensor` assumes trivial braiding (appropriate for bosonic tensors), `@planar` +correctly handles non-trivial braiding statistics. Specifically: + +1. `@tensor` ignores the order of tensor factors in a product (they commute) +2. `@planar` respects the order and handles explicit braiding tensors `τ` +3. `@planar` uses planar operations (`planaradd!`, `planarcontract!`, etc.) that + account for braiding when permuting indices + +For bosonic tensors (trivial braiding), `@planar` reduces to the same operations as +`@tensor`, but with potential performance overhead. Use `@plansor` to automatically +dispatch to the appropriate implementation based on the braiding style. + +# Example + +```julia +# Fermionic tensors +V = Vect[FermionParity](0 => 2, 1 => 2) +A = randn(V ← V) +B = randn(V ← V') +C = randn(V' ← V) + +# Correct handling of fermionic statistics +@planar D[i; j] := A[i; j] + B[i; k] * C[k; j] + +# With explicit braiding tensor +A2 = randn(V ← V ⊗ V) +B2 = randn(V ⊗ V ← V) +@planar C2[i; j] := A2[i; k l] * τ[k l; m n] * B2[m n; j] +``` + +See also: [`@plansor`](@ref), [`@tensor`](@ref) +""" macro planar(args::Vararg{Expr}) isempty(args) && throw(ArgumentError("No arguments passed to `planar`")) @@ -91,6 +158,76 @@ function planarparser(planarexpr, kwargs...) return parser end +""" + @plansor [kwargs...] expr + +Automatically dispatch between `@tensor` and `@planar` based on the braiding style +of the tensors involved in the contraction. + +This macro provides a unified interface for tensor contractions that works correctly +for both bosonic tensors (trivial braiding) and fermionic/anyonic tensors (non-trivial +braiding). At runtime, it checks `BraidingStyle(sectortype(tensor))` and: + +- If `Bosonic`: uses the standard `@tensor` implementation (efficient, ignores braiding) +- If `Fermionic` or `Anyonic`: uses the `@planar` implementation (handles braiding correctly) + +# Syntax + +The syntax is the same as `@tensor` and `@planar`: + +```julia +@plansor C[i; j] := A[i; k] * B[k; j] +@plansor C[i; j] := A[i; k] * τ[k l; m n] * B[m; n] * τ[n o; j l] +``` + +# Keyword arguments + +The same keyword arguments as `@tensor` and `@planar` are supported: +- `opt=true` or `opt=optexpr`: optimize the contraction order +- `backend=backend`: specify the backend for the contraction +- `allocator=allocator`: specify the allocator for temporary tensors +- `order=(...)`: manually specify the contraction order +- `contractcheck=true`: add runtime checks for contraction validity + +# When to use `@plansor` vs `@tensor` vs `@planar` + +- Use `@tensor` when you know your tensors have trivial braiding (bosonic, or no symmetries) + and want maximum performance +- Use `@planar` when you know your tensors have non-trivial braiding (fermionic, anyonic) + and need correct handling of braiding +- Use `@plansor` when writing generic code that should work for both cases, or when you're + unsure about the braiding properties of your tensors + +# Performance considerations + +For bosonic tensors, `@plansor` has a small runtime overhead compared to `@tensor` due to +the runtime dispatch. For performance-critical code with known bosonic tensors, prefer +`@tensor`. For fermionic/anyonic tensors, `@plansor` and `@planar` have identical performance. + +# Example + +```julia +# Works correctly for both bosonic and fermionic tensors +function my_contraction(A::AbstractTensorMap, B::AbstractTensorMap) + @plansor C[i; j] := A[i; k] * B[k; j] + return C +end + +# Bosonic case (uses @tensor internally) +V_bosonic = ℂ^2 +A_bosonic = randn(V_bosonic ← V_bosonic) +B_bosonic = randn(V_bosonic ← V_bosonic) +C_bosonic = my_contraction(A_bosonic, B_bosonic) + +# Fermionic case (uses @planar internally) +V_fermionic = Vect[FermionParity](0 => 2, 1 => 2) +A_fermionic = randn(V_fermionic ← V_fermionic) +B_fermionic = randn(V_fermionic ← V_fermionic) +C_fermionic = my_contraction(A_fermionic, B_fermionic) +``` + +See also: [`@planar`](@ref), [`@tensor`](@ref) +""" macro plansor(args::Vararg{Expr}) isempty(args) && throw(ArgumentError("No arguments passed to `planar`"))