Complete v1.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "002cf83e",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import random as rd\n",
    "from time import time\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "ef333182",
   "metadata": {},
   "outputs": [],
   "source": [
    "class site:\n",
    "    def __init__(self,lattice,x,y,z):\n",
    "        self.x = x\n",
    "        self.y = y\n",
    "        self.z = z\n",
    "        self.lattice = lattice\n",
    "\n",
    "    def setspin(self,s):\n",
    "        self.spin = s\n",
    "    \n",
    "    def setrandomspin(self): # set random spin for this site\n",
    "        self.spin = rd.choice([-1,1]) # set spin = -1 or 1\n",
    "     \n",
    "    def getspin(self):\n",
    "        return self.spin\n",
    "\n",
    "    def flipspin(self):\n",
    "        self.spin *= -1\n",
    "\n",
    "    def getposition(self):\n",
    "        return np.array([self.x,self.y,self.z])\n",
    "\n",
    "    def getzeemanenergy(self,B):\n",
    "        g = 2.0\n",
    "        muB = 5e-2 #mev/T\n",
    "        return -g*muB*self.getspin()*B\n",
    "\n",
    "    def getdeltaexchangeenergy(self,J):\n",
    "        E = 0\n",
    "        NN = self.getnearestneighbourscubic()\n",
    "        for nn in NN:\n",
    "            E += -J*self.spin*nn.spin\n",
    "        return E\n",
    "\n",
    "    def getlistposition(self,x,y,z):\n",
    "        return x*self.lattice.N2*self.lattice.N3 + y*self.lattice.N3 + z\n",
    "\n",
    "    def getnearestneighbourscubic(self):\n",
    "        #i-1 i+1, j-1 j+1, k-1 k+1 with periodic boundary conditions\n",
    "        NN1 = self.lattice.sites[self.getlistposition((self.x-1) % self.lattice.N1, self.y, self.z)]\n",
    "        NN2 = self.lattice.sites[self.getlistposition((self.x+1) % self.lattice.N1, self.y, self.z)]\n",
    "        NN3 = self.lattice.sites[self.getlistposition(self.x, (self.y-1) % self.lattice.N1, self.z)]\n",
    "        NN4 = self.lattice.sites[self.getlistposition(self.x, (self.y+1) % self.lattice.N1, self.z)]\n",
    "        NN5 = self.lattice.sites[self.getlistposition(self.x, self.y, (self.z-1) % self.lattice.N1)]\n",
    "        NN6 = self.lattice.sites[self.getlistposition(self.x, self.y, (self.z+1) % self.lattice.N1)]\n",
    "        return [NN1,NN2,NN3,NN4,NN5,NN6]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "d56f340e",
   "metadata": {},
   "outputs": [],
   "source": [
    "class lattice:\n",
    "    def __init__(self,N1,N2,N3):\n",
    "        self.N1 = N1 # Dimension of the lattice\n",
    "        self.N2 = N2\n",
    "        self.N3 = N3\n",
    "        self.sites = [] # Creating a list to take in the information of sites\n",
    "\n",
    "    def populatesites(self):\n",
    "        for i in range(self.N1):\n",
    "            for j in range(self.N2):\n",
    "                for k in range(self.N3):\n",
    "                    self.sites.append(site(lat,i,j,k)) # Assign each site a coordianate\n",
    "\n",
    "    def initialisespins(self): # Set spin for each site\n",
    "        for i in range(self.N1*self.N2*self.N3):\n",
    "            self.sites[i].setrandomspin()\n",
    "    \n",
    "    def gettotalenergy(self,J):\n",
    "        E_t = 0\n",
    "        for i in range(self.N1*self.N2*self.N3):\n",
    "            E_t = E_t + self.sites[i].getdeltaexchangeenergy(J)\n",
    "        E_t = abs(E_t/2)\n",
    "        return E_t\n",
    "    \n",
    "    def filprandomspins(self): # input number of moves\n",
    "        #self.filp = np.ones(self.N1*self.N2*self.N3)\n",
    "        #self.filp[rd.randrange(self.N1*self.N2*self.N3)] = -1 # choose a random site to flip\n",
    "        site_f = rd.randrange(self.N1*self.N2*self.N3)\n",
    "        self.sites[site_f].flipspin()\n",
    "        return site_f\n",
    "    \n",
    "    def checkfilp(self,E,J,T_pace,T,site_f):\n",
    "        deltaE = self.sites[site_f].getdeltaexchangeenergy(J)\n",
    "        #accept flip if energy change in temperature range:\n",
    "        if deltaE <= 0 or np.exp(-deltaE/(kB*T))>rd.random():\n",
    "            E += deltaE\n",
    "            T = T_pace*T # decrease temperature\n",
    "        else: # flip back\n",
    "            self.sites[site_f].flipspin()\n",
    "        return T,E\n",
    "    \n",
    "    def MC(self,E_t,MCS,J,T_pace,T_start): # Monte carlo model\n",
    "        T_plot = np.zeros(MCS) # setup \n",
    "        E_plot = np.zeros(MCS) # setup\n",
    "        T = T_start\n",
    "        E = E_t\n",
    "        i = 0\n",
    "        #for i in range(MCS):\n",
    "        while MCS > 0:\n",
    "            site_f = self.filprandomspins()\n",
    "            deltaE = self.sites[site_f].getdeltaexchangeenergy(J)\n",
    "            #print(\"deltaE = \", deltaE)\n",
    "            if deltaE <= 0 or np.exp(-deltaE/(kB*T))>rd.random():\n",
    "                E += deltaE\n",
    "                T = T_pace*T # decrease temperature\n",
    "                T_plot[i] = T\n",
    "                E_plot[i] = E\n",
    "                i = i+1\n",
    "                MCS = MCS-1\n",
    "                #print(\"E = \", E)\n",
    "                #print(\"T = \", T)\n",
    "                #print(\"MCS = \", MCS)\n",
    "            else:\n",
    "                self.sites[site_f].flipspin()\n",
    "        \n",
    "        return T_plot,E_plot"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ff112860",
   "metadata": {},
   "source": [
    "- Sites\n",
    "\n",
    "- spins\n",
    "\n",
    "- flip \n",
    "\n",
    "- energy\n",
    "\n",
    "- vary nom\n",
    "\n",
    "- energy\n",
    "\n",
    "- 100"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "357e0b99",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "12400.0\n",
      "Runtime:  0.006695\n"
     ]
    }
   ],
   "source": [
    "starttime = time()\n",
    "\n",
    "# set up some constants\n",
    "g = 2.0\n",
    "muB = 5.788e-2 # meV/T\n",
    "kB = 8.617e-2 # meV/K\n",
    "T_start = 200 # K\n",
    "T_pace = 0.92\n",
    "MCS = 100\n",
    "J = 100 # +ve is FM, -ve is AFM, units meV\n",
    "\n",
    "# set up the sites\n",
    "lat = lattice(10,10,10) # set up the size of lattice\n",
    "lat.populatesites()     # Assign each site a coordinate\n",
    "lat.initialisespins()   # Assign each site a random spin\n",
    "E_t = lat.gettotalenergy(J)    # get the total energy of the system\n",
    "print(E_t)\n",
    "T,E = lat.MC(E_t,MCS,J,T_pace,T_start)\n",
    "# Repeat\n",
    "\n",
    "#lat.fliprandomspin()    # flip random spin\n",
    "#lat.checkflip(T_step,T) # check the flip\n",
    "# Record T and E\n",
    "# stop until T falls beyond certain threshold or MC steps\n",
    "\n",
    "# Repeat\n",
    "\n",
    "print (\"Runtime:  %.6f\"  %(time()-starttime))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "9853db67",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEGCAYAAACKB4k+AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAbiklEQVR4nO3de5hkdX3n8fdnZgIyigv0tDzcuht4ABeJGWRCMIqCSARiJLrGQBqcJWxaEaNc8iSQ3gu59MZ4TR5ZdRshIDM7yHqdNaCyaGDDotjDZWYQhuv0MDDLNDMGiK3AzHz3j3OqONNUV1d3V9WpOufzep56uupXt++ZGfrD7/x+5/dTRGBmZgawIO8CzMysczgUzMysyqFgZmZVDgUzM6tyKJiZWdWivAuYjyVLlsTAwEDeZZiZdZU1a9Y8ExG9tZ7r6lAYGBhgbGws7zLMzLqKpPHpnvPpIzMzq3IomJlZlUPBzMyqHApmZlblUDAzs6pShsLKlSsZGBhgwYIFDAwMsHLlyrxLMjPrCF09JXUuVq5cydDQEJOTkwCMj48zNDQEwODgYJ6lmZnlrnQ9heHh4WogVExOTjI8PJxTRWZmnaN0obBp06ZZtZuZlUnpQqGvr29W7WZmZVK6UBgZGWHx4sW7tS1evJiRkZGcKjIz6xylC4XBwUFGR0fp7+9HEv39/YyOjnqQ2cwMUDfv0bxs2bLwgnhmZrMjaU1ELKv1XOl6CmZmNj2HgpmZVbUsFCRdI2mrpPWZtq9Kuje9bZR0b9o+IOkXmee+1Kq6zMxseq28ovla4ErgK5WGiPj9yn1JnwGezbz+0YhY2sJ6zMxsBi0LhYi4XdJAreckCfgA8I5Wfb+Zmc1eXmMKJwJPR8TDmbZDJd0j6TZJJ073RklDksYkjU1MTLS+UjOzEskrFM4GVmUebwH6IuJY4BLgf0h6ba03RsRoRCyLiGW9vTX3nTYzszlqeyhIWgS8D/hqpS0iXoiIben9NcCjwJHtrs3MrOzy6Cm8E3gwIjZXGiT1SlqY3j8MOAJ4rB3FeG8FM7OXtXJK6irgTuAoSZslnZ8+dRa7nzoCeBuwVtJ9wNeAD0fE9lbVVlHZW2F8fJyIqO6t4GAws7Iq9TIXAwMDjI+Pv6K9v7+fjRs3zqMyM7PO5WUupuG9FczMdlfqUPDeCmZmuyt1KHhvBTOz3ZU6FLy3gpnZ7ko90GxmVkYeaDYzs4Y4FMzMrMqhYGZmVQ4FMzOrciiYmVmVQ8HMzKocCmZmVuVQyPAy2mZWdi3bo7nbVJbRnpycBKguow34CmczKw33FFLDw8PVQKiYnJxkeHg4p4rMzNrPoZDyMtpmZg6FKi+jbWbmUKjyMtpmZq3do/kaSVslrc+0XSHpSUn3prczMs9dLukRSRskvatVdU3Hy2ibmbVw6WxJbwP+FfhKRByTtl0B/GtEfHrKa48GVgHHAwcC/xs4MiJ21vsOL51tZjZ7uSydHRG3A9sbfPmZwA0R8UJEPA48QhIQZmbWRnmMKXxU0tr09NK+adtBwBOZ12xO215B0pCkMUljExMTra7VzKxU2h0KXwQOB5YCW4DPpO2q8dqa57UiYjQilkXEst7e3pYUaWZWVm0NhYh4OiJ2RsQu4CpePkW0GTgk89KDgafaWZuZmbU5FCQdkHn4XqAyM2k1cJakPSUdChwB3NXO2szMrLVTUlcBdwJHSdos6Xzgk5LWSVoLnAxcDBAR9wM3Aj8FvgtcONPMo1bz4nhmVkYtm5LaDq2akjp1cTxILmTzdQtmVgS5TEntZl4cz8zKyqFQgxfHM7OycijU4MXxzKysHAo1eHE8Mysrh0INXhzPzMrKs4/MzErGs4/MzKwhDgUzM6tyKJiZWZVDwczMqhwKDfA6SGZWFovyLqDTTV0HaXx8nKGhIQBPUTWzwnFPYQZeB8nMysShMAOvg2RmZeJQmIHXQTKzMnEozMDrIJlZmTgUZuB1kMysTFq29pGka4B3A1sj4pi07VPA7wAvAo8C50XEv0gaAB4ANqRv/1FEfHim7/DaR2Zms5fX2kfXAqdNabsFOCYi3gg8BFyeee7RiFia3mYMBDMza76WhUJE3A5sn9L2/YjYkT78EXBwq77fzMxmL88xhT8Ebs48PlTSPZJuk3RiXkWZmZVZLlc0SxoGdgCV9SK2AH0RsU3SccC3JL0hIp6r8d4hYAg8LdTMrNna3lOQtJxkAHow0lHuiHghIral99eQDEIfWev9ETEaEcsiYllvb2+7yq7yOkhmVmRt7SlIOg34M+DtETGZae8FtkfETkmHAUcAj7WztkZ4HSQzK7qW9RQkrQLuBI6StFnS+cCVwN7ALZLulfSl9OVvA9ZKug/4GvDhiNhe84Nz5HWQzKzovEfzLCxYsIBaf16S2LVrV9vqMDObj3lfpyCpX9I70/t7Sdq7mQV2C6+DZGZFN2MoSPojklM6/z1tOhj4Vgtr6lheB8nMiq6RnsKFwFuA5wAi4mHgda0sqlN5HSQzK7pGZh+9EBEvSgJA0iKgewci5mlwcNAhYGaF1UhP4TZJfw7sJelU4H8C/6u1ZZmZWR4aCYU/AyaAdcCHgJuA/9jKoszMLB91Q0HSAmBdRFwVEb8XEe9P75f29FGWr242s6KpO6YQEbsk3SepLyK8KXGGr242syKa8eI1ST8Afh24C/h5pT0i3tPa0maW5yY7AwMDjI+Pv6K9v7+fjRs3tr8gM7MG1bt4rZHZR3/R5HoKYdOm2h2n6drNzLrBjKEQEbe1o5Bu09fXV7On4KubzaybNXJF8/OSnktvv5S0U9Ir9jkoG1/dbGZFNGMoRMTeEfHa9PYq4N+RrHZaar662cyKaE6rpEr6UUSc0IJ6ZiXPgWYzs241r4FmSe/LPFwALKPEy1yYmRVZI7OPfidzfwewETizJdWYmVmuGlnm4ssRcV56+6OIGCHZLtMyfHWzmRVBI6Hw+QbbSqtydfP4+DgRUb262cFgZt1m2lCQ9GZJlwK9ki7J3K4AFs70wZKukbRV0vpM236SbpH0cPpz38xzl0t6RNIGSe+a53G1lfduNrOiqNdT2AN4Dcm4w96Z23PA+xv47GuB06a0XQbcGhFHALemj5F0NHAW8Ib0PV+QNGPwdApf3WxmRTHtQHN6JfNtkq6NiFdeujuDiLhd0sCU5jOBk9L71wH/RLI095nADRHxAvC4pEeA44E7Z/u9efDVzWZWFI2MKUxK+pSkmyT9oHKb4/ftHxFbANKflW09DwKeyLxuc9rWFXx1s5kVRSOhsBJ4EDiUZHG8jcBPmlyHarTVvBZC0pCkMUljExMTTS5jbnx1s5kVRSPXKfRExNWSPp45pTTXRfKelnRARGyRdACwNW3fDBySed3BwFO1PiAiRoFRSK5onmMdTee9m82sCBrpKbyU/twi6bclHUvyS3suVgPL0/vLgW9n2s+StKekQ0mug7hrjt+RO1+zYGbdqpGewl9L+jfApSTXJ7wWuHimN0laRTKovETSZuC/AJ8AbpR0PrAJ+D2AiLhf0o3AT0mumr4wInbO/nDy5x3ZzKyb1V0QL50W+rGI+Fz7SmpcJy6I5x3ZzKzT1VsQr+7po/T/1nPfdrOb+JoFM+tmjYwp/F9JV0o6UdKbKreWV9alprs2wdcsmFk3aGRM4TfTn3+ZaQvgHc0vp/uNjIzsNqYAvmbBzLpHI3s0n9yOQoqiMpg8PDzMpk2b6OvrY2RkxIPMZtYVZtx5TdL+wH8FDoyI09N1it4cEVe3o8B6OnGg2cys0815oDl1LfA94MD08UPARU2pzMzMOkojobAkIm4EdgFExA6gK68hyIMvZDOzbtLIQPPPJfWQrkUk6QTg2ZZWVRC+kM3Muk0jYwpvIrmS+RhgPdALvD8i1ra+vPo6fUzBF7KZWSeqN6bQyOyjuyW9HTiKZDXTDRHx0gxvM3whm5l1nxlDQdKrgI8AbyU5hfR/JH0pIn7Z6uK6nTffMbNu08hA81dItsn8PHAlcDRwfSuLKgpvvmNm3aaRgeajIuLXMo9/KOm+VhVUJL6Qzcy6TSM9hXvSGUcASPoN4I7WlVQsg4ODbNy4kV27djEyMsLw8LCnp5pZx2qkp/AbwAclVUZH+4AHJK0DIiLe2LLqCsTTU82sGzQyJbW/3vMR8cqR1Dbp9CmpWZ6eamadYr5TUscl7Uuyh/KiTPvdzSux+Dw91cy6QSNTUv8K+PfAo6RXNeOls2fN01PNrBs0MtD8AeDwiDgpIk5Ob3MOBElHSbo3c3tO0kWSrpD0ZKb9jLl+Ryfy9FQz6waNhMJ6YJ9mfWFEbIiIpRGxFDgOmAS+mT79ucpzEXFTs76zEwwODjI6Okp/fz+S6O/vZ3R01IPMZtZRGgmFvyGZlvo9SasrtyZ9/ynAo3kOVreTp6eaWadrZErqdcDfAutIl89uorOAVZnHH5X0QWAMuDQifjb1DZKGgCHo3vPxnp5qZp2qkSmpt0XE25v+xdIewFPAGyLi6XSHt2dIBrH/CjggIv6w3md005TULE9PNbM8zWtKKrBG0t8Aq4EXKo1NmJJ6OnB3RDydft7TlSckXQV8Z56f37E8PdXMOlUjoXBs+vOETFszpqSeTebUkaQDImJL+vC9JAPcheTpqWbWqRq5eO3kZn+ppMXAqcCHMs2flLSUJHA2TnmuUEZGRnYbUwBPTzWzzjDj7CNJ+0u6WtLN6eOjJZ0/ny+NiMmI6ImIZzNt50bEr0bEGyPiPZleQ+FMnZ7a09PDXnvtxbnnnuuZSGaWq0ampF4LfA84MH38EHBRi+opjcr01Ouvv55f/OIXbNu2jYiozkRyMJhZHqYNBUmVU0tLIuJG0umoEbED2NmG2kpheHh4t9NIAJOTkwwPD+dUkZmVWb2ewl3pz59L6iFd9yjdW+HZad9ls+KZSGbWSeoNNCv9eQnJdNTDJd0B9ALvb3VhZeGZSGbWSer1FHolXQKcRLI20SeBm4GrgHe2vrRy8EJ5ZtZJ6oXCQuA1wN7Aq0l6FQuBxWmbNUGthfKWL1/udZHMLBfTLnMh6e6IeFOb65mVbl3mop6p6yJB0nPwiqpm1iz1lrmo11NQneesRTwbyczyVC8UTmlbFVbl2UhmlqdpQyEitrezEEtMN+vIs5HMrB0auaLZ2qjWbCRJjI+Pe9DZzFrOodBhsrORIAmEymQAL4FhZq024yY7nayIs4+yvBmPmbXCXGcfWc486Gxm7eZQ6GAedDazdnModDAPOptZuzkUOpgHnc2s3TzQ3CU86GxmzVJvoHnGPZpbQdJG4HmSzXp2RMQySfsBXwUGSPZo/kBE/CyP+jqRB53NrB3yPH10ckQszaTVZcCtEXEEcGv62FIedDazduikMYUzgevS+9cBv5tfKZ3Hg85m1g55hUIA35e0RtJQ2rZ/RGwBSH++rtYbJQ1JGpM0NjEx0aZy8+dBZzNrh1wGmiUdGBFPSXodcAvwx8DqiNgn85qfRcS+9T6nTAPNWR50NrP56LgrmiPiqfTnVpKtPo8HnpZ0AED6c2setXUDDzqbWau0PRQkvVrS3pX7wG8B64HVwPL0ZcuBb7e7tm4x3eByRHh8wczmJY+ewv7AP0u6D7gL+MeI+C7wCeBUSQ8Dp6aPrYZag84VHl8ws/nwxWtdauXKlQwPD9ccWwCPL5jZ9OqNKTgUutyCBQuo9XcoiV27duVQkZl1uo4baLbm8fiCmTWTQ6HLeXzBzJrJodDlpl7UNtXk5CTDw8NtrsrMupVDoQAGBwfZuHEjkmo+76UwzKxRDoUCqbc4nk8lmVkjHAoFUm98AXwqycxm5lAokJnGF8BLYZhZfQ6FgqmML0wXDJ6qamb1OBQKylNVzWwuHAoF5amqZjYXDoUC81RVM5sth0IJeKqqmTXKoVACjUxVPeecc9xrMDOHQhk0MlUV3GswMy+dXTrT7e+c5b0YzIrNS2db1UynksAD0GZl5lAoGZ9KMrN62h4Kkg6R9ENJD0i6X9LH0/YrJD0p6d70dka7ayuLylTVFStWeADazHazKIfv3AFcGhF3S9obWCPplvS5z0XEp3OoqZQGBwcB6u71DC/3GrLvMbNiantPISK2RMTd6f3ngQeAg9pdhyVmWiupwr0Gs3LIdUxB0gBwLPDjtOmjktZKukbSvtO8Z0jSmKSxiYmJdpVaeI0MQIPHGsyKLrdQkPQa4OvARRHxHPBF4HBgKbAF+Eyt90XEaEQsi4hlvb297Sq38BodgAb3GsyKLJdQkPQrJIGwMiK+ARART0fEzojYBVwFHJ9HbWXW6AB0hXsNZsWTx+wjAVcDD0TEZzPtB2Re9l5gfbtrs8Rsew3Lly9nwYIF7jmYFUAePYW3AOcC75gy/fSTktZJWgucDFycQ22Wmk2vYefOnUSEew5mBeBlLmxGK1eunHHa6lT9/f2MjIx4CqtZB/IyFzYvsx1rgGS84dxzz0WSTyuZdRGHgjUsO9YgiYULF9Z9faUX6tNKZt3DoWCzUuk17Nq1i+uuu67hnoOnsZp1B4eCzdlsZilV+LSSWWdzKNi8zGW8IXtayQFh1lkcCtYUU3sNyeUoM3NAmHUWh4I1TaXXEBFcf/31szqtBA4Is07gULCWmMtppSwHhFk+HArWUnM9rZTlgDBrH4eCtdx0p5UcEGadx6FgbdXqgFiyZAlLlizxAn1mc+RQsNy0IiC2bdvGtm3bqgv0nXfeedWQcGCYzcyhYB2hmQGR9dJLL1VDYmpguHdh9koOBes4rQqIqabrXTgsrMwcCtbR2hUQWbMJCweHFY1DwbpGHgGRVSssZhscDhHrdA4F60q1AkISPT099PT0VO/vscceba1rpuCYLkQGBgb4yEc+wsDAQN0gccBYq3nnNSu0yq5xmzZtYr/99gOSX9iS6OZ/+7VUjqmnpweA7du3V495vvf7+vo444wzuOmmm3b7s5zrZ3lXvnzV23mNiOioG3AasAF4BLis3muPO+64MJuLFStWRH9/f0iKnp6e6OnpCSAkBeBbi2+VP+fKn33272G+9/v7++OCCy54xd9vp9xvZn39/f2xYsWKWf/7B8Ziut/B0z2Rxw1YCDwKHAbsAdwHHD3d6x0K1mwOC9+67bZ48eJZBwN1QqHTxhSOBx6JiMci4kXgBuDMnGuyEsnuLPfMM8/wzDPPzDhuUTld064Bb7OsyclJhoeHm/Z5nRYKBwFPZB5vTtuqJA1JGpM0NjEx0dbirLxqhcVsg8MhYq2yadOmpn1Wp4VCrf9KYrcHEaMRsSwilvX29rapLLOZzRQc9UKkv7+fCy64oOFQccBYVl9fX9M+a1HTPqk5NgOHZB4fDDyVUy1mLTU4ONi0GTi1Zll14uyjos78ytPixYsZGRlp2ud11JRUSYuAh4BTgCeBnwB/EBH313q9p6SadZ9uCbBW3e+E6b31pqR2VCgASDoD+DuSmUjXRMS0EehQMDObvXqh0Gmnj4iIm4Cb8q7DzKyMOm2g2czMcuRQMDOzKoeCmZlVORTMzKyq42YfzYakCWB8Hh+xBHimSeV0izIeM5TzuH3M5THb4+6PiJpX/3Z1KMyXpLHppmUVVRmPGcp53D7m8mjmcfv0kZmZVTkUzMysquyhMJp3ATko4zFDOY/bx1weTTvuUo8pmJnZ7sreUzAzswyHgpmZVZUyFCSdJmmDpEckXZZ3Pa0g6RBJP5T0gKT7JX08bd9P0i2SHk5/7pt3ra0gaaGkeyR9J31c6OOWtI+kr0l6MP07f3PRjxlA0sXpv+/1klZJelURj1vSNZK2SlqfaZv2OCVdnv5+2yDpXbP5rtKFgqSFwH8DTgeOBs6WdHS+VbXEDuDSiPi3wAnAhelxXgbcGhFHALemj4vo48ADmcdFP+6/B74bEa8Hfo3k2At9zJIOAj4GLIuIY0iW2z+LYh73tcBpU9pqHmf63/lZwBvS93wh/b3XkNKFAnA88EhEPBYRLwI3AGfmXFPTRcSWiLg7vf88yS+Jg0iO9br0ZdcBv5tLgS0k6WDgt4EvZ5oLe9ySXgu8DbgaICJejIh/ocDHnLEI2CvdoGsxyU6NhTvuiLgd2D6lebrjPBO4ISJeiIjHgUdIfu81pIyhcBDwRObx5rStsCQNAMcCPwb2j4gtkAQH8LocS2uVvwP+FNiVaSvycR8GTAD/kJ4y+7KkV1PsYyYingQ+DWwCtgDPRsT3KfhxZ0x3nPP6HVfGUKi1w3lh5+VKeg3wdeCiiHgu73paTdK7ga0RsSbvWtpoEfAm4IsRcSzwc4pxyqSu9Bz6mcChwIHAqyWdk29VHWFev+PKGAqbgUMyjw8m6XIWjqRfIQmElRHxjbT5aUkHpM8fAGzNq74WeQvwHkkbSU4NvkPSCop93JuBzRHx4/Tx10hCosjHDPBO4PGImIiIl4BvAL9J8Y+7YrrjnNfvuDKGwk+AIyQdKmkPkgGZ1TnX1HSSRHKO+YGI+GzmqdXA8vT+cuDb7a6tlSLi8og4OCIGSP5ufxAR51Dg446I/wc8IemotOkU4KcU+JhTm4ATJC1O/72fQjJ2VvTjrpjuOFcDZ0naU9KhwBHAXQ1/akSU7gacATwEPAoM511Pi47xrSRdxrXAventDKCHZKbCw+nP/fKutYV/BicB30nvF/q4gaXAWPr3/S1g36Ifc3rcfwE8CKwHrgf2LOJxA6tIxk1eIukJnF/vOIHh9PfbBuD02XyXl7kwM7OqMp4+MjOzaTgUzMysyqFgZmZVDgUzM6tyKJiZWZVDwUpPUki6PvN4kaSJygqradvpksbSFUgflPTpBj97QNIftKJus1ZwKJgly0IcI2mv9PGpwJOVJyUdA1wJnBPJqrPHAI81+NkDgEPBuoZDwSxxM8nKqgBnk1wsVPGnwEhEPAgQETsi4gtTP0DS2yXdm97ukbQ38AngxLTt4nSfh09J+omktZI+lL73JEm3S/qmpJ9K+pKkBenrr033C1gn6eKW/ilY6S3KuwCzDnED8J/TU0ZvBK4BTkyfOwb4TAOf8SfAhRFxR7oQ4S9JFqb7k4h4N4CkIZLVPH9d0p7AHZK+n77/eJI9PsaB7wLvAx4HDopkvwAk7TPvIzWrwz0FMyAi1pKc6jkbuGmOH3MH8FlJHwP2iYgdNV7zW8AHJd1LspR5D8naNAB3RbLPx06SnspbSU5THSbp85JOAwq/0q3ly6Fg9rLVJOvzr5rSfj9w3ExvjohPAP8B2Av4kaTX13iZgD+OiKXp7dBI9gCAVy5vHBHxM5Kd1P4JuJDdNw4yazqHgtnLrgH+MiLWTWn/FPDnko4ESM/1XzL1zZIOj4h1EfG3JIvTvR54Htg787LvAReky5oj6ch0QxyA49PVexcAvw/8s6QlwIKI+Drwn0iWxDZrGY8pmKUiYjPJXsdT29dKughYJWkxyf/R/2ONj7hI0snATpKlq28m2f1th6T7SPbZ/XuS01R3p8s9T/DyNop3kgxM/ypwO/DN9P4/pEEBcPl8j9OsHq+SatYBJJ1EZkDaLC8+fWRmZlXuKZiZWZV7CmZmVuVQMDOzKoeCmZlVORTMzKzKoWBmZlX/H3En0PzQ5gr+AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "MCS_p = np.arange(0,100)\n",
    "plt.close('all')\n",
    "plt.figure()\n",
    "plt.plot(MCS_p,T,'ok')\n",
    "plt.xlabel('MC steps')\n",
    "plt.ylabel('Temperature')\n",
    "#plt.xlim([0,120])\n",
    "#plt.ylim([-0.25,0.05])\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "14da06ff",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEGCAYAAABPdROvAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAgwElEQVR4nO3de5SdVZnn8e8vSQspFYUksDAhVdGOdge8INV0vHRrizTxMobVozNxypBR2swAixZxlg2eNatXO1OrsbW90Ao9EZBA10AzqE2Ww0U6Oq3jIBhAuSNRqkIkDVWFQpax0VDP/PHugsOpc069p879nN9nrVrnPfu876n9KtTD++xn762IwMzMbKEWtbsDZmbW3RxIzMysLg4kZmZWFwcSMzOriwOJmZnVZUm7O9Bqy5cvj6GhoXZ3w8ysq9x+++1TEbGi3Gd9F0iGhobYtWtXu7thZtZVJE1U+sypLTMzq4sDiZmZ1cWBxMzM6uJAYmZmdXEgMTOzujiQ5DA2NsbQ0BCLFi1i+fLlLF++fM7x0NAQY2Nj7e6qmVnL9V35b63GxsbYunUrBw4cAGB6evrZz4qPJyYm2Lp1KwAjIyOt7aSZWRv5iWQehULh2SAynwMHDlAoFJrcIzOzztK0QCLpMkmPS7qnqO3Tkh6QdJekr0t6adFn50vaLelBSacUtZ8g6e702YWSlNoPkfQPqf1WSUPNuI89e/bUdP7ExIRTXWbWV5r5RHI5sKGk7WbguIh4DfBj4HwASeuATcCx6ZqLJC1O11wMbAXWpp/Z7zwd+HlE/DbwOeBTzbiJ1atX13xNRDyb6nIwMbNe17RAEhHfAZ4oaftmRBxMb78PrErHG4GrI+LpiHgY2A2cKOlo4LCIuCWyrRyvAE4tumZ7Or4WOGn2aaWRRkdHGRgYWNC1TnWZWT9o5xjJh4Ab0vFK4JGiz/amtpXpuLT9edek4PQksKzcL5K0VdIuSbsmJydr6uTIyAjbtm1jcHAQSSxbtoxly5bNOa6k1tSYmVm3aUsgkVQADgKzeZ9yf4mjSnu1a+Y2RmyLiOGIGF6xouzilVWNjIwwPj7OzMwMU1NTTE1NzTkeHBwse21EVCwZdvmwmfWClgcSSVuAdwMjKV0F2ZPGMUWnrQIeTe2ryrQ/7xpJS4CXUJJKa6VqKbDp6Wmmp6eJiIrHHlMxs27V0kAiaQPw58B7IqK4pnYHsClVYq0hG1S/LSL2AfslrU/jH6cB1xVdsyUdvxf4VlFgarniFNhCeUzFzLqRmvW3V9JVwFuB5cBjwF+QVWkdAszO5Pt+RPzndH6BbNzkIHBORNyQ2ofJKsCWko2pnB0RIelQ4ErgeLInkU0R8dP5+jU8PBzN3o9k0aJFLPR/V0nMzMw0uEdmZvWRdHtEDJf9rI3/Ed8WrQgkQ0NDTExU3AOmqsHBQcbHxxvbITOzOlULJJ7Z3gQLLRkeGBhgdHS0CT0yM2seB5ImyFsyXHw8ODjItm3bvE6XmXUdp7bMzGxeTm2ZmVnTOJB0iOI9Tzw50cy6ifcj6QCle554bxMz6yZ+IukA5fY88eREM+sWDiQdoNLCjl7w0cy6gQNJB6i050npgo8eOzGzTuRA0gHyLvjohR3NrBM5kHSAWhZ89NiJmXUaB5IOMbvnSZ5NHicmJpzmMrOO4UDSYfLuEe80l5l1CgeSDlPLgo9Oc5lZJ3Ag6TCVFnysZGJiwlv5mllbedHGLlHPHiezBgYGvMKwmS2IF23sAQvd46SYU2Fm1gwOJF2iEXvCw3MVX2eeeaYXiTSzhnBqqws1Is1VymkvM6vGqa0e04g0VymnvcxsoRxIulCtW/nm5UUizWwhHEi61OxM+JmZGaamppiamqp4nHdcpXSRyFqPPdZi1p88RtIHSjfOaiaPtZj1Jo+R9LnSVNjg4CBnnHFG3RVg5Xisxaz/NC2QSLpM0uOS7ilqO0LSzZIeSq+HF312vqTdkh6UdEpR+wmS7k6fXai0qqGkQyT9Q2q/VdJQs+6lFxSnwsbHx7noootyLxJZq9LZ9k55mfW2Zj6RXA5sKGk7D9gZEWuBnek9ktYBm4Bj0zUXSVqcrrkY2AqsTT+z33k68POI+G3gc8CnmnYnPSzvIpG18j4qZv2jaYEkIr4DPFHSvBHYno63A6cWtV8dEU9HxMPAbuBESUcDh0XELZEN5lxRcs3sd10LnKRm/Od1j2tGKXE5TnmZ9a5Wj5EcFRH7ANLrkal9JfBI0Xl7U9vKdFza/rxrIuIg8CRQttZV0lZJuyTtmpycbNCt9IZaS4nrKTGutsCk019m3WtJuzuQlHuSiCrt1a6Z2xixDdgGWdXWQjrYy0ZGRuqusso72356errs8Wz6a7Y/ZtY9Wv1E8lhKV5FeH0/te4Fjis5bBTya2leVaX/eNZKWAC9hbirNWsSLSpr1r1YHkh3AlnS8BbiuqH1TqsRaQzaofltKf+2XtD6Nf5xWcs3sd70X+Fb026SYDlLrPiqVeHa9WfdpZvnvVcAtwKsk7ZV0OnABcLKkh4CT03si4l7gGuA+4EbgrIh4Jn3VGcAlZAPwPwFuSO2XAssk7QbOJVWAWfuUm21f61yVZlWRmVnzeGa7NVUts+o9K96sc3lmu7VN3qqwwcFBBxGzLuVAYk2XZ4HJ0dFRCoXCs6XA3njLrHs4tWVtlyf95bSXWXs5tWUdrVAozDuG4tJgs87lQGJtl7fkt9rMeM+SN2ufTpnZbn1s9erVufegrzQz3rPkzdrHTyTWdt6D3qy7OZBY2zVr4y3PkjdrDVdtWUfLuxhkOYODg4yPjze2Q2Z9ylVb1rUWmvYaGBhgdHS0CT0ys1IOJNbRFrpfytKlS9m8eXOuKi9XfJnVx6kt6ym1rO1ViSc/ms3l1Jb1jTyTG+fjii+z2jiQWE9pVKWWK77M8nMgsZ7SqP1MIsLjJWY5OZBYT2nk5MbZGfIOJmbVOZBYT1lolVelbYE9XmI2P1dtmSWLFi2i0r8Pg4ODjI6OupLL+partsxyqDa+4jSXWWUOJGbJfOMrTnOZledAYpYUj69UMjEx4WousxIOJGZFZveXny+YOM1l9hwHErMynOYyy68tgUTSRyXdK+keSVdJOlTSEZJulvRQej286PzzJe2W9KCkU4raT5B0d/rsQklqx/1Y78mT5vLsd7NMywOJpJXAnwHDEXEcsBjYBJwH7IyItcDO9B5J69LnxwIbgIskLU5fdzGwFVibfja08Fasx82X5mrULHqzbteu1NYSYKmkJcAA8CiwEdiePt8OnJqONwJXR8TTEfEwsBs4UdLRwGERcUtkxf9XFF1j1jDl0lze78TsOS0PJBHxM+AzwB5gH/BkRHwTOCoi9qVz9gFHpktWAo8UfcXe1LYyHZe2zyFpq6RdknZNTk428nasD5TbCtjLzJs9Z0mekyQNAmsj4p8kLQWWRMT+hfzCNPaxEVgD/AL4X5I+UO2SMm1RpX1uY8Q2YBtkM9tr6a8ZZMHEgcOsvHmfSCR9GLgW+B+paRXwj3X8zrcDD0fEZET8Bvga8EbgsZSuIr0+ns7fCxxTdP0qslTY3nRc2m5mZi2UJ7V1FvAm4CmAiHiI59JOC7EHWC9pIFVZnQTcD+wAtqRztgDXpeMdwCZJh0haQzaofltKf+2XtD59z2lF15iZWYvkCSRPR8SvZ9+kAfIFp4ci4layJ5w7gLtTH7YBFwAnS3oIODm9JyLuBa4B7gNuBM6KiGfS150BXEI2AP8T4IaF9sssr7GxMYaGhrzHu1ky7+q/kv6abCzjNOBs4EzgvojoytlYXv3X6lFuT3jv8W79oNrqv3kCiYA/Bf6YbID7JuCS6NL15x1IrB5DQ0NMTEzMaR8cHGR8fLz1HTJrkQUvIy9pEXB3RHw5It4XEe9Nx10ZRMzqVWk2+8TEBIsWLWL58uUsX758zrFTYNbLqgaSiJgBfiTJU3jNqD6bPSKYnp5menp6zrEXerRelmew/WjgXkk7Je2Y/Wl2x8w6UT17wnuhR+tVeSYk/mXTe2HWJWYH1AuFAnv27Km4NW8lXujRepH3bDerQ6XB90o8KG/dqq492yXtl/RU+vlXSc9Ieqrx3TTrPrWkurzQo/WqeQNJRLw4Ig5LP4cC/xb4YvO7Ztb5Shd0XLZsGcuWLSt7vHTpUjZv3uzKLus5C0ptSfp+RKxvQn+azqkta7Vykxgr8eRG61T1Tkj8k6K3i4Bh4C0R8YbGdbF1HEis1TyOYr2gWiDJU7X1b4qODwLjZMvAm1kOtVZqubLLuk2eeSSXRMQH08+HI2KUbAVeM8uh1i15I8LjJdZV8gSSv83ZZmZlLGQSo2fCWzepmNqS9AayDadWSDq36KPDgMXN7phZryidxHjEEUcA8MQTTzx7PD09Pee62ZnwHni3TlftieQFwIvIgs2Li36eAt7b/K6Z9Y6RkRHGx8eZmZlhamqKqamp5x1ni2zPNTExUbFcOM+xU2TWCnmqtgYjIn/JSYdz1ZZ1oloru2rhkmJrhLpmtgMHJH1a0vWSvjX70+A+mvW1ehaDnI8Xi7RmyxNIxoAHgDVkCziOAz9oYp/M+k7xDPlmKLdfitNe1ih5Ulu3R8QJku6KiNektn+OiLe0pIcN5tSWdbpmprlKOe1ledWb2vpNet0n6V2SjgdWNax3ZvY8zUxzlXLayxohz8z2/y7pJcDHyOaPHAZ8tKm9MutjecqF8x7nWUvPM+mtXlUDiaTFwNqI+AbwJPBHLemVWZ8bGRlpSLopT5qs1pn3ZqXm27P9GeA9LeqLmTXYfGky75FijZBnjOT/SfqipD+Q9PrZn6b3zMzqVm2/lMHBQQ+0W0PkCSRvBI4FPgn8Tfr5TD2/VNJLJV0r6QFJ90t6g6QjJN0s6aH0enjR+edL2i3pQUmnFLWfIOnu9NmFqjQ92KyPVZpVPzo6SqFQcCmw1W3ewfaIaMa4yBeAGyPivZJeAAwAnwB2RsQFks4DzgP+XNI6YBNZMHsZ8E+SXpnSbhcDW4HvA9cDG4AbmtBfs55SutnW7CKRgJ9QrGZ59mw/StKlkm5I79dJOn2hv1DSYcAfApcCRMSvI+IXZHucbE+nbQdOTccbgasj4umIeBjYDZwo6WjgsIi4JbLSlCuKrjGzKgqFwpwdG10KbAuVJ7V1OXAT2dMAwI+Bc+r4nS8HJoGvSLpT0iWSXggcFRH7ANLrken8lcAjRdfvTW0r03Fp+xyStkraJWnX5ORkHV036w2VSn5dCmwLkSeQLI+Ia4AZgIg4CDxTx+9cArweuDgijgd+SZbGqqTcuEdUaZ/bGLEtIoYjYnjFihW19tes51Qq+XUpsC1EnkDyS0nLSH+kJa0nm1OyUHuBvRFxa3p/LVlgeSylq0ivjxedf0zR9auAR1P7qjLtZjaPcmXBLgW2hcoTSM4FdgCvkPQ9srGIsxf6CyPiX4BHJL0qNZ0E3Jd+x5bUtgW4Lh3vADZJOkTSGrJtfm9L6a/9ktanaq3Tiq4xsyrKlQUvXbqUzZs3e2FHq9m8izYCSFoCvIosnfRgRPxmnkvm+77XAZeQbZ71U+CDZEHtGmA1sAd4X0Q8kc4vAB8CDgLnRMTswP8w2RjOUrJqrbNjnhvyoo1mz1dawVXKCzsaVF+0Mc/qv4cCZwJvJktvfRf4u4j410Z3tBUcSMyeL88yKoODg4yPj7emQ9aRqgWSPIs2XgHsJ1uwEeD9wJXA+xrTPTNrpzyVWq7msmryBJJXRcRri95/W9KPmtUhM2ut1atXe2FHq0uewfY7U6UWAJJ+H/he87pkZq3khR2tXnkCye+TLdw4LmkcuAV4S1rj6q6m9s7Mms4LO1q98gy2V91EOiJasydog3iw3Sy/sbGxshtsrV69mtHRUQeYPlLXYHtETKSVeI8pPj8i7mhcF82s05SWBU9PTz/7mRd5tGJ5nkj+G/AfgZ/w3BIkERFva27XmsNPJGb5uCzYitVb/vvvgFdExK8b2y0z62QuC7a88gy23wO8tMn9MLMOk6fk12XBBvkCyV+RlQDfJGnH7E+zO2Zm7eWyYMsrTyDZDnwKuIDnttr9m2Z2yszar1pZcKVFHr3gY3/KM9j+zxHxlhb1p+k82G5Wn/kWeSzmBR97R72LNn4WeJpsOfenZ9u7tfzXgcSsPnmquYq5sqs31Fu1dXx6XV/UFkBXlv+aWX1qrdRyZVfvyzMh8Y9a0REz6w55FnksPd9627yD7ZKOknSppNnNpNZJOr35XTOzTjRfNVcxV3b1hzxVW5cDNwEvS+9/DJzTpP6YWYebr5rLCz72n4qBJG2vC7A8Iq4BZgAi4iDwTAv6ZmYdamRkhPHxcWZmZpiammJqamrO8ejoKIVCwWXBfaDaGMltwOuBX0paRlpnK+1N8mQL+mZmXcoLPvaXaqktpddzyUp/XyHpe2Rb757d7I6ZWfcqFApV55kcOHCAQqHQwh5ZM1ULJCsknQu8Ffg68NfADcCXgbc3v2tm1q3ylPxOTEw45dUjqqW2FgMv4rknk1n5yjXMrG/lLRF2yqs3VJzZLumOiHh9i/vTdJ7ZbtZ8tSyjUsoz4TtTtZntecZImkLSYkl3SvpGen+EpJslPZReDy8693xJuyU9KOmUovYT0t7xuyVdKKmpfTazfCqVCOcxMTHhVFeXqRZITmry7/4IcH/R+/OAnRGxFtiZ3iNpHbAJOBbYAFwkaXG65mJgK7A2/Wxocp/NLKdyJcKDg4O5ro2IZ1NdDiadr2IgiYgnmvVLJa0C3gVcUtS8kWzJetLrqUXtV0fE0xHxMLAbOFHS0cBhEXFLZPm5K4quMbMOVMuseHB1V7fIM7O9GT4PfJw0yTE5KiL2AaTXI1P7SuCRovP2praV6bi0fQ5JWyXtkrRrcnKyITdgZrWrNiu+Ei/62PlaHkgkvRt4PCJuz3tJmbao0j63MWJbRAxHxPCKFSty/loza4ZKs+Irpb286GPna8cTyZuA90gaB64G3ibp74HHUrqK9Pp4On8vcEzR9auAR1P7qjLtZtaFyqW9vOhjd2h5IImI8yNiVUQMkQ2ifysiPkA2e35LOm0LcF063gFsknSIpDVkg+q3pfTXfknrU7XWaUXXmFmXKU17edHH7tGuMZJyLgBOlvQQcHJ6T0TcC1wD3AfcCJwVEbOLRp5BNmC/G/gJ2cx7M+tSxWmv4kUfXQrc2ebdarfXeEKiWecrN6HR+7+310InJJqZtUW5RR9dCty5HEjMrONUKvn1Qo+dyYHEzDpOtZLf6elppqenPfu9gziQmFnHqWUGvFNe7edAYmYdp7gUOI+JiQmnudrIgcTMOtJsKXAtwcRprvZwIDGzjuY0V+dzIDGzjlbr3iallV2u8mo+T0g0s640NDSUazvfUp7YuDCekGhmPafWvU1mOf3VeA4kZtaVaq3sKuY9ThrLgcTMulatlV2zIsLjJQ3kQGJmXW8haS6XCzeOA4mZdb1qW/hWq/LyeEljOJCYWU+otIXv7HGlfeE9K75+DiRm1heqLQTpNFd9HEjMrC/MN47iNNfCOZCYWV/IUy7ssuCFcSAxs74xX7lwtfSXVeZAYmZ9p1yaa2BggNHR0Tb1qLs5kJhZ3ylXLrx06VI2b97sBR8XwIHEzPrSbJrryiuv5Fe/+tWz2/cWb+XrbX3zcSAxs75WKBQ4cOBArnNd2VVeywOJpGMkfVvS/ZLulfSR1H6EpJslPZReDy+65nxJuyU9KOmUovYTJN2dPrtQlWYcmZlVUGulliu75mrHE8lB4GMR8bvAeuAsSeuA84CdEbEW2Jnekz7bBBwLbAAukrQ4fdfFwFZgbfrZ0MobMbPuV2ul1qJFi+aMn/T7WErLA0lE7IuIO9LxfuB+YCWwEdieTtsOnJqONwJXR8TTEfEwsBs4UdLRwGERcUtku3NdUXSNmVkutS74+Mwzz3gspURbx0gkDQHHA7cCR0XEPsiCDXBkOm0l8EjRZXtT28p0XNpuZpZbngUfJbF48eL5vyzpt7GUtgUSSS8CvgqcExFPVTu1TFtUaS/3u7ZK2iVp1+TkZO2dNbOeNt+CjzMzM8zMzNT0nf20GGRbAomk3yILImMR8bXU/FhKV5FeH0/te4Fjii5fBTya2leVaZ8jIrZFxHBEDK9YsaJxN2JmfWMhs977Jc3VjqotAZcC90fEZ4s+2gFsScdbgOuK2jdJOkTSGrJB9dtS+mu/pPXpO08rusbMrKG8R3xl7XgieROwGXibpB+mn3cCFwAnS3oIODm9JyLuBa4B7gNuBM6KiGfSd50BXEI2AP8T4IaW3omZ9Y2Fbp4F2ZNJL1d5KSt46h/Dw8Oxa9eudnfDzHrQ0NAQExMTNV83MDDAtm3bGBkZaUKvGkPS7RExXO4zz2w3M2uQfk1/OZCYmTVInj1PKunmKi8HEjOzBppvz5NqurXKy4HEzKwJ+inN5UBiZtYEjazy6vSUl6u2zMzapJYqr3ZXdrlqy8ysA9WS/urklJcDiZlZm9Ra5dWpe6E4kJiZtVEtVV4LWe+rFRxIzMw6wHxproGBAUZHR1vYo/wcSMzMOkC1Kq/BwcGOXkLFgcTMrENU2hdldHSUQqHQsaXAS9rdATMzq2xsbIytW7dy4MAB4LnZ70DHPKH4icTMrIMVCoVng8isTisFdiAxM+tglUp+y81+P/PMMxkaGmp5Cswz283MOthC9ziBxs6G98x2M7MutdDFH6F1KTAHEjOzDlbPHifQmtnwDiRmZh2unj1OWjEb3oHEzKxL1JrmatVseAcSM7MuMd/s9zPOOOPZz1o5G95VW2ZmPWhsbIxCocCePXtYvXo1o6OjdQWValVbntluZtZjWj0b3qktM7Me0+rZ8F0fSCRtkPSgpN2Szmt3f8zM2q3abPhmzHjv6kAiaTHwJeAdwDrg/ZLWtbdXZmbtVa3kdzbN1chg0tWBBDgR2B0RP42IXwNXAxvb3Cczs7aar0y40Wmubg8kK4FHit7vTW3PI2mrpF2Sdk1OTrasc2Zm7ZBnNnwjZ7x3eyBRmbY59cwRsS0ihiNieMWKFS3olplZe803G76RM967PZDsBY4per8KeLRNfTEz6zjl0lyNnvHe7YHkB8BaSWskvQDYBOxoc5/MzDpG6Wz4Zsx47/qZ7ZLeCXweWAxcFhFVw6xntpuZ1a6nZ7ZHxPXA9e3uh5lZv+r21JaZmbWZA4mZmdXFgcTMzOriQGJmZnXp+qqtWkmaBCYWePlyYKqB3ekW/Xjf/XjP0J/33Y/3DLXf92BElJ3R3XeBpB6SdlUqf+tl/Xjf/XjP0J/33Y/3DI29b6e2zMysLg4kZmZWFweS2mxrdwfapB/vux/vGfrzvvvxnqGB9+0xEjMzq4ufSMzMrC4OJGZmVhcHkpwkbZD0oKTdks5rd3+aQdIxkr4t6X5J90r6SGo/QtLNkh5Kr4e3u6+NJmmxpDslfSO974d7fqmkayU9kP4/f0Ov37ekj6Z/tu+RdJWkQ3vxniVdJulxSfcUtVW8T0nnp79tD0o6pdbf50CSg6TFwJeAdwDrgPdLWtfeXjXFQeBjEfG7wHrgrHSf5wE7I2ItsDO97zUfAe4vet8P9/wF4MaI+B3gtWT337P3LWkl8GfAcEQcR7b1xCZ6854vBzaUtJW9z/Tv+Cbg2HTNRelvXm4OJPmcCOyOiJ9GxK+Bq4GNbe5Tw0XEvoi4Ix3vJ/vDspLsXren07YDp7alg00iaRXwLuCSouZev+fDgD8ELgWIiF9HxC/o8fsm2zpjqaQlwADZjqo9d88R8R3giZLmSve5Ebg6Ip6OiIeB3WR/83JzIMlnJfBI0fu9qa1nSRoCjgduBY6KiH2QBRvgyDZ2rRk+D3wcmClq6/V7fjkwCXwlpfQukfRCevi+I+JnwGeAPcA+4MmI+CY9fM8lKt1n3X/fHEjyUZm2nq2blvQi4KvAORHxVLv700yS3g08HhG3t7svLbYEeD1wcUQcD/yS3kjpVJTGBDYCa4CXAS+U9IH29qoj1P33zYEkn73AMUXvV5E9EvccSb9FFkTGIuJrqfkxSUenz48GHm9X/5rgTcB7JI2TpSzfJunv6e17huyf6b0RcWt6fy1ZYOnl+3478HBETEbEb4CvAW+kt++5WKX7rPvvmwNJPj8A1kpaI+kFZANTO9rcp4aTJLKc+f0R8dmij3YAW9LxFuC6VvetWSLi/IhYFRFDZP+/fisiPkAP3zNARPwL8IikV6Wmk4D76O373gOslzSQ/lk/iWwcsJfvuVil+9wBbJJ0iKQ1wFrgtlq+2DPbc5L0TrJc+mLgsogYbW+PGk/Sm4HvAnfz3HjBJ8jGSa4BVpP9y/i+iCgdyOt6kt4K/JeIeLekZfT4PUt6HVmBwQuAnwIfJPuPy569b0l/Cfx7sgrFO4E/BV5Ej92zpKuAt5ItFf8Y8BfAP1LhPiUVgA+R/e9yTkTcUNPvcyAxM7N6OLVlZmZ1cSAxM7O6OJCYmVldHEjMzKwuDiRmZlYXBxKzGkkKSVcWvV8iaXJ25eDU9g5Ju9Kqug9I+kyZ7zlK0jck/UjSfZKuT+1Dkv5Da+7GrH4OJGa1+yVwnKSl6f3JwM9mP5R0HPBF4ANpJeXjyOZplPokcHNEvDYiZldZBhgCHEisaziQmC3MDWQrBgO8H7iq6LOPA6MR8QBARByMiIvKfMfRZMtTkM67Kx1eAPyBpB+m/TMWS/q0pB9IukvSf4JsAqWk70j6enqi+TtJi9L5l6c9N+6W9NEG37vZ8ziQmC3M1WTLShwKvIZs9v+s44A8i0B+CbhU2WZiBUkvS+3nAd+NiNdFxOeA08lWqv094PeAD6elLCBb7vtjwKuBVwB/ArwOWBkRx0XEq4Gv1HOjZvNxIDFbgPT0MET2NHL9Ar/jJrLl3L8M/A5wp6QVZU79Y+A0ST8kC1jLyNZDArgt7ZPzDNlT0ZvJ0mgvl/S3kjYAPb2Cs7WfA4nZwu0g29/iqpL2e4ET8nxBRDwREf8zIjaTLQ76h2VOE3B2ekJ5XUSsSftowNzlviMifk624+H/Ac7i+Rt2mTWcA4nZwl0GfDIi7i5p/zTwCUmvBEjjFueWXizpbZIG0vGLyVJTe4D9wIuLTr0JOCMt8Y+kV6ZNqABOTKtSLyJbjPD/SloOLIqIrwL/lWx5eLOmWdLuDph1q4jYS7bveWn7XZLOAa5KgSKA/13mK04AvijpINl/1F0SET9IAeOgpB+R7b39BbI02h1p+fNJntsm9RaywflXA98Bvp6Ov5KCC8D5dd+sWRVe/desSxUve9/mrlifc2rLzMzq4icSMzOri59IzMysLg4kZmZWFwcSMzOriwOJmZnVxYHEzMzq8v8B+kRjT/hFXyoAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "MCS_p = np.arange(0,100)\n",
    "plt.close('all')\n",
    "plt.figure()\n",
    "plt.plot(MCS_p,E,'ok')\n",
    "plt.xlabel('MC Steps')\n",
    "plt.ylabel('Energy')\n",
    "#plt.xlim([0,120])\n",
    "#plt.ylim([-0.25,0.05])\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1d00bdc4",
   "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.9.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}