exercise_3_dmrg.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "2fe29049",
   "metadata": {},
   "source": [
    "# DMRG runs\n",
    "\n",
    "In this notebook, we use the `SimpleDMRGEngine` class from `tenpy_toycodes/d_dmrg.py` to run DMRG and find MPS ground states."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "56956a97",
   "metadata": {},
   "outputs": [],
   "source": [
    "# standard imports and cosmetics\n",
    "\n",
    "import numpy as np\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "np.set_printoptions(precision=5, suppress=True, linewidth=100, threshold=50)\n",
    "plt.rcParams['figure.dpi'] = 150"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "238d735d-7666-4181-bf50-481d8cc61adb",
   "metadata": {},
   "source": [
    "In previous notebooks, we learned how to initialize `SimpleMPS`..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8b2cfaee",
   "metadata": {},
   "outputs": [],
   "source": [
    "from tenpy_toycodes.a_mps import SimpleMPS, init_FM_MPS, init_Neel_MPS"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fd3a1a46-378d-4041-86ec-4ae992de0a7e",
   "metadata": {},
   "outputs": [],
   "source": [
    "L = 12"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a4f0ecc0",
   "metadata": {},
   "outputs": [],
   "source": [
    "psi_FM = init_FM_MPS(L=L, d=2, bc='finite')\n",
    "print(psi_FM)\n",
    "SigmaZ = np.diag([1., -1.])\n",
    "print(psi_FM.site_expectation_value(SigmaZ))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fe6d3e07-f52d-466e-913d-245e541a0ced",
   "metadata": {},
   "outputs": [],
   "source": [
    "def init_PM_MPS(L, bc='finite'):\n",
    "    \"\"\"Return a paramagnetic MPS (= product state with all spins pointing in +x direction)\"\"\"\n",
    "    d = 2\n",
    "    B = np.zeros([1, d, 1], dtype=float)\n",
    "    B[0, 0, 0] = 1./np.sqrt(2)\n",
    "    B[0, 1, 0] = 1./np.sqrt(2)\n",
    "    S = np.ones([1], dtype=float)\n",
    "    Bs = [B.copy() for i in range(L)]\n",
    "    Ss = [S.copy() for i in range(L)]\n",
    "    return SimpleMPS(Bs, Ss, bc=bc)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bcb5c97c-5f8b-4c1c-94c6-d3c2041db75c",
   "metadata": {},
   "source": [
    "... and how to initialize the `TFIModel`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d439f112",
   "metadata": {},
   "outputs": [],
   "source": [
    "from tenpy_toycodes.b_model import TFIModel\n",
    "\n",
    "g = 1.2\n",
    "model = TFIModel(L=L, J=1., g=g, bc='finite')\n",
    "\n",
    "print(\"<H_bonds> = \", psi_FM.bond_expectation_value(model.H_bonds))\n",
    "print(\"energy:\", np.sum(psi_FM.bond_expectation_value(model.H_bonds)))\n",
    "# (make sure the model and state have the same length and boundary conditions!)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "93368432-f304-4327-987c-b9c12a818c26",
   "metadata": {},
   "source": [
    "\n",
    "For small enough system size $L \\lesssim 16$, you can compare the energies to exact diagonalization in the full Hilbert space (which is exponentially expensive!):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "25871fc0-6b4f-4266-b78c-c25519210372",
   "metadata": {},
   "outputs": [],
   "source": [
    "from tenpy_toycodes.tfi_exact import finite_gs_energy\n",
    "\n",
    "if L <= 16:\n",
    "    energy_exact = finite_gs_energy(L=L, J=1., g=g)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "69be6188-dfad-4313-90dc-ca0be87793ad",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "bb7b318b",
   "metadata": {},
   "source": [
    "## The DMRG algorithm\n",
    "\n",
    "The file `tenpy_toycodes/d_dmrg.py` implements the DMRG algorithm.\n",
    "It can be called like this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3196f562",
   "metadata": {},
   "outputs": [],
   "source": [
    "from tenpy_toycodes.d_dmrg import SimpleDMRGEngine, SimpleHeff2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ffc67843",
   "metadata": {},
   "outputs": [],
   "source": [
    "chi_max = 15\n",
    "\n",
    "psi = init_FM_MPS(model.L, model.d, model.bc)\n",
    "eng = SimpleDMRGEngine(psi, model, chi_max=chi_max, eps=1.e-10)\n",
    "for i in range(10):\n",
    "    E_dmrg = eng.sweep()\n",
    "    E = np.sum(psi.bond_expectation_value(model.H_bonds))\n",
    "    print(\"sweep {i:2d}: E = {E:.13f}\".format(i=i + 1, E=E))\n",
    "print(\"final bond dimensions: \", psi.get_chi())\n",
    "mag_x = np.mean(psi.site_expectation_value(model.sigmax))\n",
    "mag_z = np.mean(psi.site_expectation_value(model.sigmaz))\n",
    "print(\"magnetization in X = {mag_x:.5f}\".format(mag_x=mag_x))\n",
    "print(\"magnetization in Z = {mag_z:.5f}\".format(mag_z=mag_z))\n",
    "if model.L <= 16:\n",
    "    E_exact = finite_gs_energy(L=model.L, J=model.J, g=model.g)\n",
    "    print(\"err in energy = {err:.3e}\".format(err=E - E_exact))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b0eee6bf",
   "metadata": {},
   "source": [
    "### Exercise: read d_dmrg.py\n",
    "\n",
    "Read the code of `tenpy_toycodes/d_dmrg.py` and try to undertstand the general structure of how it works.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "97226a24-892e-493a-88f4-9394615f240a",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "2a57a415-a81d-46f1-8366-533466ff70ee",
   "metadata": {},
   "source": [
    "### Exercise: measure correlation functions\n",
    "\n",
    "Just looking at expectation values of local operators is not enough.\n",
    "Also measure correlation functions $\\langle Z_{L/4} Z_{3L/4} \\rangle$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8fc2a15b-f45d-4073-aa56-7af10cb01bdc",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "76d1d8ed-c18d-4a84-8faf-f5390cd7e0b4",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "4799e6ed-f16d-4189-b224-8d1e972bf14d",
   "metadata": {},
   "source": [
    "### Exercise: DMRG runs\n",
    "\n",
    "Try running DMRG for various different parameters:\n",
    "\n",
    "- Change the bond dimension `chi` and truncation threashold `eps`.\n",
    "- Change the system size `L`\n",
    "- Change the model parameter `g` (at fixed $J=1$) to both the ferromagnetic phase $g<J$ \n",
    "  and the paramagnetic phase $g>1$.\n",
    "- Change the initial state.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c489cca3-3b21-4602-8329-6fb9e06137ff",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "efe930da-18c3-4091-8ea9-412b562c7857",
   "metadata": {},
   "source": [
    "### Exercise: Phase diagram\n",
    "\n",
    "To map out the phase diagram, it can be convenient to define a function that just runs DMRG for a given model. Fill in the below template. Use it obtain and plot the energy, correlations and magnetizations for different $L=16, 32, 64$ against $g$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8157051c-7bfa-407b-a889-536e821eedf6",
   "metadata": {},
   "outputs": [],
   "source": [
    "def run_DMRG(model, chi_max=50):\n",
    "    print(f\"runnning DMRG for L={model.L:d}, g={model.g:.2f}, bc={model.bc}, chi_max={chi_max:d}\")\n",
    "    \n",
    "    raise NotImplementedError(\"TODO: this is an Exercise!\")\n",
    "    \n",
    "    return psi"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b7eb8ba4-897e-495a-a22f-dc833cb246b5",
   "metadata": {},
   "outputs": [],
   "source": [
    "results_all = {}\n",
    "# TIP: comment this out after the first run to avoid overriding your data!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c1d120e7-0ecc-4800-857b-f12cd3fd23e7",
   "metadata": {},
   "outputs": [],
   "source": [
    "L = 8\n",
    "\n",
    "mag_X = []\n",
    "mag_Z = []\n",
    "E = []\n",
    "corr = []\n",
    "max_chi = []\n",
    "\n",
    "gs = [0.1, 0.5, 0.8, 0.9, 1.0, 1.1, 1.2, 1.5]\n",
    "for g in gs:\n",
    "    model = TFIModel(L=L, J=1., g=g, bc='finite')\n",
    "    psi = run_DMRG(model)\n",
    "    mag_X.append(np.mean(psi.site_expectation_value(model.sigmax)))\n",
    "    mag_Z.append(np.mean(psi.site_expectation_value(model.sigmaz)))\n",
    "    E.append(np.sum(psi.bond_expectation_value(model.H_bonds)))\n",
    "    corr.append(psi.correlation_function(model.sigmaz, model.L//4, model.sigmaz, model.L * 3 // 4))\n",
    "    max_chi.append(max(psi.get_chi()))\n",
    "\n",
    "results_all[L] = {\n",
    "    'g': gs,\n",
    "    'mag_X': mag_X,\n",
    "    'mag_Z': mag_Z,\n",
    "    'E': E,\n",
    "    'corr': corr,\n",
    "    'max_chi': max_chi\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bdc10397-269f-469d-a2c0-1bb0232dfbb0",
   "metadata": {},
   "outputs": [],
   "source": [
    "key = 'corr'\n",
    "\n",
    "plt.figure()\n",
    "for L in results_all:\n",
    "    res_L = results_all[L]\n",
    "    plt.plot(res_L['g'], res_L[key], marker='o', label=\"L={L:d}\".format(L=L))\n",
    "plt.xlabel('g')\n",
    "plt.ylabel(key)\n",
    "plt.legend(loc='best')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6d59b6c5-551f-4b2a-b7cd-7d8c1581de67",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ea8d130a-62e8-4dd6-b3d0-68849714f444",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ba570aab-8ad5-4781-927f-1c9d71906e14",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ab629c07-7576-4ec5-baef-d71f20f211e4",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "b70436f7-f7c9-49a8-bdc4-c90d322681ca",
   "metadata": {},
   "source": [
    "## Infinite DMRG"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "78ae81b0-3ade-4133-9280-dd57b9958807",
   "metadata": {},
   "source": [
    "The given DMRG code also works with `bc='infinite'` boundary conditions of the model and state.\n",
    "The given `SimpleDMRG` code also allows to run infinite DMRG, simply by replacing the `bc='finite'` for both the model and the MPS. \n",
    "\n",
    "- Look at the implementation of `d_dmrg.py` (and `a_mps.py`) to see where the differences are.\n",
    "\n",
    "Again, we can compare to analytic calculations possible for the Transverse Field Ising model.\n",
    "  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "98e49bcb-9dc9-4dc1-aff3-f8d08fff5085",
   "metadata": {},
   "outputs": [],
   "source": [
    "from tenpy_toycodes.tfi_exact import infinite_gs_energy"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f1ec34c0-68ac-4e4a-823c-f8d0ef7e7a63",
   "metadata": {},
   "source": [
    "\n",
    "The `L` parameter now just indices the number of tensors insite the unit cell of the infinite MPS.\n",
    "It has to be at least `2`, since we optimize 2 tensors at once in our DMRG code.\n",
    "Note that we now use the `mean` to calculate densities of observables instead of extensive quantities:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "545b6a75-1550-433e-8763-968281137a41",
   "metadata": {},
   "outputs": [],
   "source": [
    "model = TFIModel(L=2, J=1., g=0.8, bc='infinite')  # just change bc='infinite' here\n",
    "\n",
    "chi_max = 10\n",
    "\n",
    "psi = init_FM_MPS(model.L, model.d, model.bc)\n",
    "eng = SimpleDMRGEngine(psi, model, chi_max=chi_max, eps=1.e-7)\n",
    "for i in range(10):\n",
    "    E_dmrg = eng.sweep()\n",
    "    E = np.mean(psi.bond_expectation_value(model.H_bonds))\n",
    "    #if i % 10 == 9:\n",
    "    print(\"sweep {i:2d}: E/L = {E:.13f}\".format(i=i + 1, E=E))\n",
    "print(\"final bond dimensions: \", psi.get_chi())\n",
    "mag_x = np.mean(psi.site_expectation_value(model.sigmax))\n",
    "mag_z = np.mean(psi.site_expectation_value(model.sigmaz))\n",
    "print(\"magnetization density in X = {mag_x:.5f}\".format(mag_x=mag_x))\n",
    "print(\"magnetization density in Z = {mag_z:.5f}\".format(mag_z=mag_z))\n",
    "E_exact = infinite_gs_energy(model.J, model.g)\n",
    "print(\"err in energy = {err:.3e}\".format(err=E - E_exact))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d2ed4ee1-6b8e-448b-9a76-b2a9f2a1ef6e",
   "metadata": {},
   "source": [
    "### Exercise: Infinite DMRG\n",
    "\n",
    "Try running the infinite DMRG code. How many sweeps do you now need to perform to converge now?\n",
    "  Does it depend on `g` and `chi`?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "80476e8b-e84a-41b3-950e-1b378bd853ae",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ca2aa6c1-19d6-4a9b-b63c-8f5706a9f6e0",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "292897df-7862-46ab-a1ae-1239185d9828",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "6045388b-47b8-4a5e-bc3b-bcd2dfa07be2",
   "metadata": {},
   "source": [
    "From the exercise, you should see that you need significantly more sweeps when the correlation length gets larger, in particular at the critical point! At the critical point, the [correlation length (of infinite MPS) scales](https://arxiv.org/abs/0812.2903) as $\\xi \\propto \\chi^\\kappa$, so you again need more sweeps to converge at larger bond dimensions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "60d545ba-a00d-4461-b6db-82336a02465e",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d16df3aa",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "1e99ab49",
   "metadata": {},
   "source": [
    "### Advanced exercises for the experts (and those who want to become them ;-) )\n",
    "\n",
    "- Obtain the ground state of the transverse field ising model at the critical point with DMRG for large `L`.\n",
    "  Try to plot the corrlation function as a function of `j-i`.\n",
    "  What form does it have? Is an MPS a good ansatz for that?\n",
    "\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8adcb517-0981-4faa-a756-1b792a1ee00f",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "119170d8-d2bf-4d5b-b799-a08574dc7e37",
   "metadata": {},
   "source": [
    "- Compare running DMRG and imaginary time evolution with TEBD from `tenpy_toycodes/c_tebd.py` \n",
    "  for various parameters of `L`, `J`, `g`, and `bc`. Which one is faster? Do they always agree?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5f1aabfc-7a3d-4450-a1fa-d94f8f2d73fe",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1dfbc061",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b9d25343",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "da649492-fe07-437f-931a-7843a3e4a6a4",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "97801804-24e0-4261-912f-699fc8b93228",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1c5c22e8-1254-4073-86d6-8d4c68e1617b",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}