Add 3-site simple update (aka 3-site cluster update)

.github/workflows/Tests.yml

@@ -32,6 +32,7 @@ jobs:
           - examples
           - utility
           - bondenv
+          - timeevol
         os:
           - ubuntu-latest
           - macOS-latest

docs/make.jl

@@ -43,7 +43,8 @@ mathengine = MathJax3(
 # examples pages
 examples_optimization =
     joinpath.(["heisenberg", "bose_hubbard", "xxz", "fermi_hubbard"], Ref("index.md"))
-examples_time_evolution = joinpath.(["heisenberg_su", "hubbard_su"], Ref("index.md"))
+examples_time_evolution =
+    joinpath.(["heisenberg_su", "hubbard_su", "j1j2_su"], Ref("index.md"))
 examples_partition_functions =
     joinpath.(
         ["2d_ising_partition_function", "3d_ising_partition_function"], Ref("index.md")

docs/src/assets/pepskit.bib

@@ -134,3 +134,16 @@ @article{vanderstraeten_variational_2022
   publisher = {American Physical Society},
   doi = {10.1103/PhysRevB.105.195140}
 }
+
+@article{liu_gapless_2022,
+  title = {Gapless quantum spin liquid and global phase diagram of the spin-1/2 J1-J2 square antiferromagnetic Heisenberg model},
+  journal = {Science Bulletin},
+  volume = {67},
+  number = {10},
+  pages = {1034-1041},
+  year = {2022},
+  issn = {2095-9273},
+  doi = {10.1016/j.scib.2022.03.010},
+  url = {https://www.sciencedirect.com/science/article/pii/S2095927322001001},
+  author = {Wen-Yuan Liu and Shou-Shu Gong and Yu-Bin Li and Didier Poilblanc and Wei-Qiang Chen and Zheng-Cheng Gu},
+}

docs/src/examples/j1j2_su/main.ipynb

@@ -0,0 +1,239 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "using Markdown #hide"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Three-site simple update for the $J_1$-$J_2$ model\n",
+    "\n",
+    "In this example, we will use `SimpleUpdate` imaginary time evolution to treat\n",
+    "the two-dimensional $J_1$-$J_2$ model, which contains next-nearest-nearest neighbor interactions:\n",
+    "\n",
+    "$$\n",
+    "H = J_1 \\sum_{\\langle i,j \\rangle} \\mathbf{S}_i \\cdot \\mathbf{S}_j\n",
+    "+ J_2 \\sum_{\\langle \\langle i,j \\rangle \\rangle} \\mathbf{S}_i \\cdot \\mathbf{S}_j\n",
+    "$$\n",
+    "\n",
+    "Here we will exploit the $U(1)$ spin rotation symmetry in the $J_1$-$J_2$ model. The goal\n",
+    "will be to calculate the energy at $J_1 = 1$ and $J_2 = 1/2$, first using the simple update\n",
+    "algorithm and then, to refine the energy estimate, using AD-based variational PEPS\n",
+    "optimization.\n",
+    "\n",
+    "We first import all required modules and seed the RNG:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "using Random\n",
+    "using TensorKit, PEPSKit\n",
+    "Random.seed!(2025);"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Simple updating a challenging phase\n",
+    "\n",
+    "Let's start by initializing an `InfiniteWeightPEPS` for which we set the required parameters\n",
+    "as well as physical and virtual vector spaces. We use the minimal unit cell size\n",
+    "($2 \\times 2$) required by the simple update algorithm for Hamiltonians with\n",
+    "next-nearest-neighbour interactions:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "Dbond, χenv, symm = 4, 32, U1Irrep\n",
+    "trscheme_env = truncerr(1e-10) & truncdim(χenv)\n",
+    "Nr, Nc, J1 = 2, 2, 1.0\n",
+    "\n",
+    "# random initialization of 2x2 iPEPS with weights and CTMRGEnv (using real numbers)\n",
+    "Pspace = Vect[U1Irrep](1//2 => 1, -1//2 => 1)\n",
+    "Vspace = Vect[U1Irrep](0 => 2, 1//2 => 1, -1//2 => 1)\n",
+    "Espace = Vect[U1Irrep](0 => χenv ÷ 2, 1//2 => χenv ÷ 4, -1//2 => χenv ÷ 4)\n",
+    "wpeps = InfiniteWeightPEPS(rand, Float64, Pspace, Vspace; unitcell=(Nr, Nc));"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "The value $J_2 / J_1 = 0.5$ corresponds to a [possible spin liquid phase](@cite liu_gapless_2022),\n",
+    "which is challenging for SU to produce a relatively good state from random initialization.\n",
+    "Therefore, we shall gradually increase $J_2 / J_1$ from 0.1 to 0.5, each time initializing\n",
+    "on the previously evolved PEPS:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "dt, tol, maxiter = 1e-2, 1e-8, 30000\n",
+    "check_interval = 4000\n",
+    "trscheme_peps = truncerr(1e-10) & truncdim(Dbond)\n",
+    "alg = SimpleUpdate(dt, tol, maxiter, trscheme_peps)\n",
+    "for J2 in 0.1:0.1:0.5\n",
+    "    H = real( ## convert Hamiltonian `LocalOperator` to real floats\n",
+    "        j1_j2_model(ComplexF64, symm, InfiniteSquare(Nr, Nc); J1, J2, sublattice=false),\n",
+    "    )\n",
+    "    result = simpleupdate(wpeps, H, alg; check_interval)\n",
+    "    global wpeps = result[1]\n",
+    "end"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "After we reach $J_2 / J_1 = 0.5$, we gradually decrease the evolution time step to obtain\n",
+    "a more accurately evolved PEPS:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "dts = [1e-3, 1e-4]\n",
+    "tols = [1e-9, 1e-9]\n",
+    "J2 = 0.5\n",
+    "H = real(j1_j2_model(ComplexF64, symm, InfiniteSquare(Nr, Nc); J1, J2, sublattice=false))\n",
+    "for (dt, tol) in zip(dts, tols)\n",
+    "    alg′ = SimpleUpdate(dt, tol, maxiter, trscheme_peps)\n",
+    "    result = simpleupdate(wpeps, H, alg′; check_interval)\n",
+    "    global wpeps = result[1]\n",
+    "end"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Computing the simple update energy estimate\n",
+    "\n",
+    "Finally, we measure the ground-state energy by converging a CTMRG environment and computing\n",
+    "the expectation value, where we make sure to normalize by the unit cell size:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "peps = InfinitePEPS(wpeps)\n",
+    "normalize!.(peps.A, Inf) ## normalize PEPS with absorbed weights by largest element\n",
+    "env₀ = CTMRGEnv(rand, Float64, peps, Espace)\n",
+    "env, = leading_boundary(env₀, peps; tol=1e-10, alg=:sequential, trscheme=trscheme_env);\n",
+    "E = expectation_value(peps, H, env) / (Nr * Nc)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Let us compare that estimate with benchmark data obtained from the\n",
+    "[YASTN/peps-torch package](https://github.com/jurajHasik/j1j2_ipeps_states/blob/ea4140fbd7da0fc1b75fac2871f75bda125189a8/single-site_pg-C4v-A1_internal-U1/j20.5/state_1s_A1_U1B_j20.5_D4_chi_opt96.dat).\n",
+    "which utilizes AD-based PEPS optimization to find $E_\\text{ref}=-0.49425$:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "E_ref = -0.49425\n",
+    "@show (E - E_ref) / abs(E_ref);"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Variational PEPS optimization using AD\n",
+    "\n",
+    "As a last step, we will use the SU-evolved PEPS as a starting point for a `fixedpoint`\n",
+    "PEPS optimization. Note that we could have also used a sublattice-rotated version of `H` to\n",
+    "fit the Hamiltonian onto a single-site unit cell which would require us to optimize fewer\n",
+    "parameters and hence lead to a faster optimization. But here we instead take advantage of\n",
+    "the already evolved `peps`, thus giving us a physical initial guess for the optimization:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "peps_opt, env_opt, E_opt, = fixedpoint(\n",
+    "    H, peps, env; optimizer_alg=(; tol=1e-4, maxiter=120)\n",
+    ");"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Finally, we compare the variationally optimized energy against the reference energy. Indeed,\n",
+    "we find that the additional AD-based optimization improves the SU-evolved PEPS and leads to\n",
+    "a more accurate energy estimate."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "E_opt /= (Nr * Nc)\n",
+    "@show E_opt\n",
+    "@show (E_opt - E_ref) / abs(E_ref);"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "---\n",
+    "\n",
+    "*This notebook was generated using [Literate.jl](https://github.com/fredrikekre/Literate.jl).*"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Julia 1.11.5",
+   "language": "julia",
+   "name": "julia-1.11"
+  },
+  "language_info": {
+   "file_extension": ".jl",
+   "mimetype": "application/julia",
+   "name": "julia",
+   "version": "1.11.5"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 3
+}

examples/Cache.toml

@@ -6,4 +6,5 @@ boundary_mps = "31e14bd742c12acfc4facd5bc8f28752f025a1a523c72ffccc8581e38880d3cd
 heisenberg_su = "299b79f36b94b301dafa0d39fc54498b44dcdf0e4988e3c5416aa6c4e9e8f584"
 xxz = "b5ab0efe6b9addb4b029250d7212d8b2e771136ca4ec4a4bfd7bccc8c688b949"
 fermi_hubbard = "13d6d827a9851bef2c47cbe8894bf6eec9e1fb3222ecc24d2a876336553439ff"
+j1j2_su = "883df8a1ccf061754538620c1086e2e13f1fd7821b8f2f489f305d0cf7098908"
 heisenberg = "f6b057dacb058656ad903778fc00832c6c8a61df930118872e04006b086030c7"