{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from tenpy.models.lattice import Lattice\n",
    "from tenpy.models.model import CouplingMPOModel\n",
    "from tenpy.networks.site import SpinSite\n",
    "from tenpy.tools.params import get_parameter\n",
    "from tenpy.algorithms import dmrg\n",
    "from tenpy.networks.mps import MPS\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import time\n",
    "import random"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Model Definition"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Create a ComplingMPOModel to describe the system.\n",
    "#The parameters are passed via model_params\n",
    "class Heisenberg(CouplingMPOModel):\n",
    "\n",
    "    def __init__(self,model_params):\n",
    "        CouplingMPOModel.__init__(self,model_params)\n",
    "    \n",
    "    def init_lattice(self, model_params):\n",
    "        #Here initialize the type of lattice considered\n",
    "        lattice = get_parameter(model_params, 'lattice', self.name, False)\n",
    "        return lattice\n",
    "        \n",
    "    def init_terms(self, model_params):\n",
    "        D= get_parameter(model_params, 'D', 0., self.name, True)\n",
    "        #Try to get also the J and lambda paramter in the same way\n",
    "        #J= ...\n",
    "        #lam= ...\n",
    "        \n",
    "        #Here implement the couplings of the chain\n",
    "        \n",
    "        #D term\n",
    "        self.add_onsite(D, 0, 'Sz Sz')     \n",
    "        \n",
    "        #Heisenberg interaction\n",
    "        self.add_coupling(J, 0, 'Sx', 0, 'Sx', 1,) \n",
    "        #Implement the remaining couplings for the Heisenberg exchange\n",
    "        #Have a look at the documentation of the function add_coupling()\n",
    "        #\n",
    "        #\n",
    "        \n",
    "        #Quadratic Heisenberg term, the one with lambda\n",
    "        #Implement this term yourself\n",
    "        "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Phase diagram for $\\lambda=0$ and varying $D$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Define the model "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "L=2 #Chain length (2 for the iMPS)\n",
    "J=1 #Anitferromagnetic interaction\n",
    "lam=0 #Quadratic term\n",
    "Dmin=-1\n",
    "Dmax=2\n",
    "points=30\n",
    "d=np.linspace(Dmin,Dmax,points,endpoint=True) #Values considered for D\n",
    "\n",
    "#Definition of the lattice model\n",
    "#Define the local Hilbert space, have a look at the documentation of SpinSite\n",
    "#site=SpinSite(...)\n",
    "#Define the lattice by looking at the documentation of Lattice\n",
    "#Inifinite MPS requires periodic boundaries for the lattice\n",
    "#lat=Lattice(...) #We choose an infinite MPS  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Perform iDMRG"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "entEnt=[]\n",
    "entSp=[]\n",
    "stagMag=[]\n",
    "for D in d:\n",
    "    \n",
    "    #Define the paramters of the model\n",
    "    model_params={\n",
    "        'verbose':0,\n",
    "        'J':J,\n",
    "        'D':D,\n",
    "        'lambda':lam, \n",
    "        'lattice':lat \n",
    "    }\n",
    "\n",
    "    #Create the model\n",
    "    chain=Heisenberg(model_params)\n",
    "\n",
    "    #Define the paramters for the DMRG\n",
    "    dmrg_paramsGs = {\n",
    "        'trunc_params': {\n",
    "            'chi_max': 20 #Maximum bond dimension\n",
    "        },\n",
    "        'verbose': 1\n",
    "    }\n",
    "    \n",
    "    #Initialize the MPS with an arbitrary product state\n",
    "    product_state=['up',0] \n",
    "    #Have a look at the documentation MPS.from_product_state()\n",
    "    #psiGs=MPS.from_product_state(...)\n",
    "    \n",
    "    #Perform iDMRG\n",
    "    #Check documentation of dmrg.run\n",
    "    #info=dmrg.run(...)\n",
    "   \n",
    "    #Check documentation for expectation_value(). Which operator do you want to measure for magnetization?\n",
    "    #mag=psiGs.expectation_value(...) #Measure magnetization\n",
    "    #Implement function to compute staggered magnetization from mag\n",
    "    #staggered=...\n",
    "    stagMag.append(staggered) #Compute staggered magnetization\n",
    "    \n",
    "    #Check documentation of entanglement_spectrum() and compute it\n",
    "    #spectrum=...\n",
    "    entSp.append(spectrum[0])\n",
    "    \n",
    "    #Check documentation of entanglement_entropy() and compute it\n",
    "    #entropy=...\n",
    "    entEnt.append(entropy[0])\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Study staggered magnetization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize=(14,5))\n",
    "plt.plot(d,stagMag,'o',ls='-',c='k')\n",
    "plt.xlabel('D');\n",
    "plt.ylabel('Staggered magnetization');\n",
    "plt.axvspan(-1, -0.33, alpha=0.2, color='red')\n",
    "plt.axvspan(-0.33, 1, alpha=0.2, color='b')\n",
    "plt.axvspan(1, 2, alpha=0.2, color='g')\n",
    "plt.xlim(-1,2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Study of the entanglement entropy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize=(14,5))\n",
    "plt.tight_layout()\n",
    "plt.plot(d,entEnt,'o',ls='-',c='k')\n",
    "plt.xlabel('D');\n",
    "plt.ylabel('Entanglement Entropy');\n",
    "plt.axvspan(-1, -0.33, alpha=0.2, color='red')\n",
    "plt.axvspan(-0.33, 1, alpha=0.2, color='b')\n",
    "plt.axvspan(1, 2, alpha=0.2, color='g')\n",
    "plt.xlim(-1,2);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Study of the entanglement spectrum"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize=(14,5))\n",
    "plt.tight_layout()\n",
    "\n",
    "for i in range(len(d)): \n",
    "    for j in range(int(len(entSp[i])/2)):\n",
    "        if abs(entSp[i][2*j]-entSp[i][2*j+1])<1e-4:\n",
    "            plt.scatter(d[i]-0.01,entSp[i][2*j],c='k',s=50,marker='.')\n",
    "            plt.scatter(d[i]+0.01,entSp[i][2*j+1],c='k',s=50,marker='.')\n",
    "            plt.plot([d[i]-0.03,d[i]+0.03],[entSp[i][2*j],entSp[i][2*j]],c='k')\n",
    "        else:\n",
    "            plt.scatter(d[i],entSp[i][2*j],c='k',s=50,marker='.')\n",
    "            plt.scatter(d[i],entSp[i][2*j+1],c='k',s=50,marker='.')\n",
    "            plt.plot([d[i]-0.03,d[i]+0.03],[entSp[i][2*j],entSp[i][2*j]],c='k')\n",
    "            plt.plot([d[i]-0.03,d[i]+0.03],[entSp[i][2*j+1],entSp[i][2*j+1]],c='k')\n",
    "        \n",
    "plt.xlabel('D');\n",
    "plt.ylabel('Entanglement Spectrum');\n",
    "plt.axvspan(-1.1, -0.33, alpha=0.2, color='red')\n",
    "plt.axvspan(-0.33, 1, alpha=0.2, color='b')\n",
    "plt.axvspan(1, 2.1, alpha=0.2, color='g')\n",
    "plt.xlim(-1.1,2.1);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Phase diagram for $D=0$ and varying $\\lambda$ "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Do the same studies for this new situation."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.7.0"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}