diff --git a/.github/workflows/PR_comment.yml b/.github/workflows/PR_comment.yml new file mode 100644 index 0000000..ea00f15 --- /dev/null +++ b/.github/workflows/PR_comment.yml @@ -0,0 +1,17 @@ +name: Docs preview comment +on: + pull_request: + types: [labeled] + +permissions: + pull-requests: write +jobs: + pr_comment: + runs-on: ubuntu-latest + steps: + - name: Create PR comment + if: github.event_name == 'pull_request' && github.repository == github.event.pull_request.head.repo.full_name && github.event.label.name == 'documentation' + uses: thollander/actions-comment-pull-request@v3 + with: + message: 'After the build completes, the updated documentation will be available [here](https://quantumkithub.github.io/MultiTensorKit.jl/previews/PR${{ github.event.number }}/)' + comment-tag: 'preview-doc' \ No newline at end of file diff --git a/.gitignore b/.gitignore index 10593a9..5f9fc64 100644 --- a/.gitignore +++ b/.gitignore @@ -11,4 +11,3 @@ Manifest.toml benchmark/*.json docs/Manifest.toml docs/build/ -docs/src/index.md diff --git a/Project.toml b/Project.toml index a542085..b40fb6a 100644 --- a/Project.toml +++ b/Project.toml @@ -5,15 +5,22 @@ version = "0.1.0" [deps] Artifacts = "56f22d72-fd6d-98f1-02f0-08ddc0907c33" +BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf" +BlockTensorKit = "5f87ffc2-9cf1-4a46-8172-465d160bd8cd" +DelimitedFiles = "8bb1440f-4735-579b-a4ab-409b98df4dab" JSON3 = "0f8b85d8-7281-11e9-16c2-39a750bddbf1" +MPSKit = "bb1c41ca-d63c-52ed-829e-0820dda26502" +MPSKitModels = "ca635005-6f8c-4cd1-b51d-8491250ef2ab" +TensorKit = "07d1fe3e-3e46-537d-9eac-e9e13d0d4cec" TensorKitSectors = "13a9c161-d5da-41f0-bcbd-e1a08ae0647f" +TestExtras = "5ed8adda-3752-4e41-b88a-e8b09835ee3a" +TupleTools = "9d95972d-f1c8-5527-a6e0-b4b365fa01f6" [compat] Aqua = "0.8.9" Artifacts = "1.10, 1" JSON3 = "1.14.1" SafeTestsets = "0.1" -TensorKitSectors = "0.1.2" Test = "1.10" TestExtras = "0.3" julia = "1.10" diff --git a/docs/Project.toml b/docs/Project.toml new file mode 100644 index 0000000..a36ea53 --- /dev/null +++ b/docs/Project.toml @@ -0,0 +1,4 @@ +[deps] +Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" +DocumenterCitations = "daee34ce-89f3-4625-b898-19384cb65244" +MultiTensorKit = "f0555a46-f681-4ef1-bed5-d64870568636" diff --git a/docs/make.jl b/docs/make.jl new file mode 100644 index 0000000..4e93ef2 --- /dev/null +++ b/docs/make.jl @@ -0,0 +1,26 @@ +using Documenter +using DocumenterCitations +using MultiTensorKit + +pages = ["Home" => "index.md", + "Manual" => ["man/fusioncats.md", "man/multifusioncats.md", "man/implementation.md"], + "Library" => "lib/library.md", + "References" => "references.md"] + +# bibliography +bibpath = joinpath(@__DIR__, "src", "assets", "MTKrefs.bib") +bib = CitationBibliography(bibpath; style=:authoryear) + +mathengine = MathJax3(Dict(:loader => Dict("load" => ["[tex]/physics"]), + :tex => Dict("inlineMath" => [["\$", "\$"], ["\\(", "\\)"]], + "tags" => "ams", + "packages" => ["base", "ams", "autoload", "physics", + "mathtools"]))) + +makedocs(; sitename="MultiTensorKit.jl", modules=[MultiTensorKit], + authors="Boris De Vos, Laurens Lootens and Lukas Devos", + pages=pages, pagesonly=true, plugins=[bib], + format=Documenter.HTML(; prettyurls=true, mathengine=mathengine, + assets=["assets/custom.css"])) + +deploydocs(; repo="https://github.com/QuantumKitHub/MultiTensorKit.jl", push_preview = true) \ No newline at end of file diff --git a/docs/src/assets/MTKrefs.bib b/docs/src/assets/MTKrefs.bib new file mode 100644 index 0000000..7a95ca2 --- /dev/null +++ b/docs/src/assets/MTKrefs.bib @@ -0,0 +1,53 @@ +@book{etingof2016tensor, + title={Tensor Categories}, + author={Etingof, P. and Gelaki, S. and Nikshych, D. and Ostrik, V.}, + isbn={9781470434410}, + lccn={2015006773}, + series={Mathematical Surveys and Monographs}, + url={https://books.google.be/books?id=Z6XLDAAAQBAJ}, + year={2016}, + publisher={American Mathematical Society} +} + +@article{Lootens_2023, + title={Dualities in One-Dimensional Quantum Lattice Models: Symmetric Hamiltonians and Matrix Product Operator Intertwiners}, + volume={4}, + ISSN={2691-3399}, + url={http://dx.doi.org/10.1103/PRXQuantum.4.020357}, + DOI={10.1103/prxquantum.4.020357}, + number={2}, + journal={PRX Quantum}, + publisher={American Physical Society (APS)}, + author={Lootens, Laurens and Delcamp, Clement and Ortiz, Gerardo and Verstraete, Frank}, + year={2023}, + month=jun } + +@misc{etingof2009, + title={Fusion categories and homotopy theory}, + author={Pavel Etingof and Dmitri Nikshych and Victor Ostrik and with an appendix by Ehud Meir}, + year={2009}, + eprint={0909.3140}, + archivePrefix={arXiv}, + primaryClass={math.QA}, + url={https://arxiv.org/abs/0909.3140}, +} + +@misc{henriques2020, + title={Representations of fusion categories and their commutants}, + author={AndrΓ© Henriques and David Penneys}, + year={2020}, + eprint={2004.08271}, + archivePrefix={arXiv}, + primaryClass={math.OA}, + url={https://arxiv.org/abs/2004.08271}, +} + +@misc{Lootens_2024, + title={Entanglement and the density matrix renormalisation group in the generalised Landau paradigm}, + author={Laurens Lootens and Clement Delcamp and Frank Verstraete}, + year={2024}, + eprint={2408.06334}, + archivePrefix={arXiv}, + primaryClass={quant-ph}, + url={https://arxiv.org/abs/2408.06334}, +} \ No newline at end of file diff --git a/docs/src/assets/custom.css b/docs/src/assets/custom.css new file mode 100644 index 0000000..4e38a20 --- /dev/null +++ b/docs/src/assets/custom.css @@ -0,0 +1,47 @@ +.center { + display: block; + margin-left: auto; + margin-right: auto; +} + +/* set fixed non-trivial inversion and hue rotation */ +:root { + --inversionFraction: 100%; + --hueRotation: -180deg; +} + +/* color-invert images */ +.color-invertible { + filter: invert(var(--inversionFraction)) hue-rotate(var(--hueRotation)) !important; + -ms-filter: invert(var(--inversionFraction)) !important; + -webkit-filter: invert(var(--inversionFraction)) hue-rotate(var(--hueRotation)) !important; +} +/* including the logo when we make one +.docs-logo { + filter: invert(var(--inversionFraction)) hue-rotate(var(--hueRotation)) !important; + -ms-filter: invert(var(--inversionFraction)) !important; + -webkit-filter: invert(var(--inversionFraction)) hue-rotate(var(--hueRotation)) !important; +} */ + +/* but disable for the two light themes */ +.theme--documenter-light .color-invertible { + filter: invert(0%) hue-rotate(0deg) !important; + -ms-filter: invert(var(0%)) !important; + -webkit-filter: invert(var(0%)) hue-rotate(0deg) !important; +} +.theme--catppuccin-latte .color-invertible { + filter: invert(0%) hue-rotate(0deg) !important; + -ms-filter: invert(var(0%)) !important; + -webkit-filter: invert(var(0%)) hue-rotate(0deg) !important; +} +/* for the logo as well */ +/* .theme--documenter-light .docs-logo { + filter: invert(0%) hue-rotate(0deg) !important; + -ms-filter: invert(var(0%)) !important; + -webkit-filter: invert(var(0%)) hue-rotate(0deg) !important; +} +.theme--catppuccin-latte .docs-logo { + filter: invert(0%) hue-rotate(0deg) !important; + -ms-filter: invert(var(0%)) !important; + -webkit-filter: invert(var(0%)) hue-rotate(0deg) !important; +} */ \ No newline at end of file diff --git a/docs/src/index.md b/docs/src/index.md new file mode 100644 index 0000000..63c1b2d --- /dev/null +++ b/docs/src/index.md @@ -0,0 +1,34 @@ +# MultiTensorKit + +**TensorKit extension to multifusion categories** + +## Package summary +MultiTensorKit.jl provides the user a package to work with multifusion categories, the extension of regular fusion categories where the unit is no longer simple and unique. +Multifusion categories naturally embed the structure of module categories over fusion categories. Hence, MultiTensorKit.jl allows not only the fusion of objects within the same +fusion category (as TensorKit.jl), but also the fusion with and between module categories over these fusion categories. + +MultiTensorKit.jl is built to be compatible with TensorKit, thus allowing the construction of symmetric tensors with new symmetries due to the module structure. Through this, +tensor network simulations of quantum many-body systems with aid of [MPSKit.jl](https://github.com/QuantumKitHub/MPSKit.jl) can be performed. + +## Table of contents + +```@contents +Pages = ["man/fusioncats.md","man/multifusioncats.md","lib/library.md", "references.md"] +Depth = 2 +``` + +## Installation + +MultiTensorKit.jl is currently not registered to the Julia General Registry. You can install the package as +``` +pkg> add https://github.com/QuantumKitHub/MultiTensorKit.jl.git +``` + +## Usage + +As the name suggests, MultiTensorKit is an extension of [TensorKit.jl](https://github.com/Jutho/TensorKit.jl) and +[TensorKitSectors.jl](https://github.com/QuantumKitHub/TensorKitSectors.jl). Therefore, we recommend including TensorKit +to your project. Additionally, MultiTensorKit was made to be functional with [MPSKit.jl](https://github.com/QuantumKitHub/MPSKit.jl) +and [MPSKitModels.jl](https://github.com/QuantumKitHub/MPSKitModels.jl) for Matrix Product State (MPS) calculations, supporting symmetries +which go beyond TensorKit. All these packages are registered in JuliaRegistries and can be added through the package manager. + diff --git a/docs/src/lib/library.md b/docs/src/lib/library.md new file mode 100644 index 0000000..86a8f6c --- /dev/null +++ b/docs/src/lib/library.md @@ -0,0 +1,5 @@ +# Library + +```@autodocs +Modules = [MultiTensorKit] +``` \ No newline at end of file diff --git a/docs/src/man/fusioncats.md b/docs/src/man/fusioncats.md new file mode 100644 index 0000000..f584678 --- /dev/null +++ b/docs/src/man/fusioncats.md @@ -0,0 +1,59 @@ +# Introduction + +The manual has been divided into different sections in an attempt to break down the information the user requires to use MultiTensorKit.jl. +We start off with a short summary of fusion category theory. Users familiar with TensorKit.jl may have read the [Optional introduction to category theory](https://jutho.github.io/TensorKit.jl/stable/man/categories/) in the documentation of TensorKit; this section can then largely be skipped. Be aware that notation may differ from the literature. + +Afterwards, the extension to multifusion categories is explained, and its relation to (bi)module categories over fusion categories is shown. + +# Fusion category theory + +The aim of this section is to explain the bare minimum required to proceed to the next section on multifusion category theory and bimodule categories. More details can be found in the [TensorKit](https://jutho.github.io/TensorKit.jl/stable/man/categories/) documentation or the book Tensor Categories [etingof2016tensor](@cite). + +Let us start simple and introduce the **fusion ring** $\mathcal{C}$ in a black-box manner. This ring +* consists of finitely many simple objects $\{ X_1, X_2, ..., X_R \}$, with $R$ the rank of the fusion ring, +* which can be fused with one another: $X_i \otimes X_j \cong \sum_k N_{ij}^k X_k,$ introducing the **N-symbol** $N_{ij}^k \in \mathbb{N}$ in the fusion rules, +* contains a *unique* unit object $1_\mathcal{C}$ which satisfies $1_\mathcal{C} \otimes X_i \cong X_i \otimes 1_\mathcal{C} \cong X_i$ for all objects $X_i \in \mathcal{C}$, +* has a dual object $\overline{X}$ for every object $X$ such that $X \otimes \overline{X} \cong \overline{X} \otimes X \cong 1_\mathcal{C} \oplus ...$, generalising the notion of the inverse element. + +To extend the fusion ring to the **fusion category**, we need to add the following structure: +* Consider only the representatives of isomorphism classes of simple objects $\mathcal{I}_\mathcal{C}$, +* The associator $F^{X_iX_jX_k}: (X_i \otimes X_j) \otimes X_k \xrightarrow{\sim} X_i \otimes (X_j \otimes X_k)$ + which fulfills the famous pentagon equation, +* Morphisms between (simple) objects $\text{Hom}_\mathcal{C}(X_i, X_j)$, which are empty vector spaces unless the objects are isomorphic, the latter then giving $\mathbb{C}$, +* More general morphisms $\text{Hom}_\mathcal{C}(X_i \otimes X_j, X_k) \cong \mathbb{C}^{N_{ij}^k}$. + +This way, we can describe fusion categories by a triple $(\otimes, 1_\mathcal{C}, F)$ of $\mathcal{C}$ defining its monoidal product, unit object and monoidal associator, the latter also commonly called the **F-symbol**. In particular, the simple objects have their respective quantum dimensions $d_i = \dim(X_i)$ which form their own one-dimensional representation of the fusion algebra: $d_i d_j = \sum_k N_{ij}^k d_k$. In particular, the unit object always has quantum dimension 1, and all other quantum dimensions are larger or equal to one. These quantum dimensions are encoded in the F-symbol. The isomorphisms instead of the equalities are a technical detail, so we drop that notation. + +Vectors in these hom-spaces are graphically denoted as living in the trivalent junction +```@raw html + +``` + + +With the F-symbol, we can perform F-moves: + +```@raw html + +``` + +TensorKit requires the F-symbols to be unitary. This way, we can interpret the F-symbol $F^{ijk}_l$ as a unitary matrix, and the F-move as a unitary basis transformation. Unitarity is also useful from a diagrammatic point of view because the category is then equipped with a pivotal and spherical structure. This essentially means that morphisms can be drawn and moved around freely on a 2-sphere, such that vector spaces can be moved freely from domain (codomain) to codomain (domain). + +## Examples + +### $\mathsf{Vec_G}$ and $\mathsf{Rep(G)}$ +Colloquially speaking, category theory attempts to generalise mathematical structures and their relations in a way that different structures can be treated in an equal manner. This is noted in particular as fusion category theory encompasses not only finite and compact groups, but also their representations. We show a table sketching how these are put on equal footing categorically. + +|$\mathsf{Vec_G}$|$\mathsf{Rep(G)}$|Categorical generalisation| +|:---:|:---:|:---:| +|$G$-graded vector spaces $V_1, V_2, ...$|Representations of $G$ $(V_1, \pi_1), (V_2, \pi_2), ...$|Objects| +|$G$-graded preserving linear maps $\phi: V \rightarrow W$|Intertwiners $f: V_1 \rightarrow V_2$, $f \circ \pi_1 = \pi_2 \circ f$|Morphisms $\text{Hom}_\mathcal{C}$| +1d subspaces $\mathbb{C}_{g_1}, \mathbb{C}_{g_2}$: $\text{Hom}_{\mathsf{Vec_G}}(\mathbb{C}_{g_1},\mathbb{C}_{g_2}) = \delta_{g_1g_2}$|Irreps: $\text{Hom}_{\mathsf{Rep(G)}}(\rho_i,\rho_j) = \delta_{ij} \mathbb{C}$ (Schur)|Simple objects: $\text{Hom}_{\mathcal{C}}(a,b) = \delta_{ab}\mathbb{C}$| +$G$-graded tensor product $(V \otimes W)_g = \oplus_{hk=g} V_h \otimes W_k$| $\pi_i \otimes \pi_j \simeq \oplus_i N_{ij}^k\rho_k$ | Direct sum, monoidal product, fusion rules, multiplicity| +$\mathbb{C}_1 \otimes W \simeq W \simeq W \otimes \mathbb{C}_1$ | Trivial rep 1: $1 \otimes \rho = \rho = \rho \otimes 1$ | Monoidal unit $1_\mathcal{C}$ +$\mathbb{C}_g \otimes \mathbb{C}_{g^{-1}} = \mathbb{C}_1 = \mathbb{C}_{g^{-1}} \otimes \mathbb{C}_g$ | $\pi \otimes \overline{\pi} = 1 \oplus ...$ | Dual object +$F:(V \otimes W) \otimes U \xrightarrow{\sim}V \otimes (W \otimes U)$|$F: (\pi_1 \otimes \pi_2) \otimes \pi_3 \xrightarrow{\sim} \pi_1 \otimes (\pi_2 \otimes \pi_3)$|F-symbol| + +### $\mathsf{Fib}$ and $\mathsf{Ising}$ +Arguably the simplest fusion category besides the familiar groups or representations of groups is the Fibonacci fusion category. This contain 2 simple objects $1$ and $\tau$, with non-trivial fusion rule $\tau \otimes \tau = 1 \oplus \tau$. This fusion category is in fact **braided** as well, and actually **modular**. + +Another simple fusion category is the Ising category, commonly denoted $\mathsf{Ising}$. 3 simple objects form this category, namely $\{1, \psi, \sigma\}$, where $1$ and $\psi$ behave like the trivial charged representations of $\mathbb{Z}_2$, while $\sigma$ is the $\mathbb{Z}_2$ extension of this. The fusion rules reflect these: $1 \otimes \psi = \psi = \psi \otimes 1, \sigma \otimes X = \sigma = X \otimes \sigma$ for $X = 1, \psi$, and $\sigma \otimes \sigma = 1 \oplus \psi$. This fusion category is also **modular**. Both these modular fusion categories are already implemented in [TensorKitSectors](https://github.com/QuantumKitHub/TensorKitSectors.jl). \ No newline at end of file diff --git a/docs/src/man/img/A4_sym_entanglement_spectrum.svg b/docs/src/man/img/A4_sym_entanglement_spectrum.svg new file mode 100644 index 0000000..c7fcead --- /dev/null +++ b/docs/src/man/img/A4_sym_entanglement_spectrum.svg @@ -0,0 +1,559 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/docs/src/man/img/Bmove_MF.svg b/docs/src/man/img/Bmove_MF.svg new file mode 100644 index 0000000..27513af --- /dev/null +++ b/docs/src/man/img/Bmove_MF.svg @@ -0,0 +1,963 @@ + + + + diff --git a/docs/src/man/img/Bmove_lit.svg b/docs/src/man/img/Bmove_lit.svg new file mode 100644 index 0000000..17667aa --- /dev/null +++ b/docs/src/man/img/Bmove_lit.svg @@ -0,0 +1,724 @@ + + + + diff --git a/docs/src/man/img/Fmove.svg b/docs/src/man/img/Fmove.svg new file mode 100644 index 0000000..0ed24e1 --- /dev/null +++ b/docs/src/man/img/Fmove.svg @@ -0,0 +1,1071 @@ + + + + diff --git a/docs/src/man/img/Fmove_CM.svg b/docs/src/man/img/Fmove_CM.svg new file mode 100644 index 0000000..b2faaa9 --- /dev/null +++ b/docs/src/man/img/Fmove_CM.svg @@ -0,0 +1,654 @@ + + + + diff --git a/docs/src/man/img/Fmove_CMD.svg b/docs/src/man/img/Fmove_CMD.svg new file mode 100644 index 0000000..27792b6 --- /dev/null +++ b/docs/src/man/img/Fmove_CMD.svg @@ -0,0 +1,678 @@ + + + + diff --git a/docs/src/man/img/Fmove_D.svg b/docs/src/man/img/Fmove_D.svg new file mode 100644 index 0000000..1838d4c --- /dev/null +++ b/docs/src/man/img/Fmove_D.svg @@ -0,0 +1,651 @@ + + + + diff --git a/docs/src/man/img/Fmove_MD.svg b/docs/src/man/img/Fmove_MD.svg new file mode 100644 index 0000000..e360755 --- /dev/null +++ b/docs/src/man/img/Fmove_MD.svg @@ -0,0 +1,668 @@ + + + + diff --git a/docs/src/man/img/Fmove_coloring.svg b/docs/src/man/img/Fmove_coloring.svg new file mode 100644 index 0000000..804e15d --- /dev/null +++ b/docs/src/man/img/Fmove_coloring.svg @@ -0,0 +1,2336 @@ + + + + diff --git a/docs/src/man/img/Nsymbol_coloring.svg b/docs/src/man/img/Nsymbol_coloring.svg new file mode 100644 index 0000000..1fcf989 --- /dev/null +++ b/docs/src/man/img/Nsymbol_coloring.svg @@ -0,0 +1,1208 @@ + + + + diff --git a/docs/src/man/img/homvector.svg b/docs/src/man/img/homvector.svg new file mode 100644 index 0000000..d019e1a --- /dev/null +++ b/docs/src/man/img/homvector.svg @@ -0,0 +1,342 @@ + + + + diff --git a/docs/src/man/img/pentagon_colored.svg b/docs/src/man/img/pentagon_colored.svg new file mode 100644 index 0000000..15b59e1 --- /dev/null +++ b/docs/src/man/img/pentagon_colored.svg @@ -0,0 +1,3478 @@ + + + + diff --git a/docs/src/man/img/qdim_fs_MF.svg b/docs/src/man/img/qdim_fs_MF.svg new file mode 100644 index 0000000..9ec6610 --- /dev/null +++ b/docs/src/man/img/qdim_fs_MF.svg @@ -0,0 +1,382 @@ + + + + diff --git a/docs/src/man/implementation.md b/docs/src/man/implementation.md new file mode 100644 index 0000000..28abc23 --- /dev/null +++ b/docs/src/man/implementation.md @@ -0,0 +1,166 @@ +# MultiTensorKit implementation: $\mathsf{Rep(A_4)}$ as a guiding example +This tutorial is dedicated to explaining how MultiTensorKit was implemented to be compatible with with TensorKit and MPSKit for matrix product state simulations. In particular, we will be making a generalised anyonic spin chain. We will demonstrate how to reproduce the entanglement spectra found in [Lootens_2024](@cite). The model considered there is a spin-1 Heisenberg model with additional terms to break the usual $\mathsf{U_1}$ symmetry to $\mathsf{Rep(A_4)}$, while having a non-trivial phase diagram and relatively easy Hamiltonian to write down. + +This will be done with the `A4Object = BimoduleSector{A4}` `Sector`, which is the multifusion category which contains the structure of the module categories over $\mathsf{Rep(A_4)}$. Since there are 7 module categories, `A4Object` is a $r=7$ multifusion category. There are 3 fusion categories up to equivalence: +- ``\mathsf{Vec_{A_4}}``: the category of $\mathsf{A_4}$-graded vector spaces. The group $\mathsf{A}_4$ is order $4!/2 = 12$. It has thus 12 objects. +- ``\mathsf{Rep(A_4)}``: the irreducible representations of the group $\mathsf{A}_4$, of which there are 4. One is the trivial representation, two are one-dimensional non-trivial and the last is three-dimensional. +- ``\mathsf{Rep(H)}``: the representation category of some Hopf algebra which does not have a name. It has 6 simple objects. + +For this example, we will require the following packages: +````julia +using TensorKit, MultiTensorKit, MPSKit, MPSKitModels, Plots +```` + +## Identifying the simple objects +We first need to select which fusion category we wish to use to grade the physical Hilbert space, and which fusion category to represent e.g. the symmetry category. In our case, we are interested in selecting $\mathcal{D} = \mathsf{Rep(A_4)}$ for the physical Hilbert space. We know the module categories over $\mathsf{Rep(G)}$ to be $\mathsf{Rep^\psi(H)}$ for a subgroup $\mathsf{H}$ and 2-cocycle $\psi$. Thus, the 7 module categories $\mathcal{M}$ one can choose over $\mathsf{Rep(A_4)}$ are + +- ``\mathsf{Rep(A_4)}`` itself as the regular module category, +- ``\mathsf{Vec}``: the category of vector spaces, +- ``\mathsf{Rep(\mathbb{Z}_2)}``, +- ``\mathsf{Rep(\mathbb{Z}_3)}``, +- ``\mathsf{Rep(\mathbb{Z}_2 \times \mathbb{Z}_2)}``, +- ``\mathsf{Rep^\psi(\mathbb{Z}_2 \times \mathbb{Z}_2)}``, +- ``\mathsf{Rep^\psi(A_4)}``. + +When referring to specific fusion and module categories, we will use this non-multifusion notation. + +The easiest way to identify which elements of the multifusion category correspond to the subcategories we wish to use is ... #TODO + +Now that we have identified the fusion and module categories, we want to select the relevant objects we wish to place in our graded spaces. Unfortunately, due to the nature of how the N-symbol and F-symbol data are generated, the objects of the fusion subcategories are not ordered such that `label=1` corresponds to the unit object. Hence, the simplest way to find the unit object of a fusion subcategory is + +````julia +one(A4Object(i,i,1)) +```` + +Left and right units of subcategories are uniquely specified by their fusion rules. For example, the left unit of a subcategory $\mathcal{C}_{ij}$ is the simple object in $\mathcal{C}_i$ for which + +$$^{}_a \mathbb{1} \times a = a \quad \forall a \in \mathcal{C}_i.$$ + +A similar condition uniquely defines the right unit of a subcategory. For fusion subcategories, a necessary condition is that the left and right units coincide. + +Identifying the other simple objects of a (not necessarily fusion) category requires more work. We recommend a combination of the following to uniquely determine any simple object `a`: +- Check the dimension of the simple object: `dim(a)` +- Check the dual of the simple object: `dual(a)` +- Check how this simple object fuses with other (simple) objects: `Nsymbol(a,b,c)` + +The dual object of some simple object $a$ of an arbitrary subcategory $\mathcal{C}_{ij}$ is defined as the unique object $a^* \in \mathcal{C}_{ji}$ satisfying + +$$^{}_a \mathbb{1} \in a \times a^* \quad \text{and} \quad \mathbb{1}_a \in a^* \times a$$ + +with multiplicity 1. +## Constructing the Hamiltonian and matrix product state +TensorKit has been made compatible with the multifusion structure by keeping track of the relevant units in the fusion tree manipulations. With this, we can make `GradedSpace`s whose objects are in `A4Object`: + +````julia +D1 = A4Object(6, 6, 1) # unit in this case +D2 = A4Object(6, 6, 2) # non-trivial 1d irrep +D3 = A4Object(6, 6, 3) # non-trivial 1d irrep +D4 = A4Object(6, 6, 4) # 3d irrep +```` +Since we want to replicate a spin-1 Heisenberg model, it makes sense to use the 3-dimensional irrep to grade the physical space, and thus construct our Hamiltonian. We don't illustrate here how to derive the considered Hamiltonian in a $\mathsf{Rep(A_4)}$ basis, but simply give it. For now, we construct a finite-size spin chain; below we will repeat the calculation for an infinite system. + +````julia +P = Vect[A4Object](D4 => 1) # physical space +T = ComplexF64 +# usual Heisenberg part +h1_L = TensorMap(zeros, T, P βŠ— P ← P) +h1_R = TensorMap(zeros, T, P ← P βŠ— P) +block(h1_L, D4) .= [0; 1] +block(h1_R, D4) .= [0 1;] +@plansor h1[-1 -2; -3 -4] := h1_L[-1 1; -3] * h1_R[-2; 1 -4] + +# biquadratic term +h2_L = TensorMap(zeros, T, P βŠ— Vect[A4Object](D1 => 1, D2 => 1, D3 => 1) ← P) +h2_R = TensorMap(ones, T, P ← Vect[A4Object](D1 => 1, D2 => 1, D3 => 1) βŠ— P) +block(h2_L, D4) .= [4 / 3; 1 / 3; 1 / 3] +@plansor h2[-1 -2; -3 -4] := h2_L[-1 1; -3] * h2_R[-2; 1 -4] + +# anti-commutation term +h3_L = TensorMap(zeros, T, P βŠ— P ← P) +h3_R = TensorMap(zeros, T, P ← P βŠ— P) +block(h3_L, D4) .= [1; 0] +block(h3_R, D4) .= [0 1;] +@plansor h3[-1 -2; -3 -4] := h3_L[-1 1; -3] * h3_R[-2; 1 -4] + +L = 60 +J1 = -2.0 # probing the A4 symmetric phase first +J2 = -5.0 +lattice = FiniteChain(L) +H1 = @mpoham sum(-2 * h1{i,j} for (i, j) in nearest_neighbours(lattice)) +H2 = @mpoham sum(h2{i,j} for (i, j) in nearest_neighbours(lattice)) +H3 = @mpoham sum(2im * h3{i,j} for (i, j) in nearest_neighbours(lattice)) + +H = H1 + J1 * H2 + J3 * H3 +```` + +For the matrix product state, we will select $\mathsf{Vec}$ as the module category for now: +````julia +M = A4Object(1, 6, 1) # Vec +```` +and construct the finite MPS: +````julia +D = 40 # bond dimension +V = Vect[A4Object](M => D) +Vb = Vect[A4Object](M => 1) # non-degenerate boundary virtual space +init_mps = FiniteMPS(L, P, V; left=Vb, right=Vb) +```` + +!!! warning "Important" + We must pass on a left and right virtual space to the keyword arguments `left` and `right` of the `FiniteMPS` constructor, since these would by default try to place a trivial space of the `Sector`, which does not exist for any `BimoduleSector` due to the semisimple unit. + +## DMRG2 and the entanglement spectrum +We can now look to find the ground state of the Hamiltonian with two-site DMRG. We use this instead of the "usual" one-site DMRG because the two-site one will smartly fill up the blocks of the local tensor during the sweep, allowing one to initialise as a product state in one block and more likely avoid local minima, a common occurence in symmetric tensor network simulations. +````julia +dmrg2alg = DMRG2(;verbosity=2, tol=1e-7, trscheme=truncbelow(1e-4)) +ψ, _ = find_groundstate(init_mps, H, dmrg2alg) +```` +The truncation scheme keyword argument is mandatory when calling `DMRG2` in MPSKit. Here, we choose to truncate such that all singular values are larger than $10^{-4}$, while setting the default tolerance for convergence to $10^{-7}$. More information on this can be found in the [MPSKit](https://github.com/QuantumKitHub/MPSKit.jl) documentation. To run one-site DMRG anyway, use `DMRG` which does not require a truncation scheme. + +Now that we've found the ground state, we can compute the entanglement spectrum in the middle of the chain. +````julia +spec = entanglement_spectrum(ψ, round(Int, L/2)) +```` +This returns a dictionary which maps the objects grading the virtual space to the singular values. In this case, there is one key corresponding to $\mathsf{Vec}$. We can also immediately return a plot of this data by the following: +````julia +entanglementplot(ψ;site=round(Int, L/2)) +```` +This plot will show the singular values per object, as well as include the "effective" bond dimension, which is simply the dimension of the virtual space where we cut the system. The next section will show this plot along with those when selecting the other module categories. + +## Search for the correct dual model + +Consider a quantum lattice model with its symmetries determing the phase diagram. For every phase in the phase diagram, the dual model for which the ground state maximally breaks all symmetries spontaneously is the one where the entanglement is minimised and the tensor network is represented most efficiently [Lootens_2024](@cite). Let us confirm this result, starting with the $\mathsf{Rep(A_4)}$ spontaneous symmetry breaking phase. The code will look exactly the same as above, except the virtual space of the MPS will change to be graded by the other module categories: + +````julia +module_numlabels(i) = MultiTensorKit._numlabels(A4Object, i, 6) +V = Vect[A4Object](A4Object(i, 6, label) => D for label in 1:module_numlabels(i)) +Vb = Vect[A4Object](first(sectors(V)) => 1) # not all charges on boundary, play around with what is there +```` + +```@raw html + +``` + +The plot shows the entanglement spectra of the various dual models in the middle of the ground state for the original Hamiltonian probing the $\mathsf{Rep(A_4)}$-symmetric phase, along with the ground state without symmetries. Every subtitle also mentions the memory required to store the same middle ground state tensor. This plot should be compared to [Lootens_2024; Figure 2](@cite). + +This can be repeated with other parameter values for $J_1$ and $J_2$ in the Hamiltonian to probe the $\mathsf{Rep(\mathbb{Z}_2 \times \mathbb{Z}_2)}$-symmetric or $\mathsf{Rep^\psi(A_4)}$ SPT phase. + +!!! note "Additional functions and keyword arguments" + Certain commonly used functions within MPSKit require extra keyword arguments to be compatible with multifusion MPS simulations. In particular, the keyword argument `sector` (note the lowercase "s") appears in + - `excitations` with `QuasiparticleAnsatz`: the sector is selected by adding an auxiliary space to the *domain* of each eigenvector of the transfer matrix. Since in a full contraction the domain of the eigenvector lies in the opposite side of the physical space (labeled by objects in $\mathcal{D} = \mathsf{Rep(A_4)}$), the charged excitations lie in the symmetry category $\mathcal{C} = \mathcal{D^*_M}$. + - `exact_diagonalization`: the `sector` keyword argument now requires an object in $\mathcal{D}$, since this is the fusion category which specifies the bond algebra from which the Hamiltonian is constructed. This is equivalent to adding a charged leg on the leftmost (or rightmost) virtual space of the MPS in conventional MPS cases. + +## Differences with the infinite case +We can repeat the above calcalations also for an infinite system. The `lattice` variable will change, as well as the MPS constructor and the algorithm: +````julia +lattice = InfiniteChain(1) +init = InfiniteMPS([P], [V]) +inf_alg = VUMPS(; verbosity=2, tol=1e-7) +```` + +Besides `VUMPS`, `IDMRG` and `IDMRG2` are as easy to run with the `A4Object` `Sector`. It is also clear that boundary terms do not play a role in this case. + +!!! note "More functions for infinite systems" + When dealing with an infinite system, additional information can be retrieved from the matrix product state. These also require a keyword argument `sector` to be specified. These are + - `transfer_spectrum`: similar to `excitations`, the (partial) transfer matrix spectrum is selected by adding a charged auxiliary space to the transfer matrix eigenvectors. + - `correlation_length`: since this function calls `transfer_spectrum`, the same logic applies. + - `excitations` in the infinite case also requires the keyword argument. \ No newline at end of file diff --git a/docs/src/man/multifusioncats.md b/docs/src/man/multifusioncats.md new file mode 100644 index 0000000..414eea1 --- /dev/null +++ b/docs/src/man/multifusioncats.md @@ -0,0 +1,201 @@ +# Extending to multifusion category theory + +This section will explain how to go from a fusion category to a multifusion category, as well as why one would want to consider the latter. Multifusion categories naturally embed the structure of **bimodule categories**. To explain this, we must start by explaining module categories over fusion categories, following this up with (invertible) bimodule categories, and finishing off with the multifusion structure. + +## Module categories +We start from a fusion category $\mathcal{D}$ defined by the triple $(\otimes_\mathcal{D}, \mathbb{1}_\mathcal{D}, {}^\mathcal{D}\!F)$ with simple objects $\alpha, \beta, ... \in \mathcal{I}_\mathcal{D}$. We drop the $\mathcal{D}$ subscript when there is no ambiguity concerning the fusion category. We call its associator + +$${}^\mathcal{D}\!F^{\alpha \beta \gamma}: \alpha \otimes (\beta \otimes \gamma) \rightarrow (\alpha \otimes \beta) \otimes \gamma$$ + + the **monoidal associator**. An F-move is now graphically portrayed as: + +```@raw html + +``` + +We can consider the **right module category** $\mathcal{M}$ over $\mathcal{D}$, which is a category (not necessarily fusion!) with (isomorphism classes of) simple objects $\mathcal{I}_\mathcal{M} = \{A,B,...\}$, a right action + +$$\triangleleft: \mathcal{M} \times \mathcal{D} \rightarrow \mathcal{M}$$ + +and the **right module associator** + +$${}^\triangleleft\!F^{A\alpha\beta}: A \triangleleft (\alpha \otimes \beta) \rightarrow (A \triangleleft \alpha) \triangleleft \beta.$$ + +An F-move with this module associator can be expressed as: + +```@raw html + +``` + +The module structure of $\mathcal{M}$ is now defined as the triple $(\mathcal{M}, \triangleleft, {}^\triangleleft\!F)$. The right module associator now satisfies a mixed pentagon equation with ${}^\mathcal{D}\!F$. + +Similarly, we can define a **left module category** $(\mathcal{M}, \triangleright, {}^\triangleright\!F)$ over a fusion category $(\otimes_\mathcal{C}, 1_\mathcal{C}, {}^\mathcal{C}\!F)$. The functor is now a left action of $\mathcal{C}$ on $\mathcal{M}$ + +$$\triangleright: \mathcal{C} \times \mathcal{M} \rightarrow \mathcal{M},$$ + +while the **left module associator** is a natural isomorphism that acts as + +$${}^\triangleright\!F^{abA}: (a \otimes b) \triangleright A \rightarrow a \triangleright (b \triangleright A)$$ + +for $\mathcal{I}_\mathcal{C} = \{a,b,...\}$. The left module associator also fulfills a mixed pentagon equation with ${}^\mathcal{C}\!F$. An F-move with ${}^\mathcal{C}\!F$ takes on the form: + +```@raw html + +``` + +We can combine the concepts of left and right module categories as follows. Say there are two fusion categories $\mathcal{C}$ and $\mathcal{D}$. A $(\mathcal{C}, \mathcal{D})$-**bimodule category** is a module category, now defined through a sextuple $(\mathcal{M}, \triangleright, \triangleleft, {}^\triangleright\!F, {}^\triangleleft\!F, {}^{{\triangleright \hspace{-1.2mu}\triangleleft}}\!F)$ such that $(\mathcal{M}, \triangleright, {}^\triangleright\!F)$ is a left $\mathcal{C}$-module category, $(\mathcal{M}, \triangleleft, {}^\triangleleft\!F)$ is a right $\mathcal{D}$-module category, and with additional structure such that the **bimodule associator** acts as + +$${}^{{\triangleright \hspace{-1.2mu}\triangleleft}}\!F^{aA\alpha}: (a \triangleright A) \triangleleft \alpha \rightarrow a \triangleright (A \triangleleft \alpha)$$ + +for $a \in \mathcal{I}_\mathcal{C}, \alpha \in \mathcal{I}_\mathcal{D}, A \in \mathcal{I}_\mathcal{M}$. The bimodule associator fulfills a mixed pentagon equation with the module associators. An F-move with ${}^{{\triangleright \hspace{-1.2mu}\triangleleft}}\!F$ is given by: + +```@raw html + +``` +## Opposite module categories +Consider a fusion category $\mathcal{D}$ and a *right* module category $\mathcal{M}$ over $\mathcal{D}$. We can define $\mathcal{M}^{\text{op}}$ to be the opposite category of $\mathcal{M}$ [etingof2009](@cite). Then, $\mathcal{M}^{\text{op}}$ is a *left* module category over $\mathcal{D}$. A similar statement can be made starting from a left module category and getting an opposite right module category. In particular, given a $(\mathcal{C}, \mathcal{D})$-bimodule category $\mathcal{M}$ over the fusion categories $\mathcal{C}$ and $\mathcal{D}$, we see immediately that $\mathcal{M}^{\text{op}}$ is a $(\mathcal{D}, \mathcal{C})$-bimodule category. + +Interestingly, due to the opposite actions + +$$\triangleleft^{\text{op}}: \mathcal{D} \times \mathcal{M}^{\text{op}} \rightarrow \mathcal{M}^{\text{op}}$$ + +and + +$$\triangleright^{\text{op}}: \mathcal{M}^{\text{op}} \times \mathcal{C} \rightarrow \mathcal{M}^{\text{op}},$$ + +there is a valid notion of multiplying a module category with its opposite: + +$$\mathcal{M} \times \mathcal{M}^{\text{op}} \rightarrow \mathcal{C}, \quad \mathcal{M}^{\text{op}} \times \mathcal{M} \rightarrow \mathcal{D}.$$ + +## Invertible bimodule categories and Morita equivalence + +The bimodule categories we consider are restricted to be **invertible**. This means that this bimodule category fulfills the condition $\mathcal{C} \equiv \mathcal{D}^*_\mathcal{M}$, or in other words the two fusion categories $\mathcal{C}$ and $\mathcal{D}$ are **Morita equivalent** (or each other's **Morita dual** with respect to $\mathcal{M}$). Morita equivalence between $\mathcal{C}$ and $\mathcal{D}$ can be shown to hold if and only if the bimodule category fulfills + +$$\mathcal{M} \boxtimes_\mathcal{D} \mathcal{M}^\text{op} \simeq \mathcal{C}, \quad \text{and} \quad \mathcal{M}^\text{op} \boxtimes_\mathcal{C} \mathcal{M} \simeq \mathcal{D},$$ + +where $\boxtimes_\mathcal{C}$ and $\boxtimes_\mathcal{D}$ denote the Deligne product relative to $\mathcal{C}$ and $\mathcal{D}$ respectively. This is precisely the invertibility property of the bimodule category. + +## Multifusion categories +A fusion category has the important condition that $\mathsf{End}_\mathcal{C}(1_\mathcal{C}) \cong \mathbb{C}$, i.e., the unit of the fusion category is simple. If we drop this condition, then we consider a **multifusion category**. We assume the multifusion category itself to be indecomposable, meaning it is not the direct sum of two non-trivial multifusion categories. Let us call this multifusion category $\mathcal{C}$. It will be clear in a moment that this will not be ambiguous. Its unit can then be decomposed as + +$$1_\mathcal{C} = \bigoplus_{i=1}^r 1_r,$$ + +i.e., it is decomposed into simple unit objects of the subcategories $\mathcal{C}_{ij} \coloneqq 1_i \otimes \mathcal{C} \otimes 1_j$. With this, we see that the multifusion category itself can be decomposed into its subcategories + +$$\mathcal{C} = \bigoplus_{i,j=1}^r \mathcal{C}_{ij}.$$ + +We call this an $r \times r$ multifusion category. Due to this structure, we represent a simple object in the multifusion category "Name" by + +````julia +struct BimoduleSector{Name} <: Sector + i::Int + j::Int + label::Int +end +```` + +`i` and `j` specify which subcategory $\mathcal{C}_{ij}$ we are considering, and `label` selects a particular simple object within that subcategory. + +We want to consider multifusion categories because **their structure encapsulates that of (bi)module categories**. Every diagonal category $\mathcal{C}_{ii} \coloneqq \mathcal{C}_i$ (also known as a component category) is a fusion category, and every off-diagonal category $\mathcal{C}_{ij}$ is an invertible $(\mathcal{C}_{i}, \mathcal{C}_{j})$-bimodule category. That way, as long as we know how the simple objects of the fusion and module categories fuse with one another, and we can determine all the monoidal and module associators, we can treat the multifusion category as one large fusion category with limited fusion rules. In particular, the tensor product + +$$\otimes_\mathcal{C}: \mathcal{C}_{ij} \times \mathcal{C}_{kl} \rightarrow \delta_{jk}\mathcal{C}_{il}$$ + +takes on the same structure as the product of two matrix units. This is not a coincidence; there is a deep relation between multifusion categories and matrix algebras [etingof2016tensor; Section 4.3](@cite). + +Given a subcategory $\mathcal{C}_{ij}$, we can define the **left** unit as the unit object of the fusion category $\mathcal{C}_i$, while the **right** unit is the unit object of the fusion category $\mathcal{C}_j$. In other words, the left unit of $\mathcal{C}_{ij}$ is the unique object of the multifusion category $\mathcal{C}$ for which + +$$1_i \otimes_\mathcal{C} m_{ij} = m_{ij} \quad \forall m_{ij} \in \mathcal{C}_{ij},$$ + +and similarly for the right unit of $\mathcal{C}_{ij}$, + +$$m_{ij} \otimes_\mathcal{C} 1_j = m_{ij} \quad \forall m_{ij} \in \mathcal{C}_{ij}.$$ + +We can also immediately see that for a bimodule subcategory $\mathcal{C}_{ij}$, the opposite (bi)module subcategory $\mathcal{C}_{ij}^{\text{op}} \equiv \mathcal{C}_{ji}$, and as expected, + +$$\mathcal{C}_{ij} \times \mathcal{C}_{ji} \rightarrow \mathcal{C}_i,$$ + +just like what we concluded when considering opposite module categories outside of the multifusion structure. +### 2-category and coloring +Multifusion categories can also be interpreted as 2-categories. We still interpret the objects of this 2-category the same way. The 1-morphisms are the subcategories themselves, and the 2-morphisms the morphisms of the multifusion category. The graphical calculus of monoidal 1-categories can be extended to 2-categories by use of *colorings* [henriques2020](@cite). We have previously differed between module strands and fusion strands by the color of the strand itself. However, in 2-categories the strands (1-morphisms) separate regions which are colored based on the fusion categories they are representing. Since we draw the strands vertically, a single strand results in a left and right region, and the colorings will determine the fusion category which fuses from the left or right with that single strand. In particular, fusion strands necessarily have the same coloring on the left and right, while module strands have a mismatching coloring. + +The simplest non-trivial fusion diagram is a trivalent junction: + +```@raw html + +``` + +The most general case is the top left figure, where all three regions have a different coloring. The top middle region having the same coloring from the top left and top right strands follow from the delta function in the tensor product definition. However, this most general trivalent junction with three colorings will never be needed. In short, we will always be considering a single bimodule category $\mathcal{C}_{ij}$ at a time, and the only other non-diagonal subcategory which fuses with this is its opposite $\mathcal{C}_{ji}$. This is displayed in the top middle and right. Similarly, two colorings are required when considering the fusion between a fusion and module strand, shown in the bottom left and middle figure. The simplest trivalent junctions boil down to fusions within fusion categories, which is obviously drawn with just one color. This is shown in the bottom right. + +With this coloring system, we can specify which associator must be called to perform a particular F-move. Such an F-move would look like + +```@raw html + +``` + +### Why opposite module categories end up being necessary in MultiTensorKit + +One of the common manipulations which can act on a tensor map is the transposition of vector spaces. We will refer to this as the bending of legs. One of the elementary bends is the **right bend**, where one of the tensor legs is bent along the right from the codomain to the domain, or vice versa. At the level of the tensor, a covariant index becomes contravariant, or vice versa. Similarly, a **left bend** can also be performed, bending the leg along the left. This guarantees that legs will not cross, preventing braidings which require extra data known as R-symbols. + +Linear algebra tells us that given a (finite-dimensional) vector space $V$ with a basis denoted $\{|e_i\rangle\}$, one can consider the **dual** vector space $V^*$, whose dual basis $\{\langle e_i^*|\}$ satisfies the property $\langle e_j^* | e_i \rangle = \delta_{ij}$. In the diagrammatic calculus, specifying whether a tensor map leg represents a vector space or its dual is done with an arrow. Following the TensorKit convention, legs with arrows pointing downwards are vector spaces, and arrows pointing upwards state that we are considering its dual. In particular, at the level of fusion trees we can also draw arrows on the strands to denote whether we are considering morphisms between objects or dual objects. + +In principle, choosing to bend, e.g., codomain legs to the right and domain legs to left is an arbitrary choice, but would require to distinguish between left and right transposes. However, TensorKit.jl is implemented in a way that does not differentiate the two. In particular, we do not have to worry about this when considering categorical symmetries where, in principle, the left and right dual of an object are not equivalent. This is because this left-right symmetry is guaranteed when considering unitary fusion categories, which is what is done in TensorKit and necessarily in MultiTensorKit. + +For this reason, at the level of the fusion trees the topological move that is performed to bend the legs along the right is called a **B-move**. Graphically, one can show that this bend boils down to a particular F-move. The typical equation found in the literature is the following: +```@raw html + +``` + +The reason to only consider B-moves is rooted in the choice of canonical form of fusion trees within TensorKit, where fusions are iterated over from left to right and splittings from right to left. + +Importantly, we identify the dual vector space labeled by a module category with a vector space labeled by the opposite module category. Consequently, + +$$\mathcal{M}^* \simeq \mathcal{M}^\text{op}.$$ + +In the multifusion setting, this can also be seen graphically. By keeping track of the colorings and the directions of the arrows of the legs, one can see that we need to slightly modify the expression for the B-move to the following: + +```@raw html + +``` + +where by $\mathbb{1}_j$ we mean the unit of $\mathcal{C}_j$. + +### More on the topological data: gauge choices and distilling properties of the subcategories + +An important property of the F-symbols is that they must satisfy the **triangle identities**. In fusion category theory, this states that isomorphisms between (simple) objects $a$ and the tensor product between $a$ and the unit $\mathbb{1}$ exists, and that + +$$(a \otimes \mathbb{1}) \otimes b \cong a \otimes (\mathbb{1} \otimes b)$$ + +for $b$ in the same fusion category. This can be straightforwardly generalised to multifusion categories. This requires a particular partial gauging of these trivalent vertices. + +Besides the triangle identities, the (multi)fusion category must also fulfill the **pentagon equations**. These encapsulate the two identical manners to evaluating the fusion of four objects in the (multi)fusion category. Every fusion category's F-symbols must satisfy these individually, but also the (bi)module associators between bimodule and fusion categories. One can check that, for every pair of fusion categories, their bimodule category and opposite bimodule category, there are 32 pentagon equations to be satisfied. In the multifusion notation, they can be represented by + +```@raw html + +``` + +The most generic F-move contains 4 colors. For that reason, MultiTensorKit requires the F-symbol data to be provided as some data file (currently .txt) with 4 + 6 + 4 + 2 = 16 columns. The first 4 refer to the colors, the following 6 label the simple objects of the corresponding subcategories, the next 4 label multiplicities, and the final 2 provide the real and imaginary value of the F-symbol itself. + +In a similar manner, the N-symbols contain maximally three colors, so these data must provide 3 columns labeling the colors, 3 columns labeling the simple objects and a final column with the dimension of the corresponding vector space. + +Besides the B-move (and closely related A-move, which we do not illustrate), we can also see how the quantum dimension and Frobenius-Schur indicator expressions get modified. We already know that an F-move of the form $F^{c \bar{c} c}_{c}$ needs to be evaluated for these topological data, which is in fact another gauge fixing. Graphically, we find + +```@raw html + +``` +In principle, a gauge fixing can be done to set the Frobenius-Schur indicator to $\pm 1$. However, this assumption is no longer required within TensorKit and can be relaxed to just be a phase. + +The above gauge fixing is a property of the unitary gauge. By choosing appropriate bases, one can transform the F-symbols of a (multi)fusion category to be unitary matrices. More details on the importance of unitary topological data can be found in the [TensorKit](https://jutho.github.io/TensorKit.jl/stable/man/categories/#ss_topologicalfusion) documentation. + +### Braiding +A very important aspect of MultiTensorKit is that all `BimoduleSector`s are defined to *not* support braiding: `TensorKitSectors.BraidingStyle(::Type{<:BimoduleSector}) = NoBraiding()`. We do this for two reasons. On the one hand, there is no natural 1-categorical way of defining braidings between the components of the multifusion category. It is possible that the diagonal fusion categories themselves are braided, but a "componentwise" braiding is unwise to support. On the other hand, it is entirely possible to write matrix product state manipulations in a planar manner (which has been done in [MPSKit](https://github.com/QuantumKitHub/MPSKit.jl)), thus avoiding the need of a braiding tensor. + +## Examples of multifusion categories +Without specifying any of the categories, the simplest non-trivial multifusion category is a $2\times 2$ one, and the categories can be organised in a matrix as + +$$\mathcal{C} = \begin{pmatrix} \mathcal{C}_1 & \mathcal{M} \\ \mathcal{M}^{\text{op}} & \mathcal{C}_2\end{pmatrix}.$$ + +We already identified the off-diagonal elements with module categories over the fusion categories on the diagonal. Accordingly, $\mathcal{M}$ is a $(\mathcal{C}_1, \mathcal{C}_2)$-bimodule category, and $\mathcal{M}^{\text{op}}$ is the opposite module category and a $(\mathcal{C}_2, \mathcal{C}_1)$-bimodule category. + +If we take $\mathcal{C}_1 = \mathcal{C}_2 = \mathsf{Rep}(\mathbb{Z}_2)$ and $\mathcal{M} = \mathsf{Vec}$, then the entire multifusion category is Morita equivalent to $\mathsf{Rep}(\mathbb{Z}_2)$, and we view the $\mathbb{Z}_2$-extension of $\mathsf{Rep}(\mathbb{Z}_2)$ to be precisely the $\mathsf{Ising}$ category [etingof2009; Section 9](@cite) [etingof2016tensor; Example 4.10.5](@cite). We identify the trivial representation of $\mathsf{Rep}(\mathbb{Z}_2)$ with the unit of $\mathsf{Ising}$, the sign representation with $\psi$ and the unique object of $\mathsf{Vec}$ with the duality object $\sigma$. One can easily check that the fusion rules of $\mathsf{Ising}$ match with those we expect within $\mathsf{Rep}(\mathbb{Z}_2)$ and with its module category $\mathsf{Vec}$. Additionally, the fusion between $\mathsf{Vec}$ and $\mathsf{Vec}^\text{op}$ (and vice-versa) giving every object in $\mathcal{C}_1$ ($\mathcal{C}_2$) is consistent with $\sigma \times \sigma^* = 1 + \psi$. This particular example can be found in [TensorKitSectors](https://github.com/QuantumKitHub/TensorKitSectors.jl). + +This construction can be generalised to $\mathcal{C}_1 = \mathcal{C}_2 = \mathsf{Rep(G)}$ with $\mathsf{G}$ a finite abelian group, such that the entire multifusion category is Morita equivalent to $\mathsf{Rep(G)}$ and can be evaluated as the Tambara-Yamagami category $\mathsf{TY}(\mathsf{G})$ (with positive Frobenius-Schur indicator for our purposes), and $\mathsf{Vec}$ will represent the duality object which squares to all invertible objects of the original group. To be exact, one of the diagonal fusion categories should be $\mathsf{Vec_G}$ for the correct Morita dual relation, but it is known for abelian groups that this is isomorphic to $\mathsf{Rep(G)}$. \ No newline at end of file diff --git a/docs/src/references.md b/docs/src/references.md new file mode 100644 index 0000000..4b47677 --- /dev/null +++ b/docs/src/references.md @@ -0,0 +1,4 @@ +# References + +```@bibliography +``` \ No newline at end of file diff --git a/src/MultiTensorKit.jl b/src/MultiTensorKit.jl index 34d66bc..c321bff 100644 --- a/src/MultiTensorKit.jl +++ b/src/MultiTensorKit.jl @@ -3,9 +3,19 @@ module MultiTensorKit export BimoduleSector, A4Object using JSON3 +using DelimitedFiles using Artifacts using TensorKitSectors +using TupleTools +import TupleTools: insertafter + +using BlockTensorKit +import BlockTensorKit: SumSpace + +using TensorKit +import TensorKit: hasblock, dim + include("bimodulesector.jl") end diff --git a/src/bimodulesector.jl b/src/bimodulesector.jl index d9473c3..1cf4673 100644 --- a/src/bimodulesector.jl +++ b/src/bimodulesector.jl @@ -1,7 +1,24 @@ +""" + struct BimoduleSector{Name} <: Sector + BimoduleSector{Name}(i::Int, j::Int, label::Int) + +Represents objects in the component subcategory ``π’žα΅’β±Ό`` of the multifusion category ``π’ž = ⨁ᡒⱼ π’žα΅’β±Ό``, +where ``π’ž`` is identified as Name. + +## Fields +- `i::Int`: The row index of the object in the matrix representation of the multifusion category. +- `j::Int`: The column index of the object in the matrix representation of the multifusion category. +- `label::Int`: The label of the object within the component subcategory. +""" struct BimoduleSector{Name} <: Sector i::Int j::Int label::Int + function BimoduleSector{:A4}(i::Int, j::Int, label::Int) + i <= 12 && j <= 12 || throw(DomainError("object outside the matrix A4")) + return label <= _numlabels(BimoduleSector{:A4}, i, j) ? new{:A4}(i, j, label) : + throw(DomainError("label outside category A4($i, $j)")) + end end BimoduleSector{Name}(data::NTuple{3,Int}) where {Name} = BimoduleSector{Name}(data...) const A4Object = BimoduleSector{:A4} @@ -19,11 +36,10 @@ end Base.IteratorSize(::Type{SectorValues{<:BimoduleSector}}) = Base.SizeUnknown() # TODO: generalize? -# TODO: numlabels? function Base.iterate(iter::SectorValues{A4Object}, (I, label)=(1, 1)) I > 12 * 12 && return nothing i, j = CartesianIndices((12, 12))[I].I - maxlabel = numlabels(A4Object, i, j) + maxlabel = _numlabels(A4Object, i, j) return if label > maxlabel iterate(iter, (I + 1, 1)) else @@ -31,8 +47,13 @@ function Base.iterate(iter::SectorValues{A4Object}, (I, label)=(1, 1)) end end +function Base.length(::SectorValues{A4Object}) + return sum(_numlabels(A4Object, i, j) for i in 1:12, j in 1:12) +end + TensorKitSectors.FusionStyle(::Type{A4Object}) = GenericFusion() TensorKitSectors.BraidingStyle(::Type{A4Object}) = NoBraiding() +TensorKitSectors.sectorscalartype(::Type{A4Object}) = ComplexF64 function TensorKitSectors.:βŠ—(a::A4Object, b::A4Object) @assert a.j == b.i @@ -42,8 +63,8 @@ function TensorKitSectors.:βŠ—(a::A4Object, b::A4Object) if (a_l == a.label && b_l == b.label)] end -function _numlabels(::Type{A4Object}, i, j) - return Ncache = _get_Ncache(A4Object) +function _numlabels(::Type{T}, i, j) where {T<:BimoduleSector} + return length(_get_dual_cache(T)[2][i, j]) end # Data from files @@ -70,6 +91,7 @@ const Ncache = IdDict{Type{<:BimoduleSector},Array{Dict{NTuple{3,Int},Int},3}}() function _get_Ncache(::Type{T}) where {T<:BimoduleSector} global Ncache return get!(Ncache, T) do + @debug "loading Nsymbol cache for $T" return extract_Nsymbol(T) end end @@ -82,12 +104,12 @@ function TensorKitSectors.Nsymbol(a::I, b::I, c::I) where {I<:A4Object} return get(_get_Ncache(I)[i, j, k], (a.label, b.label, c.label), 0) end -# TODO: can we define dual for modules? -const Dualcache = IdDict{Type{<:BimoduleSector},Vector{Tuple{Int,Vector{Int}}}}() +const Dualcache = IdDict{Type{<:BimoduleSector},Tuple{Vector{Int64},Matrix{Vector{Int64}}}}() function _get_dual_cache(::Type{T}) where {T<:BimoduleSector} global Dualcache return get!(Dualcache, T) do + @debug "loading dual cache for $T" return extract_dual(T) end end @@ -95,94 +117,224 @@ end function extract_dual(::Type{A4Object}) N = _get_Ncache(A4Object) ncats = size(N, 1) - return map(1:ncats) do i + Is = zeros(Int, ncats) + + map(1:ncats) do i Niii = N[i, i, i] - nobj = maximum(first, keys(Niii)) - - # find identity object: - # I x I -> I, a x I -> a, I x a -> a - I = findfirst(1:nobj) do u - get(Niii, (u, u, u), 0) == 1 || return false - for j in 1:nobj - get(Niii, (j, u, j), 0) == 1 || return false - get(Niii, (u, j, j), 0) == 1 || return false + nobji = maximum(first, keys(N[i, i, i])) + # want to return a leftone and rightone for each entry in multifusion cat + # leftone/rightone needs to at least be the unit object within a fusion cat + Is[i] = findfirst(1:nobji) do a + get(Niii, (a, a, a), 0) == 1 || return false # I x I -> I + for othera in 1:nobji + get(Niii, (othera, a, othera), 0) == 1 || return false # a x I -> a + get(Niii, (a, othera, othera), 0) == 1 || return false # I x a -> a + end + + # check leftone + map(1:ncats) do j + nobjj = maximum(first, keys(N[j, j, j])) + for b in 1:nobjj + get(N[i, j, j], (a, b, b), 0) == 1 || return false # I = leftone(b) + end + end + + # check rightone + map(1:ncats) do k + nobjk = maximum(first, keys(N[k, k, k])) + for c in 1:nobjk + get(N[k, i, k], (c, a, c), 0) == 1 || return false # I = rightone(c) + end end return true end + end + + allduals = Matrix{Vector{Int}}(undef, ncats, ncats) # ncats square matrix of vectors + for i in 1:ncats + nobji = maximum(first, keys(N[i, i, i])) + for j in 1:ncats + allduals[i, j] = Int[] - # find duals - # a x abar -> I - duals = map(1:nobj) do j - return findfirst(1:nobj) do k - return get(Niii, (j, k, I), 0) == 1 + nobjj = maximum(first, keys(N[j, j, j])) + # the nested vectors contain the duals of the objects in π’ž_ij, which are in C_ji + Niji = N[i, j, i] # π’ž_ij x π’ž_ji -> C_ii + Njij = N[j, i, j] # π’ž_ji x π’ž_ij -> C_jj + for i_ob in 1:nobji, j_ob in 1:nobjj + get(Niji, (i_ob, j_ob, Is[i]), 0) == 1 || continue # leftone(c_ij) ∈ c_ij x c_ji + get(Njij, (j_ob, i_ob, Is[j]), 0) == 1 || continue # rightone(c_ij) ∈ c_ji x c_ij + push!(allduals[i, j], j_ob) end end - - return I, duals end + return Is, allduals end function Base.one(a::BimoduleSector) a.i == a.j || error("don't know how to define one for modules") - return A4Object(a.i, a.i, _get_dual_cache(typeof(a))[a.i][1]) + return A4Object(a.i, a.i, _get_dual_cache(typeof(a))[1][a.i]) +end + +function TensorKitSectors.leftone(a::BimoduleSector) + return A4Object(a.i, a.i, _get_dual_cache(typeof(a))[1][a.i]) +end + +function TensorKitSectors.rightone(a::BimoduleSector) + return A4Object(a.j, a.j, _get_dual_cache(typeof(a))[1][a.j]) end function Base.conj(a::BimoduleSector) - a.i == a.j || error("don't know how to define dual for modules") - return A4Object(a.i, a.i, _get_dual_cache(typeof(a))[a.i][2][a.label]) + return A4Object(a.j, a.i, _get_dual_cache(typeof(a))[2][a.i, a.j][a.label]) end function extract_Fsymbol(::Type{A4Object}) - return mapreduce((x, y) -> cat(x, y; dims=1), 1:12) do i - filename = joinpath(artifact_path, "A4", "Fsymbol_$i.json") + result = Dict{NTuple{4,Int},Dict{NTuple{6,Int},Array{ComplexF64,4}}}() + for i in 1:12 + filename = joinpath(artifact_path, "A4", "Fsymbol_$i.txt") + @debug "loading $filename" @assert isfile(filename) "cannot find $filename" - json_string = read(filename, String) - Farray_part = copy(JSON3.read(json_string)) - return map(enumerate(Farray_part[Symbol(i)])) do (I, x) - j, k, l = Tuple(CartesianIndices((12, 12, 12))[I]) - y = Dict{NTuple{6,Int},Array{ComplexF64,4}}() - for (key, v) in x - a, b, c, d, e, f = parse.(Int, split(string(key)[2:(end - 1)], ", ")) + Farray_part = readdlm(filename) + for ((i, j, k, l), colordict) in convert_Fs(Farray_part) + result[(i, j, k, l)] = Dict{NTuple{6,Int},Array{ComplexF64,4}}() + for ((a, b, c, d, e, f), Fvals) in colordict a_ob, b_ob, c_ob, d_ob, e_ob, f_ob = A4Object.(((i, j, a), (j, k, b), (k, l, c), (i, l, d), (i, k, e), (j, l, f))) - result = Array{ComplexF64,4}(undef, - (Nsymbol(a_ob, b_ob, e_ob), - Nsymbol(e_ob, c_ob, d_ob), - Nsymbol(b_ob, c_ob, f_ob), - Nsymbol(a_ob, f_ob, d_ob))) - map!(result, reshape(v, size(result))) do cmplxdict - return complex(cmplxdict[:re], cmplxdict[:im]) + result[(i, j, k, l)][(a, b, c, d, e, f)] = zeros(ComplexF64, + Nsymbol(a_ob, b_ob, e_ob), + Nsymbol(e_ob, c_ob, d_ob), + Nsymbol(b_ob, c_ob, f_ob), + Nsymbol(a_ob, f_ob, d_ob)) + for (I, v) in Fvals + result[(i, j, k, l)][(a, b, c, d, e, f)][I] = v end - - y[(a, b, c, d, e, f)] = result end end end + return result +end + +function convert_Fs(Farray_part::Matrix{Float64}) # Farray_part is a matrix with 16 columns + data_dict = Dict{NTuple{4,Int}, + Dict{NTuple{6,Int},Vector{Pair{CartesianIndex{4},ComplexF64}}}}() + # want to make a Dict with keys (i,j,k,l) and vals + # a Dict with keys (a,b,c,d,e,f) and vals + # a pair of (mu, nu, rho, sigma) and the F value + for row in eachrow(Farray_part) + i, j, k, l, a, b, c, d, e, f, mu, nu, rho, sigma = Int.(@view(row[1:14])) + v = complex(row[15], row[16]) + colordict = get!(data_dict, (i, j, k, l), + Dict{NTuple{6,Int},Vector{Pair{CartesianIndex{4},ComplexF64}}}()) + Fdict = get!(colordict, (a, b, c, d, e, f), + Vector{Pair{CartesianIndex{4},ComplexF64}}()) + push!(Fdict, CartesianIndex(mu, nu, rho, sigma) => v) + end + return data_dict end const Fcache = IdDict{Type{<:BimoduleSector}, - Array{Dict{NTuple{6,Int},Array{ComplexF64,4}},4}}() + Dict{NTuple{4,Int64},Dict{NTuple{6,Int64},Array{ComplexF64,4}}}}() function _get_Fcache(::Type{T}) where {T<:BimoduleSector} global Fcache return get!(Fcache, T) do + @debug "loading Fsymbol cache for $T" return extract_Fsymbol(T) end end function TensorKitSectors.Fsymbol(a::I, b::I, c::I, d::I, e::I, f::I) where {I<:A4Object} - # TODO: should this error or return 0? - (a.j == b.i && b.j == c.i && a.i == d.i && c.j == d.j && - a.i == e.i && b.j == e.j && f.i == a.j && f.j == c.j) || - throw(ArgumentError("invalid fusion channel")) + # required to keep track of multiplicities where F-move is partially unallowed + # also deals with invalid fusion channels + Nabe = Nsymbol(a, b, e) + Necd = Nsymbol(e, c, d) + Nbcf = Nsymbol(b, c, f) + Nafd = Nsymbol(a, f, d) + + Nabe > 0 && Necd > 0 && Nbcf > 0 && Nafd > 0 || + return zeros(sectorscalartype(I), Nabe, Necd, Nbcf, Nafd) i, j, k, l = a.i, a.j, b.j, c.j - return get(_get_Fcache(I)[i, j, k, l], - (a.label, b.label, c.label, d.label, e.label, f.label)) do - return zeros(sectorscalartype(A4Object), - (Nsymbol(a, b, e), Nsymbol(e, c, d), Nsymbol(b, c, f), - Nsymbol(a, f, d))) + colordict = _get_Fcache(I)[i, j, k, l] + return colordict[(a.label, b.label, c.label, d.label, e.label, f.label)] +end + +# interface with TensorKit where necessary +#----------------------------------------- + +function TensorKit.blocksectors(W::TensorMapSpace{S,N₁,Nβ‚‚}) where + {S<:Union{GradedSpace{A4Object,NTuple{486,Int64}}, + SumSpace{GradedSpace{A4Object,NTuple{486,Int64}}}},N₁,Nβ‚‚} + codom = codomain(W) + dom = domain(W) + # @info "in the correct blocksectors" + if N₁ == 0 && Nβ‚‚ == 0 # 0x0-dimensional TensorMap is just a scalar, return all units + # this is a problem in full contractions where the coloring outside is π’ž + return NTuple{12,A4Object}(one(A4Object(i, i, 1)) for i in 1:12) # have to return all units b/c no info on W in this case + elseif N₁ == 0 + @assert Nβ‚‚ != 0 "one of Type A4Object doesn't exist" + return filter!(isone, collect(blocksectors(dom))) + elseif Nβ‚‚ == 0 + @assert N₁ != 0 "one of Type A4Object doesn't exist" + return filter!(isone, collect(blocksectors(codom))) + elseif Nβ‚‚ <= N₁ # keep intersection + return filter!(c -> hasblock(codom, c), collect(blocksectors(dom))) + else + return filter!(c -> hasblock(dom, c), collect(blocksectors(codom))) + end +end + +# TODO: definition for zero of GradedSpace? + +function dim(V::GradedSpace{<:BimoduleSector}) + T = Base.promote_op(*, Int, real(sectorscalartype(sectortype(V)))) + return reduce(+, dim(V, c) * dim(c) for c in sectors(V); init=zero(T)) +end + +# limited oneunit +function Base.oneunit(S::GradedSpace{<:BimoduleSector}) + allequal(a.i for a in sectors(S)) && allequal(a.j for a in sectors(S)) || + throw(ArgumentError("sectors of $S are not all equal")) + first(sectors(S)).i == first(sectors(S)).j || + throw(ArgumentError("sectors of $S are non-diagonal")) + sector = one(first(sectors(S))) + return β„‚[A4Object](sector => 1) +end + +function Base.oneunit(S::SumSpace{<:GradedSpace{<:BimoduleSector}}) + allequal(a.i for a in sectors(S)) && allequal(a.j for a in sectors(S)) || + throw(ArgumentError("sectors of $S are not all equal")) + first(sectors(S)).i == first(sectors(S)).j || + throw(ArgumentError("sectors of $S are non-diagonal")) + sector = one(first(sectors(S))) + return SumSpace(β„‚[A4Object](sector => 1)) +end + +# maybe from the homspace +function TensorKit.insertrightunit(P::ProductSpace{V,N}, ::Val{i}=Val(length(P)); + conj::Bool=false, + dual::Bool=false) where {i,V<:GradedSpace{<:I},N} where {I<:BimoduleSector} + #possible change to rightone of correct space for N = 0 + u = N > 0 ? oneunit(P[1]) : error("no unit object in this space") + if dual + u = TensorKit.dual(u) + end + if conj + u = TensorKit.conj(u) + end + return ProductSpace(TupleTools.insertafter(P.spaces, i, (u,))) +end + +function TensorKit.insertleftunit(P::ProductSpace{V,N}, ::Val{i}=Val(length(P) + 1); + conj::Bool=false, + dual::Bool=false) where {i,V<:GradedSpace{<:I},N} where {I<:BimoduleSector} + u = N > 0 ? oneunit(P[1]) : error("no unit object in this space") + if dual + u = TensorKit.dual(u) + end + if conj + u = TensorKit.conj(u) end + return ProductSpace(TupleTools.insertafter(P.spaces, i - 1, (u,))) end diff --git a/test/caching.jl b/test/caching.jl new file mode 100644 index 0000000..c385c86 --- /dev/null +++ b/test/caching.jl @@ -0,0 +1,60 @@ +using MultiTensorKit +using TensorKit +using MPSKit, MPSKitModels + +C1 = A4Object(1,1,1) +C0 = A4Object(1,1,4) # unit +M = A4Object(1,2,1) +D0 = A4Object(2,2,12) # unit +D1 = A4Object(2,2,2) # self-dual object + +P = Vect[A4Object](D0 => 1, D1 => 1) +h = TensorMap(ones, ComplexF64, P βŠ— P ← P βŠ— P) + +lattice = InfiniteChain(1) +H = @mpoham -sum(h{i,j} for (i,j) in nearest_neighbours(lattice)); + +D = 4 +V = Vect[A4Object](M => D); +inf_init = InfiniteMPS([P], [V]); + +ψ, envs = find_groundstate(inf_init, H, VUMPS(verbosity=3, tol=1e-12, maxiter=20)); + +# basic TensorKit tests + +# diagonal +obj = A4Object(2,2,1) +obj2 = A4Object(2,2,2) +sp = Vect[A4Object](obj=>1, obj2=>1) +A = TensorMap(ones, ComplexF64, sp βŠ— sp ← sp βŠ— sp) +transpose(A, (2,4,), (1,3,)) + +blocksectors(sp βŠ— sp) +@plansor fullcont[] := A[a b;a b] # 12 fusiontrees + +# π’ž x β„³ +obj = A4Object(1,1,1) +obj2 = A4Object(1,2,1) + +sp = Vect[A4Object](obj=>1) +sp2 = Vect[A4Object](obj2=>1) +TensorMap(rand, ComplexF64, sp βŠ— sp2 ← sp) # should throw ArgumentError +homspace = sp βŠ— sp2 ← sp2 +A = TensorMap(ones, ComplexF64, homspace) +permute(space(A),((1,),(3,2))) +transpose(A, (1,2,), (3,)) == A +transpose(A, (3,1,), (2,)) + +Aop = TensorMap(ones, ComplexF64, conj(sp2) βŠ— sp ← conj(sp2)) +transpose(Aop, (1,2,), (3,)) == Aop +transpose(Aop, (1,), (3,2)) + +@plansor Acont[a] := A[a b;b] # should not have data bc sp isn't the unit + +spfix = Vect[A4Object](one(obj)=>1) +Afix = TensorMap(ones, ComplexF64, spfix βŠ— sp2 ← sp2) +@plansor Acontfix[a] := Afix[a b;b] # should have a fusion tree + +blocksectors(sp βŠ— sp2) +A = TensorMap(ones, ComplexF64, sp βŠ— sp2 ← sp βŠ— sp2) +@plansor fullcont[] := A[a b;a b] # 12 fusiontrees \ No newline at end of file diff --git a/test/localtests.jl b/test/localtests.jl new file mode 100644 index 0000000..f8ba2bb --- /dev/null +++ b/test/localtests.jl @@ -0,0 +1,549 @@ +using TensorKitSectors +using MultiTensorKit +using Revise + +testobj = A4Object(1, 1, 1) # fusion cat object +unit = one(testobj) +collect(testobj βŠ— unit) +@assert unit == leftone(testobj) == rightone(testobj) + +testobj2 = A4Object(2, 2, 1) +unit2 = one(testobj2) +collect(testobj2 βŠ— unit2) +@assert unit2 == leftone(testobj2) == rightone(testobj2) + +testmodobj = A4Object(1, 2, 1) +one(testmodobj) +leftone(testmodobj) +rightone(testmodobj) + +Fsymbol(testobj, testobj, A4Object(1, 1, 3), testobj, A4Object(1, 1, 3), A4Object(1, 1, 4)) + +# using Artifacts +# using DelimitedFiles + +# artifact_path = joinpath(artifact"fusiondata", "MultiTensorKit.jl-data-v0.1.1") +# filename = joinpath(artifact_path, "A4", "Fsymbol_4.txt") +# txt_string = read(filename, String) +# F_arraypart = copy(readdlm(IOBuffer(txt_string))); + +# F_arraypart = MultiTensorKit.convert_Fs(F_arraypart) +# i,j,k,l = (4,12,12,2) # 5,2,8,10 +# a,b,c,d,e,f = (1, 11, 3, 1, 1, 3) #(2,1,1,2,2,3) +# testF = F_arraypart[(i,j,k,l)][(a,b,c,d,e,f)] +# a_ob, b_ob, c_ob, d_ob, e_ob, f_ob = A4Object.(((i, j, a), (j, k, b), +# (k, l, c), (i, l, d), +# (i, k, e), (j, l, f))) +# result = Array{ComplexF64,4}(undef, +# (Nsymbol(a_ob, b_ob, e_ob), +# Nsymbol(e_ob, c_ob, d_ob), +# Nsymbol(b_ob, c_ob, f_ob), +# Nsymbol(a_ob, f_ob, d_ob))) + +# map!(result, reshape(testF, size(result))) do pair +# return pair[2] +# end + +N = MultiTensorKit._get_Ncache(A4Object); + +duals = MultiTensorKit._get_dual_cache(A4Object)[2] +# checking duals is correct +for i in 1:12, j in 1:12 + for (index, a) in enumerate(duals[i, j]) + aob = A4Object(i, j, index) + bob = A4Object(j, i, a) + leftone(aob) ∈ aob βŠ— bob && rightone(aob) ∈ bob βŠ— aob || @show i, j, aob, bob + end +end + +A = MultiTensorKit._get_Fcache(A4Object) +a, b, c, d, e, f = A4Object(1, 1, 3), A4Object(1, 1, 2), A4Object(1, 1, 2), + A4Object(1, 1, 2), A4Object(1, 1, 2), A4Object(1, 1, 2) +a, b, c, d, e, f = A4Object(1, 1, 1), A4Object(1, 1, 2), A4Object(1, 1, 2), + A4Object(1, 1, 1), A4Object(1, 1, 2), A4Object(1, 1, 4) +coldict = A[a.i, a.j, b.j, c.j] +bla = get(coldict, (a.label, b.label, c.label, d.label, e.label, f.label)) do + return coldict[(a.label, b.label, c.label, d.label, e.label, f.label)] +end + +###### TensorKit stuff ###### +using MultiTensorKit +using TensorKit +using Test + +# π’ž x π’ž example + +obj = A4Object(2, 2, 1) +obj2 = A4Object(2, 2, 2) +sp = Vect[A4Object](obj => 1, obj2 => 1) +A = TensorMap(ones, ComplexF64, sp βŠ— sp ← sp βŠ— sp) +transpose(A, (2, 4), (1, 3)) + +blocksectors(sp βŠ— sp) +@plansor fullcont[] := A[a b; a b] # problem here is that fusiontrees for all 12 units are given + +# π’ž x β„³ example +obj = A4Object(1, 1, 1) +obj2 = A4Object(1, 2, 1) + +sp = Vect[A4Object](obj => 1) +sp2 = Vect[A4Object](obj2 => 1) +@test_throws ArgumentError("invalid fusion channel") TensorMap(rand, ComplexF64, + sp βŠ— sp2 ← sp) +homspace = sp βŠ— sp2 ← sp2 +A = TensorMap(ones, ComplexF64, homspace) +fusiontrees(A) +permute(space(A), ((1,), (3, 2))) +transpose(A, (1, 2), (3,)) == A +transpose(A, (3, 1), (2,)) + +Aop = TensorMap(ones, ComplexF64, conj(sp2) βŠ— sp ← conj(sp2)) +transpose(Aop, (1, 2), (3,)) == Aop +transpose(Aop, (1,), (3, 2)) + +@plansor Acont[a] := A[a b; b] # should not have data bc sp isn't the unit + +spfix = Vect[A4Object](one(obj) => 1) +Afix = TensorMap(ones, ComplexF64, spfix βŠ— sp2 ← sp2) +@plansor Acontfix[a] := Afix[a b; b] # should have a fusion tree + +blocksectors(sp βŠ— sp2) +A = TensorMap(ones, ComplexF64, sp βŠ— sp2 ← sp βŠ— sp2) +@plansor fullcont[] := A[a b; a b] # same 12 fusiontrees problem + +# completely off-diagonal example + +obj = A4Object(5, 4, 1) +obj2 = A4Object(4, 5, 1) +sp = Vect[A4Object](obj => 1) +sp2 = Vect[A4Object](obj2 => 1) +conj(sp) == sp2 + +A = TensorMap(ones, ComplexF64, sp βŠ— sp2 ← sp βŠ— sp2) +Aop = TensorMap(ones, ComplexF64, sp2 βŠ— sp ← sp2 βŠ— sp) + +At = transpose(A, (2, 4), (1, 3)) +Aopt = transpose(Aop, (2, 4), (1, 3)) + +blocksectors(At) == blocksectors(Aop) +blocksectors(Aopt) == blocksectors(A) + +@plansor Acont[] := A[a b; a b] # ignore this error for now +@plansor Acont2[] := A[b a; b a] + +testsp = SU2Space(0 => 1, 1 => 1) +Atest = TensorMap(ones, ComplexF64, testsp βŠ— testsp ← testsp βŠ— testsp) +@plansor Aconttest[] := Atest[a b; a b] + +# π’ž x β„³ ← β„³ x π’Ÿ +c = A4Object(1, 1, 1) +m = A4Object(1, 2, 1) +d = A4Object(2, 2, 1) +W = Vect[A4Object](c => 1) βŠ— Vect[A4Object](m => 1) ← + Vect[A4Object](m => 1) βŠ— Vect[A4Object](d => 1) + +# bram stuff +using TensorKitSectors +for i in 1:12, j in 1:12 + for a in A4Object.(i, j, MultiTensorKit._get_dual_cache(A4Object)[2][i, j]) + F = Fsymbol(a, dual(a), a, a, leftone(a), rightone(a))[1, 1, 1, 1] + isapprox(F, frobeniusschur(a) / dim(a); atol=1e-15) || + @show a, F, frobeniusschur(a) / dim(a) # check real + isreal(frobeniusschur(a)) || isapprox(abs(frobeniusschur(a)), 1.0; atol=1e-15) || + @show a, frobeniusschur(a), abs(frobeniusschur(a)) + end +end + +for i in 1:12, j in 1:12 # 18a + i != j || continue + objsij = A4Object.(i, j, MultiTensorKit._get_dual_cache(A4Object)[2][i, j]) + @assert all(dim(m) > 0 for m in objsij) +end + +for i in 1:12, j in 1:12 # 18b + objsii = A4Object.(i, i, MultiTensorKit._get_dual_cache(A4Object)[2][i, i]) + objsij = A4Object.(i, j, MultiTensorKit._get_dual_cache(A4Object)[2][i, j]) + + Ndict = Dict{Tuple{A4Object,A4Object,A4Object},Int}() + for a in objsii, m in objsij + for n in a βŠ— m + Ndict[(a, m, n)] = Nsymbol(a, m, n) + end + end + + for a in objsii, m in objsij + isapprox(dim(a) * dim(m), sum(Ndict[(a, m, n)] * dim(n) for n in a βŠ— m); + atol=2e-9) || @show a, m + end +end + +for i in 1:12, j in 1:12 # 18c + objsii = A4Object.(i, i, MultiTensorKit._get_dual_cache(A4Object)[2][i, i]) + objsij = A4Object.(i, j, MultiTensorKit._get_dual_cache(A4Object)[2][i, j]) + m_dimsum = sum(dim(m)^2 for m in objsij) + c_dimsum = sum(dim(c)^2 for c in objsii) + isapprox(m_dimsum, c_dimsum; atol=1e-8) || @show i, j, c_dimsum, m_dimsum +end + +(a, b, c, d, e, f) = (A4Object(2, 1, 1), A4Object(1, 2, 1), A4Object(2, 2, 11), + A4Object(2, 2, 11), A4Object(2, 2, 9), A4Object(1, 2, 1)) +Fsymbol(a, b, c, d, e, f) +zeros(ComplexF64, Nsymbol(a, b, e), Nsymbol(e, c, d), Nsymbol(b, c, f), Nsymbol(a, f, d)) == +Fsymbol(a, b, c, d, e, f) + +# testing blocksectors +using BlockTensorKit +W = Vect[A4Object](A4Object(2, 2, 12) => 1) ← + ProductSpace{GradedSpace{A4Object,NTuple{486,Int64}},0}() +W isa TensorMapSpace{GradedSpace{A4Object,NTuple{486,Int64}}} +W isa HomSpace +W isa HomSpace{S} where {S<:SumSpace{GradedSpace{A4Object,NTuple{486,Int64}}}} +W isa TensorMapSpace{SumSpace{GradedSpace{A4Object,NTuple{486,Int64}}}} + +# this appears as well (N1=N2=0) +W = ProductSpace{SumSpace{GradedSpace{A4Object,NTuple{486,Int64}}},0}() ← + ProductSpace{SumSpace{GradedSpace{A4Object,NTuple{486,Int64}}},0}() +typeof(W) +W isa TensorMapSpace{S} where {S<:GradedSpace{A4Object,NTuple{486,Int64}}} +W isa HomSpace{S} where {S<:GradedSpace{A4Object,NTuple{486,Int64}}} +W isa HomSpace +W isa TensorMapSpace + +W isa TensorMapSpace{SumSpace{GradedSpace{A4Object,NTuple{486,Int64}}},0,0} +W isa TensorMapSpace{GradedSpace{A4Object,NTuple{486,Int64}},0,0} +TensorMapSpace{SumSpace{GradedSpace{A4Object,NTuple{486,Int64}}},0,0} <: +TensorMapSpace{GradedSpace{A4Object,NTuple{486,Int64}},N₁,Nβ‚‚} where {N₁,Nβ‚‚} +W isa TensorSpace{SumSpace{GradedSpace{A4Object,NTuple{486,Int64}}}} + +GradedSpace{A4Object,NTuple{486,Int64}} <: +BlockTensorKit.SumSpace{GradedSpace{A4Object,NTuple{486,Int64}}} + +############ MPSKit wow ############ +using MultiTensorKit +using TensorKit +using MPSKit, MPSKitModels +C1 = A4Object(1, 1, 1) +C0 = A4Object(1, 1, 4) # unit +M = A4Object(1, 2, 1) +D0 = A4Object(2, 2, 12) # unit +D1 = A4Object(2, 2, 2) # self-dual object +collect(D0 βŠ— D1) +collect(D1 βŠ— D1) + +P = Vect[A4Object](D0 => 1, D1 => 1) +h = TensorMap(ones, ComplexF64, P βŠ— P ← P βŠ— P) + +using Profile +Profile.init(; delay=0.1) +lattice = InfiniteChain(1) +t = time() +H = @mpoham -sum(h{i,j} for (i, j) in nearest_neighbours(lattice)); # 15min, 10.5min +dt = time() - t +println("Time to create Hamiltonian: ", dt, " seconds") + +# testing insertleft/rightunit +sp = SU2Space(0 => 1, 1 => 1) +ht = TensorMap(ones, ComplexF64, P ← P) +htl = TensorMap(ones, ComplexF64, P ← one(P)) +htr = TensorMap(ones, ComplexF64, one(P) ← P) +htnone = TensorMap(ones, ComplexF64, one(P) ← one(P)) +insertrightunit(htr) # adding to empty space +insertleftunit(htr) +insertrightunit(htl) +insertleftunit(htl) # adding to empty space +insertrightunit(htnone) +insertleftunit(htnone) + +D = 4 +V = Vect[A4Object](M => D); +t = time() +inf_init = InfiniteMPS([P], [V]); # 8min, 6.5min +dt = time() - t +println("Time to create InfiniteMPS: ", dt, " seconds") + +# VUMPS +t = time() +ψ, envs = find_groundstate(inf_init, H, VUMPS(; verbosity=3, tol=1e-12, maxiter=20)); # 40 min, 30min +dt = time() - t +println("Time to find groundstate: ", dt, " seconds") +expectation_value(ψ, H, envs) + +using JLD2 +# jldsave("C:/Boris/Unief/PhD/Code/MultiTensorKit/Saves/testH.jld2"; H) +# jldsave("C:/Boris/Unief/PhD/Code/MultiTensorKit/Saves/testinf_init.jld2"; inf_init) +# jldsave("C:/Boris/Unief/PhD/Code/MultiTensorKit/Saves/testgs.jld2"; ψ) +# jldsave("C:/Boris/Unief/PhD/Code/MultiTensorKit/Saves/testenvs.jld2"; envs) + +ψ = jldopen("C:/Boris/Unief/PhD/Code/MultiTensorKit/Saves/testgs.jld2", "r") do file + return read(file, "ψ") +end +H = jldopen("C:/Boris/Unief/PhD/Code/MultiTensorKit/Saves/testH.jld2", "r") do file + return read(file, "H") +end +inf_init = jldopen("C:/Boris/Unief/PhD/Code/MultiTensorKit/Saves/testinf_init.jld2", + "r") do file + return read(file, "inf_init") +end +envs = jldopen("C:/Boris/Unief/PhD/Code/MultiTensorKit/Saves/testenvs.jld2", "r") do file + return read(file, "envs") +end + +entropy(ψ) +entanglement_spectrum(ψ) +transfer_spectrum(ψ; sector=C0) +correlation_length(ψ; sector=C0) +norm(ψ) + +# test correlator + +#IDMRG + +ψ, envs = find_groundstate(inf_init, H, IDMRG(; verbosity=3, tol=1e-8, maxiter=15)); +expectation_value(ψ, H, envs) + +#IDMRG2 +inf_init2 = InfiniteMPS([P, P], [V, V]) +H2 = @mpoham -sum(h{i,j} for (i, j) in nearest_neighbours(InfiniteChain(2))); +idmrg2alg = IDMRG2(; verbosity=3, tol=1e-8, maxiter=15, trscheme=truncdim(10)) +ψ2, envs2 = find_groundstate(inf_init2, H2, idmrg2alg); +expectation_value(ψ2, H2, envs2) + +#QuasiParticleAnsatz + +momenta = range(0, 2Ο€, 3) +t = time() +excE, excqp = excitations(H, QuasiparticleAnsatz(; ishermitian=true), momenta, ψ, envs; + sector=C0, num=1, parallel=false); # not working for some reason +dt = time() - t +# 15min, 100min +println("Time to create excitations: ", dt, " seconds") + +####################################################### +# debugging QPA for infinite systems +momentum = momenta[begin] +Ο•β‚€ = LeftGaugedQP(rand, ψ, ψ; sector=C0, momentum=momenta[begin]); +E = MPSKit.effective_excitation_renormalization_energy(H, Ο•β‚€, envs, envs) +H_eff = MPSKit.EffectiveExcitationHamiltonian(H, envs, envs, E); # function that acts on QP +alg = QuasiparticleAnsatz(; ishermitian=true) + +# do block with eigsolve effectively doing this +my_operator(Ο•) = H_eff(Ο•; alg.alg_environments...) +Es, Ο•s, convhist = eigsolve(my_operator, Ο•β‚€, 1, :SR, alg.alg) + +# goes back in the eigsolve during KrylovKit.initialise +iter = LanczosIterator(my_operator, Ο•β‚€, alg.alg.orth); +fact = initialize(iter; verbosity=alg.alg.verbosity) + +# goes wrong during apply of KrylovKit.initialise +Axβ‚€ = KrylovKit.apply(iter.operator, iter.xβ‚€) + +# apply(f,x) = f(x) calculates EffectiveExcitationHamiltonian(H, envs, envs, E) of QP +iter.operator(iter.xβ‚€) +MPSKit.EffectiveExcitationHamiltonian(H, envs, envs, E)(Ο•β‚€; alg.alg_environments...) +H_eff(Ο•β‚€; alg.alg_environments...) + +# error occurs when calculating the environment +qp_envs = environments(Ο•β‚€, H, envs, envs; alg.alg_environments...) + +# goes wrong when calculating environments of QP +lBs = PeriodicVector([MPSKit.allocate_GBL(Ο•β‚€, H, Ο•β‚€, i) for i in 1:length(Ο•β‚€)]); +MPSKit.left_excitation_transfer_system(lBs[1], H, Ο•β‚€; solver=MPSKit.Defaults.linearsolver) + +# problem occurs at the linsolve which calls GMRES +found = zerovector(lBs[1]) +H_partial = map(h -> getindex(h, 1:1, 1, 1, 1:1), parent(H)) +T = TransferMatrix(exci.right_gs.AR, H_partial, exci.left_gs.AL) +start = scale!(last(found[1:1] * T), cis(-momenta[begin] * 1)) +if exci.trivial && isidentitylevel(H, i) + # not using braiding tensors here, leads to extra leg + util = similar(exci.left_gs.AL[1], space(parent(H)[1], 1)[1]) + fill_data!(util, one) + @plansor start[-1 -2; -3 -4] -= start[2 1; -3 3] * + util[1] * + r_RL(exci.right_gs)[3; 2] * + l_RL(exci.right_gs)[-1; -4] * + conj(util[-2]) +end +found[1] = add!(start, lBs[1]) + +T = TransferMatrix(exci.right_gs.AR, exci.left_gs.AL) +if exci.trivial + # deal with extra leg + @plansor lRL_util[-1 -2; -3] := l_RL(exci.right_gs)[-1; -3] * conj(util[-2]) + @plansor rRL_util[-1 -2; -3] := r_RL(exci.right_gs)[-1; -3] * util[-2] + T = regularize(T, lRL_util, rRL_util) +end + +found[1], convhist = linsolve(flip(T), found[1], found[1], MPSKit.Defaults.linearsolver, 1, + -cis(-momenta[begin] * 1)) + +############################################################################## +# quick test on complex f symbols and dimensions +testp = Vect[A4Object](one(A4Object(i, i, 1)) => 1 for i in 1:12) +dim(testp) +oneunit(testp) + +# finite stuff +L = 10 +lattice = FiniteChain(L) +P = Vect[A4Object](D0 => 1, D1 => 1) +D = 4 +V = Vect[A4Object](M => D) + +dmrgalg = DMRG(; verbosity=3, tol=1e-8, maxiter=100, + alg_eigsolve=MPSKit.Defaults.alg_eigsolve(; ishermitian=false)) +fin_init = FiniteMPS(L, P, V; left=V, right=V) +Hfin = @mpoham -sum(h{i,j} for (i, j) in nearest_neighbours(lattice)); +open_boundary_conditions(H, L) == Hfin +ψfin, envsfin = find_groundstate(fin_init, Hfin, dmrgalg); +expectation_value(ψfin, Hfin, envsfin) / (L - 1) + +entropy(ψfin, round(Int, L / 2)) +entanglement_spectrum(ψfin, round(Int, L / 2)) +Es, states, convhist = exact_diagonalization(Hfin; sector=D0); +Es / (L - 1) + +dmrg2alg = DMRG2(; verbosity=3, tol=1e-8, maxiter=15, + alg_eigsolve=MPSKit.Defaults.alg_eigsolve(; ishermitian=false), + trscheme=truncdim(10)) +ψfin2, envsfin2 = find_groundstate(fin_init, Hfin, dmrg2alg); +0expectation_value(ψfin2, Hfin, envsfin2) / (L - 1) + +entropy(ψfin2, round(Int, L / 2)) +entanglement_spectrum(ψfin2, round(Int, L / 2)) + +S = left_virtualspace(Hfin, 1) +oneunit(S) +eltype(S) +oneunit(eltype(S)) # should error + +# excitations +qpaalg = QuasiparticleAnsatz(; ishermitian=false, tol=1e-8, maxiter=100) +excEfin, excqpfin = excitations(Hfin, qpaalg, ψfin, envsfin; sector=C0, num=1); +excEfin + +excFIN, excqpFIN = excitations(Hfin, FiniteExcited(; gsalg=dmrg2alg), ψfin; num=1); +excFIN + +# changebonds test +dim(left_virtualspace(ψ, 1)) +ψch, envsch = changebonds(ψ, H, OptimalExpand(; trscheme=truncerr(1e-3)), envs) +dim(left_virtualspace(ψch, 1)) + +# time evolution + +ψt, envst = timestep(ψ, H, 3, 0, + TDVP(; integrator=MPSKit.Defaults.alg_expsolve(; ishermitian=true)), + envs); +et = expectation_value(ψt, H, envst) +e = expectation_value(ψ, H, envs) +isapprox(et, e; atol=1e-12) # hermitian + +tdvpalg = TDVP(; integrator=MPSKit.Defaults.alg_expsolve(; ishermitian=true, krylovdim=10)) +ψt, envst = time_evolve(ψ, H, range(0, 1, 10), tdvpalg, envs) + +timealg = WII(; tol=1e-8, maxiter=100) +make_time_mpo(H, 0.1, timealg) + +# stat-mech stuff +mpo = InfiniteMPO([h]) +ψ, envs = leading_boundary(inf_init, mpo, + VUMPS(; verbosity=3, tol=1e-8, maxiter=15, + alg_eigsolve=MPSKit.Defaults.alg_eigsolve(; + ishermitian=false))); + +# addition, substraction, multiplication + +# finite +Hfin * ψfin + +# infinite +H * ψ # currently doesn't work +MPSKit.DenseMPO(H) * ψ # does work + +# approximate +ψa, _ = approximate(ψ, (mpo, inf_init), IDMRG()); + +mpo2 = InfiniteMPO([h, h]) +inf_init2 = InfiniteMPS([P, P], [V, V]) +H2 = @mpoham -sum(h{i,j} for (i, j) in nearest_neighbours(InfiniteChain(2))); +ψ2, envs2 = find_groundstate(inf_init2, H2, VUMPS(; verbosity=3, tol=1e-10, maxiter=15)); +ψa, _ = approximate(ψ2, (mpo2, inf_init2), IDMRG2(; trscheme=truncdim(10))); + +# # testing InfiniteMPOHamiltonian and FiniteMPOHamiltonian constructor not relying on MPSKitModels + +sp = Vect[FibonacciAnyon](:I => 1, :Ο„ => 1) +t = TensorMap(ones, ComplexF64, sp ← sp) +InfiniteMPOHamiltonian(PeriodicArray([sp]), i => t for i in 1:1) +H = @mpoham -sum(h{i,j} for (i, j) in nearest_neighbours(InfiniteChain(1))); + +############################################## +# diagonal case + +D = 4 +V = Vect[A4Object](D0 => D, D1 => D); +inf_init = InfiniteMPS([P], [V]); + +ψ, envs = find_groundstate(inf_init, H, VUMPS(; verbosity=3, tol=1e-10, maxiter=500)); + +# other module categories +# 4,5,7,8,9,10,11,12 gives poorly converged envs, +k = 6 # 6 gives lapackexception(22) +V = Vect[A4Object](A4Object(k, 2, i) => 2 + for i in 1:MultiTensorKit._numlabels(A4Object, k, 2)) +inf_init = InfiniteMPS([P], [V]); +# expectation_value(inf_init, H) + +ψ, envs = find_groundstate(inf_init, H, VUMPS(; verbosity=3, tol=1e-10, maxiter=10)); +expectation_value(ψ, H, envs) + +# checking multiplicity +function obs(i::Int) + return A4Object.(i, i, MultiTensorKit._get_dual_cache(A4Object)[2][i, i]) +end + +for i in 1:12 + !any(Nsymbol(a, b, c) > 1 for a in obs(i), b in obs(i), c in obs(i)) && @show i +end + +# trying to make heisenberg +using MultiTensorKit, TensorKit +using MPSKit, MPSKitModels + +D1 = A4Object(1, 1, 2) # 3-dimensional irrep of A4 +M = A4Object(2, 1, 1) # Vec + +P = Vect[A4Object](D1 => 1) +h_aux1 = TensorMap(ones, ComplexF64, P ← P βŠ— P) +h_aux2 = TensorMap(ones, ComplexF64, P βŠ— P ← P) + +@plansor h[-1 -2; -3 -4] := h_aux2[-1 1; -3] * h_aux1[-2; 1 -4] # different basis +lattice = InfiniteChain(1) +t = time() +H = @mpoham -sum(h{i,j} for (i, j) in nearest_neighbours(lattice)); +dt = time() - t +println("Time to create Hamiltonian: ", dt, " seconds") + +D = 4 +V = Vect[A4Object](M => D) +t = time() +inf_init = InfiniteMPS([P], [V]); +dt = time() - t +println("Time to create InfiniteMPS: ", dt, " seconds") + +# VUMPS +t = time() +ψ, envs = find_groundstate(inf_init, H, VUMPS(; verbosity=3, tol=1e-12, maxiter=200)); +dt = time() - t +println("Time to find groundstate: ", dt, " seconds") +expectation_value(ψ, H, envs) # this gives 0 oopsie + +### caching checks +length(TensorKit.GLOBAL_FUSIONBLOCKSTRUCTURE_CACHE) + +length(TensorKit.treepermutercache) # tensor stuff +length(TensorKit.treetransposercache) +length(TensorKit.treebraidercache) + +length(TensorKit.transposecache) # fusion tree stuff +length(TensorKit.braidcache) diff --git a/test/test_A4.jl b/test/test_A4.jl index 696ca5b..81cb8bd 100644 --- a/test/test_A4.jl +++ b/test/test_A4.jl @@ -11,7 +11,7 @@ I = A4Object end @testset "Fusion Category $i" for i in 1:12 - objects = A4Object.(i, i, MultiTensorKit._get_dual_cache(I)[i][2]) + objects = A4Object.(i, i, MultiTensorKit._get_dual_cache(I)[2][i,i]) @testset "Basic properties" begin s = rand(objects, 3) @@ -19,9 +19,9 @@ end @test @constinferred(hash(s[1])) == hash(deepcopy(s[1])) @test isone(@constinferred(one(s[1]))) @constinferred dual(s[1]) - @constinferred dim(s[1]) - @constinferred frobeniusschur(s[1]) - @constinferred Bsymbol(s...) + @constinferred dim(s[1]) + @constinferred frobeniusschur(s[1]) + @constinferred Bsymbol(s...) # ill-defined test, doesn't necessarily exist and will error at dictionary keys @constinferred Fsymbol(s..., s...) end @@ -41,14 +41,36 @@ end end end F = hvcat(length(fs), Fblocks...) - @test isapprox(F' * F, one(F); atol=1e-12, rtol=1e-12) + @test isapprox(F' * F, one(F); atol=1e-9, rtol=1e-9) # some are simply not unitary? end end end +end - @testset "Pentagon equation" begin - for a in objects, b in objects, c in objects, d in objects - @test pentagon_equation(a, b, c, d; atol=1e-12, rtol=1e-12) +@testset "Pentagon equation" begin + objects = collect(values(A4Object)) + for a in objects + for b in objects + a.j == b.i || continue # skip if not compatible + for c in objects + b.j == c.i || continue # skip if not compatible + for d in objects + c.j == d.i || continue # skip if not compatible + @test pentagon_equation(a, b, c, d; atol=1e-9, rtol=1e-9) # ill-defined for same reason + end + end end end end + +@testset "A4 Category ($i, $j) units and duals" for i in 1:12, j in 1:12 + Cij_obs = A4Object.(i, j, MultiTensorKit._get_dual_cache(I)[2][i,j]) + + s = rand(Cij_obs, 1)[1] + @test eval(Meta.parse(sprint(show, s))) == s + @test @constinferred(hash(s)) == hash(deepcopy(s)) + @test i == j ? isone(@constinferred(one(s))) : (isone(@constinferred(leftone(s))) && isone(@constinferred(rightone(s)))) + @constinferred dual(s) + @test dual(s) == A4Object(j, i, MultiTensorKit._get_dual_cache(I)[2][i,j][s.label]) + @test dual(dual(s)) == s +end