{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# 7-MonteCarlo Sampling" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "view-in-github" }, "source": [ "[![Open in Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/github/OpenJij/OpenJijTutorial/blob/master/source/en/007-MonteCarloSampling.ipynb)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "OpenJij implements Simulated Annealing (SA). If we keep tempereture constant, it is possible to sample spin sequences from the canonical distribution at this temperature.\n", "\n", "$$\n", "p(\\{\\sigma\\}) = \\frac{\\exp(-\\beta E(\\{\\sigma\\}))}{Z}, \\ Z = \\sum_{\\{\\sigma\\}}\\exp(-\\beta E(\\{\\sigma\\}))\n", "$$\n", "\n", "In the following, we deal with the fully-coupled ferromagnetic ising model.\n", "\n", "$$\n", "E(\\{\\sigma\\}) = \\frac{J}{N} \\sum_{i" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# visualize result\n", "plt.errorbar(temp_list, mag, yerr=mag_std)\n", "plt.plot(temp_list, mag)\n", "plt.xlabel('temperature', fontsize=15)\n", "plt.ylabel(r'$|m|$', fontsize=15)\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This phenomenon, in which the value changes significantly from one temperature to another, is called a **phase transition**. In this model (when the system size is close to infinity), a phase transition occurs at a temperature of 1.0. It has been theoretically proven. \n", "However, it is often not possible to calculate the temperature at which the phase transition occurs theoretically in actual models. For this reason, MonteCarlo simulations are often used to study the properties of phase transitions numerically." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Binder cumulant\n", "\n", "Now, let us assume that we don't know the temperature of the phase transition, try to find the temperature of the phase transition as accurately as possible from the numerical simulations. \n", "Looking at the figure above, we can see the magnetization approach 0 as the temperature increases. However, it is not clear which temperature is the phase transition point. Phase transitions are theoretically phenomena that occurs in systems of infinite size, but simulations can only deal with finite size, which results in an error from the theory. This is called the **finite size effect**. Numerical analysis of system of infinite size is a seemingly impossible. However, in numerical simulation in statistical mechanics, the method of obtaining information of infinite system size from finite system size has been developed. The one of those methods it to use a quantity called **Binder cumulant**\n", "\n", "$$U_4 \\equiv \\frac{1}{2}\\left( 3- \\frac{\\langle m^4\\rangle}{\\langle m^2\\rangle^2} \\right)$$" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "# calculation of U4\n", "def u_4(states):\n", " m = np.array([np.mean(state) for state in states])\n", " return 0.5 * (3-np.mean(m**4)/(np.mean(m**2)**2))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We defer to the statistical mechanics textbook for details. This quantity is close to 1 for ferromagnetism, where the magnetization approaches 1, and 0 for paramagnetism, where the magnetization approaches 0. Furthermore, the phase transition point is known to take a value independent of the system size. Therefore, we can perform the numerical experiments as described above for several system sizes and find that the point where the graph of $U_4$ intersects at a single point is the phase transition point. " ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "# set a list of system size\n", "n_list = [40, 80, 120, 160]\n", "# set a list of temperature\n", "temp_list = np.linspace(0.5, 1.5, 30)\n", "\n", "# set sampler\n", "sampler = oj.SASampler(num_reads=300)\n", "\n", "u4_list_n = []\n", "for n in n_list:\n", " # make instance\n", " h, J = fully_connected(n)\n", " u4_temp = []\n", " for temp in temp_list:\n", " beta = 1.0/temp\n", " schedule = [[beta, 100 if temp < 0.9 else 300]]\n", " response = sampler.sample_ising(h, J, \n", " schedule=schedule, reinitialize_state=False,\n", " num_reads=100 if temp < 0.9 else 1000\n", " )\n", " u4_temp.append(u_4(response.states))\n", " u4_list_n.append(u4_temp)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAELCAYAAADHksFtAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAABewklEQVR4nO3dd1hURxfA4d9so/ciCHZFsYu9xB5j7xpNTDGWLzGaYkxvJiYmGluMmsQ0NTEx9l5ijbH33hVUUEBRkA67O98fiwoKCAos4LzPsw/s3bl3zyVmz96ZuWeElBJFURRFyYrG2gEoiqIohZtKFIqiKEq2VKJQFEVRsqUShaIoipItlSgURVGUbOmsHUB+8PT0lGXLlrV2GIqiKEXK/v37r0spve7dXiwTRdmyZdm3b5+1w1AURSlShBAXM9uuup4URVGUbKlEoSiKomRLJQpFURQlW8VyjEJRFOW21NRUQkNDSUpKsnYohYatrS3+/v7o9foctVeJQlGUYi00NBQnJyfKli2LEMLa4VidlJKoqChCQ0MpV65cjvaxateTEOJXIUSkEOJYFq8LIcRUIcQ5IcQRIURQQceoKErRlpSUhIeHh0oSaYQQeHh45OoKy9pjFLOA9tm83gGolPYYCnxfADEpilLMqCSRUW7/HlbtepJSbhVClM2mSTdgjrTUQt8lhHAVQvhKKa/mRzx7p75LakQEwmCHxmCL5s5PA1q9Aa1Oj0ZvQKu3QavXI/QGhEGP0OkRNgY0egMagwH0+ju/C4MBvY0dDrYu6G3t1D9YRVGKnMI+RuEHXE73PDRt232JQggxFMtVB6VLl36oN4tavoIyoQ9en8OU9ngYRg2YtWDSajBrBWadBqnVIHVapE6HWa/FrNci9VrMeh3SoEPq9Xd+kvZT7+CAg5s3Tq7euLqXxNHNG62TExoHRzSODmidnBC2tioxKUoxcOnSJapWrcro0aMZNWoUAGvXruX111/HZDIxePBg3nvvvXx7/8KeKHJMSjkTmAlQr169h1qNaXufYSy9FozeFI/BlICNTMKWJGxkMgaZgl6mYCNT0JGKzmxCYzajMUt0UqI1g9Zs+akxSzRmiTCBxiSRJhMmMxilwGQWmKXAZAYptZjRIs1mMJsQphS0JonWBLpkic4IOqNEZwKDUaI3gt4EhlTQmdPOG7iZ9rjvb6LRIB1sEY6O6Jyc0Du7onNyRuPkiNbREY2jk+V3Jyc0jk5o3VzReXmh8/JC6+qqkoyiFBIjR46kQ4cOd56bTCZeffVV1q9fj7+/P/Xr16dr165UrVo1X96/sCeKMKBUuuf+advyxZj/Db/zu9ksSTaaSUo1kZhqIinVRFKqmSSj5feEZBPX4pKJvJVMZGwSkbHJRMYmc+2W5Xej+W6u0mCmrqeJ5r4m6nmmEOiYgKvpBsRGQFz43Z/x18GYDDLz6xUzlisZoxDESA0RUsc17LhpcCMGO+KMWhJSzCQnpWBKSMQ2WWKfnJj2uIZ9pMQpVItDssA+WWKTZEl0mdLr0Xl6WhLH7Z+3HyW8sa1cGZ2Pj0omipIDISEhdOjQgWbNmrFjxw78/PxYtmwZdnZ2D9x36dKllCtXDgcHhzvb9uzZQ8WKFSlfvjwA/fr1Y9myZY9tolgODBdCzAMaAjH5NT5xL41GYGfQYmfQ4pbLfc1myc2EFCJjk4m4lcTJq7HsDo7ixzM3mXjUCDhT2r08Dcq506CCOw3LuVPa3f7uh67ZBKZUMKWk+5mCxpSKxpSC3pSMXfx1fG6GQIZHMMg40IPRCa5rtYQ7eRPl4ssNRw+u2zpxxmBDlJDcMCUTlRRFXOwNjLdisE8G5wRwjZe4xoFvog6f5CQ8E67gfDYU+/3J6GLiM5yn1s0N26pVLY9qlp/6UqVU8lAKrc9WHOfElVt5esyqJZ35tEu1B7Y7e/Ysf/31Fz/99BN9+/Zl0aJFXL16lblz597Xtnnz5kydOpW4uDjGjRvH+vXrmTBhwp3Xw8LCKFXq7ndof39/du/enTcnlAmrJgohxF9AS8BTCBEKfAroAaSUPwCrgY7AOSABGGidSHNHoxF4ONrg4WhDoK8zLSt780rLCpjMkpNXb7E7+AZ7gqPYeDKChftDAfBxtiWojCtONnr0OoFeq0l73P3doNWg19qh1zlQxacMQfXbZvxQlhISouBmCLqbIfjcDMbnRjBEnYML+yHxRrog9eBeDjwqkRpQnhsuJYlwLUmY3kBY/BXC4sI4GRtGWFwYV+KvYDQb0Zq0uMRDiVgtLRNLUzvaBa/L14ifNQtSUy2HdXLCNjDwTvKwb9gQvbd3Af71FaVwKleuHLVr1wagbt26hISE8NFHH/H2229nuc/o0aN58803cXR0LKAoM2ftWU/9H/C6BF4toHDynVYjqO7nQnU/FwY1K4fZLDl3LS4tcdzgaGg0SalmUk23H5JUkzlDN1Z6VX2deaFJGbrV9sNWrwUhwMHT8vCvd/8OCTfg+lmIOpv28xxcP4v+7D+UMKdSAqjpUgoqtYOA9lDvCdDbYTKbuJZ4jbA4S+I4c+MMKy5t5Pu4w2jqa2jkUY+O1CAo2hXduUsknTjBzb/+QiYng16PS6dOuA8ciG3lgPz9AyvKA+Tkm39+sbGxufO7VqslMTGRb775Jtsrit27d7Nw4ULeeecdoqOj0Wg02NraUrduXS5fvjvPJzQ0FD8/v3yLXVg+i4uXevXqyeJUZtxslqSa0xKH0Uyy0cymU5HM3hHC6YhYXO31PF2/FAMalqGUu33u38BkhOiLELINzv4D5zdDajzo7KB8i7TE8RS4+N/ZRUrJ6Zun+SfkH9ZfXE/IrRAEgjredXiyzJO08WuF29VYohcuInrRImRiIg7NmuEx6CXsGzVS3VNKgTl58iSBgYFWjSEkJITOnTtz7Jjl3uIJEyYQFxfH6NGjc3yM0aNH4+joyKhRozAajQQEBLBx40b8/PyoX78+f/75J9Wq5TwRZvZ3EULsl1Le9y2zsI9RKFi6smw0Wmx0QNqXkmcalqZ/g1LsDr7B7B0h/PxfMD9tvUCbwBK82KQsTSrk4k5UrQ48KlgedV+wDKiHbIMz6+DMWstjFVCiuiVhBHRA+NejinsVqrhXYUSdEZyPPs/6S+tZf3E94/aOY9zecdT0rEmvp3vRedgGbs1fwI0//uDSwJewCQzE46WBOLdvj8hhrRlFUe7S6XRMmzaNp556CpPJxEsvvZSrJJFb6oqimLgSncjc3Rf5a89lbsSnUNHbkRcal6FnkD8ONo/wfUBKuH4mLWmsg0s7LbOyvKtCw/9Bjb5gyHgVExITwoZLG1gbvJbTN09T2a0y7zZ4l7puNbm1YgVRv80i5fx5dD4+uD//PK59+6C1ch+sUnwVhiuKwig3VxQqURQzSakmVh25yuydIRwJjcFWr6G0uz0+Lnb4Otvi42KLr4vlp4+LLb7Odjjb6XJ+9ZF4E06uhD0/QvhRsHWFoOeh/mBwK5OhqZSSdRfXMWnfJK7GX6Vt6baMrDcSfwc/4rZu5cavv5GwZw8aR0dce/fGqU1r7GrVQhgMef+HUR5bKlFkTiWKxzhR3Cal5NDlaFYcvkpYdALhMUmEp93jce9/cju9Fl8XW5oHePFWuwCcbHPQHSQlXNoFu3+AkysACQEdLFcZ5ZpbBtbTJBmTmH18Nr8c+wWj2cjzVZ9nSM0hOOgdSDx6jBu//catdevAZELY2WFfty4OjRpi36gxtoFVEFpt3v5xlMeKShSZU4lCJYospZrMXItN5mpMEuExSVyNSSQ8JolLNxLYcDICH2dbvuxRg1ZVcjGlNSYM9v0C+2dZpud6BULDoVDzaTDcvUkoIj6Cbw98y4oLK/C08+S1Oq/RrWI3NEKD6dYtEvbsIX7XbuJ37STl3HkANC4uODSoj32jRjg0aoShfHk1EK7kikoUmVOJQiWKh3Lw0k3eXXSEMxFx9Kjjxyedq+LmkItuoNQkOLbIcpURfgRsXaDJCGj6pmXAPM2Ra0cYt3ccR64dIdA9kPcavEdQiYwV5FMjI0nYvZv4XbtI2LmL1CtXANB5eeE5fDhuT/fNk3NWij+VKDKnEoVKFA8t2Whi+ubzzNh8Dhc7PZ91q0anGr65+xYvJVzeDTu+g1MroXRj6DkTXEunayJZHbyayfsnE5EQQfuy7Xmr3lv4OPhkcjhJ6uXLxO/aRcyixSSdPEn51asw+Pvf11ZR7qUSReZykyisvR6FUsjY6LSMfDKAFSOa4edmx/A/D/K/3/cTcSsXy0gKAaUbQb+50PMnCD8G3zeD40vSNRF0Kt+J5d2X80qtV9h8eTPdlnZjzvE5GM3Gew4nMJQujVvfvvhN/Ra0WiK/mXDvuyqKkk9UolAyFejrzOJXmvBBxyr8e+YabSf9y/y9l8n1FWjNvvDyVvCsCAtehGXDIeVuzSh7vT3Dag9jSbcl1C1Rl2/2fUO/lf04fO1wpofTlyiBx5DBxK5bR8LevY9whopS+KWmpvLCCy9Qo0YNAgMD+eqrr+68tnbtWipXrkzFihX5+uuv8zUOlSiULOm0GoY2r8DaN5oT6OvMO4uO8Nwve7h8IyF3B3IvDy+tgyfegoN/wI/N4cqhDE1KOZViepvpTG45mejkaJ5b/Ryf7fyMmOSY+w7nMXAgOl9fwr/6Cml62JVBFKXwW7BgAcnJyRw9epT9+/fz448/EhIScqfM+Jo1azhx4gR//fUXJ06cyLc4VKJQHqicpwPzhjRiTPfqHLx0k6embOVY2P0f4NnS6qHNJ/DCckhJgJ/bwo5pYDbfaSKEoG2Ztizrvoznqj7HkrNL6Lq0K8vOLctwJaOxs8P7rbdIPnGSmKVL8+gsFSX/hISEEBgYyJAhQ6hWrRrt2rUjMTHxgfsJIYiPj8doNJKYmIjBYMDZ2TlDmXGDwXCnzHh+USU8lBzRaATPNSpD6yre9JyxnTf/PsSKEc0sxQhzo1xzeGW7pQvqnw/h/Cbo/j04lbjTxEHvwNv136Zrha6M2TWGj7Z/xNJzS/mo0UdUcK0AgHOnjtz84w8iJ0/B6an2aB0dsnpHRblrzXuWG0Xzkk8N6PDgrp+HKTPeu3dvli1bhq+vLwkJCUyePBl3d/cCLzOuriiUXPFztWN871qcjYxj4j+nH+4g9u6Wge5Ok+DidvihKZxdf1+zyu6VmdNhDp82/pQzN8/Qe3lvpuyfQqIxESEEJT54H9P160TNnPmIZ6Uo+S+zMuNvv/02hw4duu8xdepUwLJAkVar5cqVKwQHBzNx4kQuXLhQ4LGrKwol11oEeDGgUWl+3hZMm8ASNCrvkfuDCAH1B0GZJrBwEMztDdV6wpOfZZhGqxEaegf0pnXp1kzcN5Ffjv3CmZtnmNF2BnY1a+LSrSs3Zs3CtW8fNV1WebAcfPPPLw9TZvzPP/+kffv26PV6vL29adq0Kfv27aNUqVIFWmZcXVEoD+WDjoGUdrfnrfmHiU1KffgDeQfCkE3Q4l04vRqm1YdNX0ByXIZm7rbufNnsS94IeoP/wv5jX7jlPhmvN99U02WVIutBVxSlS5dm06ZNAMTHx7Nr1y6qVKlC/fr1OXv2LMHBwaSkpDBv3jy6du2ab3GqRKE8FHuDjkl9a3M1JpExKx9xtoXeFlp9AMP3QWAX2PoNfFcXDv2ZYbAb4JnAZ/Cy8+K7g98hpUTv44PH4EFquqxSLL366qvExcVRrVo16tevz8CBA6lZs2aGMuOBgYH07dtXlRnPLXVndsH5Zt0ppm8+z0/P1+PJqiUevENOXN4Da9+DsP1Qsg60/9pyA1+av079xdjdY/mx7Y808WuCOTGR8x07oXVzpdyCBaqIoJKBujM7c+rObKXAvN4mgEBfZ95ffISouOS8OWipBjBoA/SYCbER8OtTsGAgRF8CoFelXvg6+N65qsg4XTb/pggqyuNKJQrlkRh0GiY/XYtbiUY+WHI093duZ0WjgVpPw4h90OI9OL3GMn6xcQwGYwqv1HqFY1HH2Hx5M2CZLmtXqxaRUyZjiot/wMEVRckNlSiUR1bFx5m32gWw7ngEiw+E5e3BDQ7Q6n1LwgjsAv9NgPnP0aVCF0o7lWbaoWmYpfnudNlrarqsouQ1lSiUPDH4ifI0KOvO6OXHCYt+8B2nuebiD71+hifHwPlN6C7uZFjtYZy9eZZ/Qv4BwK5WLZy7duHGrFmkhIbmfQyK8phSiULJE1qNYEKfWpil5O0FhzGb82mSRIOh4OQLm7+kfZmnqOhakemHpt+pOOs9cqRluuyEifnz/oryGFKJQskzpT3s+bhzVXacj2L2zpBc7ZtqMj+4EVim0jYfBZd2og3ewvDawwm5FcLKCystL9+eLrt2rZouqyh5RCUKJU89Xb8Urat48/WaU5yLjM20TarJzNHQGGZtD+a1vw7SbNwmqn2yju3nrufsTeo8Dy6lYdMXtC7VikD3QH44/AOpJsuNfx4vvYTOx0dVl1WKvKioKFq1aoWjoyPDhw+/sz0hIYFOnTpRpUoVqlWrxnvvvXfnteTkZJ5++mkqVqxIw4YNCQkJeeQ4VKJQ8pQQgq971cDeoGXk/MOkmsxcj0vmn+PhfL3mFH1/3EmN0evoMm0bo1ecYHdwFDX8XPB3s+O1vw4SHpODBZJ0BmjxDlw5iDizlhF1RhAWF8aSc5aFkdJPl42ePz+fz1hR8o+trS1jxoxhwoT7Kw+MGjWKU6dOcfDgQbZv386aNWsA+OWXX3Bzc+PcuXO8+eabvPvuu48ch0oUSp7zdrJlbI8aHAmNofFXG6n3xQaG/r6fn/+7QLLRTP8Gpfmufx22v9eaXe+34fsBdZn5fF0SU00M//NAzrqhavW3rHOx+Uua+Tahjncdfjz8I0lGS6Jx7twJ+8aNiJw4idTIyHw+Y0XJ3sOWGXdwcKBZs2bY2tpm2G5vb0+rVq0AMBgMBAUFEZo2gWPZsmW88MILAPTu3ZuNGzc+8rR1VRRQyRcdavgyvFVFzkbGElTajaAybtTwc8myLHlFbye+7lWT1/46yLg1p/ioc9Xs30Crg5bvw+IhiJPLGFFnBC+te4n5p+fzfLXnEULg++mnXOjajYixX+E/ZXI+nKVS1IzbM45TN07l6TGruFfh3QYP/tb+MGXGcyI6OpoVK1bw+uuvA2QoQa7T6XBxcSEqKgpPT89cnFVGKlEo+WbUU5Vz1b5rrZLsD7nBz9uCqVvGjQ41fLPfoXov+G8ibPmK+sN20dC3Ib8c+4XeAb2x19tjKFsWz1de5tq3U4nd0g2nli0f/mQU5RFlVmb8o48+4u23337oYxqNRvr3789rr71G+fLl8yjS+6lEoRQqH3QK5FBoDG8vPEIVX2fKeWazIJFGa7mqWPACHF3AiDojGLB6AHNPzmVIzSEAeAwaRMzKVUR8PgaHBg3Q2NsX0JkohVFOvvnnl4cpM/4gQ4cOpVKlSrzxxht3tvn5+XH58mX8/f0xGo3ExMTg4fEQSwGko8YolELFRqdlxrNB6LSCV/7YT2LKA2YtBXa1rDC25StquVelhX8Lfjv+G7dSbgEgDAZ8P/+M1CtXuDZtegGcgaLk3IPKjGfno48+IiYmhilTpmTY3rVrV2bPng3AwoULad26NUKIR4rT6olCCNFeCHFaCHFOCPFeJq+XFkJsFkIcFEIcEUJ0tEacSsHxc7VjytO1OR0Ry8fLjmU/EKfRQKuP4GYIHPqTV2u/SmxKLHOOz7nTxL5uXVz79OHG7NkknTyZ/yegKHmobNmyjBw5klmzZuHv78+JEycIDQ3lyy+/5MSJEwQFBVG7dm1+/vlnAAYNGkRUVBQVK1Zk0qRJfP31oy/WZNUy40IILXAGeBIIBfYC/aWUJ9K1mQkclFJ+L4SoCqyWUpbN7riqzHjxMOmf00zddI5xvWrwdP3SWTeUEn5uC7Hh8NoBRm57n+1h21nTaw3utu4AmGJiON+xE/qSJSk77y9VivwxosqMZ64olRlvAJyTUl6QUqYA84Bu97SRgHPa7y7AlQKMT7Gi19sG0KyiJx8vO86xsJisGwoBrT+EW6GwfzbDaw8nyZTEb8d+u9NE6+JCifffJ+noUW7++VcBRK8oxYe1E4UfcDnd89C0bemNBgYIIUKB1cCIzA4khBgqhNgnhNh37dq1/IhVKWBajeDbfrVxtzcwbO4BYhKzWXK1fCso0xT+m0B5ex86levEX6f+IjLh7j0Uzp064tC0KdemTCE1PLwAzkBRigdrJ4qc6A/MklL6Ax2B34UQ98UtpZwppawnpazn5eVV4EEq+cPD0Ybpz9bhSnQioxYcznq8Qgho9SHERcC+X3il1isYzUZ+PvpzuiYCn9GfIo1GIr78soDOQFGKPmsnijCgVLrn/mnb0hsEzAeQUu4EbIGHv3NEKXLqlnHn/Y6BrD8RwcytF7JuWLap5cpi22RK2bjSvWJ3FpxZwJW4u72VhlKl8Hz1VWLXbyB248YCiF5Rij5rJ4q9QCUhRDkhhAHoByy/p80loA2AECIQS6JQfUuPmZealqVjDR/GrzvN7gtRWTds/REkRMHuH3m51ssIBDOPZFzIyGPgi9hUqkT4mC/UaniKkgNWTRRSSiMwHFgHnATmSymPCyE+F0J0TWv2FjBECHEY+At4UVpzqpZiFUIIxvWqiZ+rHWNXZzPF1b8eBLSHHVPx0djSt3Jflp5byqVbl+4eS6/H5/PPMEZEcP27nJVJUJTHmbWvKJBSrpZSBkgpK0gpv0zb9omUcnna7yeklE2llLWklLWllP9YN2LFWpxs9TzfuAyHQ2OyLGEOQKsPICkGdk5ncI3B6DV6vj/8fYYm9nXq4NrvaW78/geJx47nc+SK8nCyKjMOkJKSwtChQwkICKBKlSosWrQIUGXGFYVutf3QagSLslub27eW5Y7tXTPwTE2lf2B/Vl1Yxfno8xmaeY8cidbDnfBPPkEajfkcuaLkXnZlxr/88ku8vb05c+YMJ06coEWLFoAqM64oeDnZ0CLAiyUHwjBlt9xqm09AmuHvAbxU+Rns9fZMP5SxhIfWyQmfDz8k6cQJbv75Zz5HrjzO8rrMOMCvv/7K+++/D4BGo7lTHVaVGVcUoFeQP5tOHWDn+SiaVcpiApxnJejxA8x/HtcNn/Nc4AB+OPIjJ6JOUNXjbglzp6eewr5hQ6J+m4Xbs8+qO7aLufCxY0k+mbdlxm0Cq+DzwQcPbJeXZcajo6MB+Pjjj9myZQsVKlRg2rRplChRIl/KjKsrCqXIaRPojbOtjkUHQrNvWLUbNH8HDs3lufhUnA3O911VCCFw6/c0xqtXid+5Kx+jVh53mZUZf9iigEajkdDQUJo0acKBAwdo3Lgxo0aNyrfY1RWFUuTY6rV0rlWSJQfCGNPdiKNNNv+MW74PkSdw3vg5A1u8zLcXV3Io8hC1vWvfaeLYpg1aFxdiFi/CsVnT/D8BxWpy8s0/v+RlmXEPDw/s7e3p2bMnAH369OGXX34BVJlxRbmjV5Afiakm1hy9mn1DjcbSBeVZmWd2/4m7wYVph6ZlbGIw4Ny1K7HrN2BKu6RXlILwsFcUQgi6dOnCli1bANi4cSNVq1q6VItlmXFFeRhBpd0o62HP4uxmP91m4wT9/8IeweDYBHZf3c2eq3syNHHt1ROZmkrMipX5FLGiPJzMyowDjBs3jtGjR1OzZk1+//13Jk6cCBTDMuP5RZUZfzxM3XiWSevPsO3dVvi75WDlugtbSP6jJx3LlMHPqzqzO8zJ8E0ruGcvpJSUX7I4H6NWCpoqM565olRmXFEeWo86lkLDSw/m4KoCoHxLbNqN5X/XIzl47RDbr2zP8LJLr54knzxJ0okTWRxAUR5PKlEoRVYpd3salnNn0YGwnM8Tb/g/elToil+qke92jMmwn0vnzgiDgehF6opCUdJTiUIp0noF+RN8PZ6Dl6NztoMQ6DtP4WXhzomEK2w6knFxI6cnnyRm5UrMycn5E7BiFcWxi/1R5PbvoRKFUqR1qOGDrV7Dov0PuKciPZ0NnXsvoKxRMm3fJMxxdxc3cu3VE3NMDLEbNuRDtIo12NraEhUVpZJFGiklUVFRmd7tnRV1H4VSpDnZ6nmqmg8rDl/hky5VsdHl7M5qnUtJhtUexjvHvmfdgr50eO4f0Bmwb9QIXUlfYhYtxqVTp3yOXikI/v7+hIaGola+vMvW1hZ/f/8ct1eJQinyegX5s+zQFTadjKRDDd8c7/dU0MvMPLeIGbfCeHLDp+jaf4XQaHDt0ZPrM2aQGhaG3u/elXmVokav11OuXDlrh1Gkqa4npchrWtGTEs42Dy7pcQ+N0DC88QeEGPSsOP4HRFvWrHDp0QOA6KVL8zpURSmSVKJQijytRtC9jh9bTl/jelzuBqFbl2pNNddK/OjqSOq/4wAw+Pvh0LgRMYuXIM3m/AhZUYoUlSiUYqFXkD9Gs2T5oSsPbpyOEIJhdd8gTKdjxbllcMOyJrdLz16khoWRsHt3foSrKEWKShRKsRBQwokafi4sPpi77ieAJ/yeoLpbZWa6OpO6xVLuwKltGzTOzuqeCkVBJQqlGOkZ5MexsFucDs9mmdRMCCF4Jeg1wnRalgevgmtn0Nja4tK5E7H//IMpJiafIlaUokElCqXY6FqrJDqNYHEuB7XBclVRw60KP7m6kLplLGDpfpIpKdxavTqvQ1WUIkUlCqXY8HC0oWVlb5YcDMNoyt0gdIariovrIeIEttWqYlOlCtELF+VTxIpSNKhEoRQrvYL8iIxNZvv5qFzv28yvGTXdqzLTzYXUzV8ghMC1Z0+Sjh8n6VTeLp+pKEWJShRKsdI60BsXO33uSnqksVxVjOCKTsuysC1w9TDOXToj9HqiF6tBbeXxpRKFUqzY6LR0qeXLuuPhxCal5nr/piWbWq4qXN1I3fQlOjc3HNu04dbyFZhTUvIhYkUp/FSiUIqdXkH+JBvNrH7QMqmZuH1VcVWnYWn4Ngjdh2uvnpiio4nbtDkfolWUwk8lCqXYqV3KlfKeDizKyTKpmWhasik1Parxk5sbqZu+wKFJE3Q+PkQvUoPayuNJJQql2BFC0KuuP3uCbzB6+XHORuT+vophdUZwVathybW9iNA9uPToTvy2baRezf1ViqIUdSpRKMXSc43L0LVWSebuvsiTk7fS98edLDsURrLRlKP9m5RsQi3PGparis1f4NqjB0hJzLJl+Ry5ohQ+KlEoxZKzrZ6p/euw8/02vNehCuExSbw+7xCNv9rE2NUnCb4en+3+QgiG1R5OuFawJOoQBmMw9g0aEL1osSoUqDx2RHFc9alevXpy37591g5DKUTMZsn289eZu+sS609GYDJLmlX05JmGpXmyagn02vu/M0kpeW71s0REHGGV9CXRawRX3nkX77dH4TFokBXOQlHylxBiv5Sy3r3b1RWF8ljQaARPVPLih+fqsuO91rz1ZADB1+MZNvcATb7exPdbzt/XLWUZq0i7qog+gXNlW5zatyfymwncWvePlc5EUQqe1ROFEKK9EOK0EOKcEOK9LNr0FUKcEEIcF0L8WdAxKsVLCWdbRrSpxNZ3WvHri/Wo6uvMuLWnaD/lP/49k3G5zMa+jantWYuf3N1J/fdLSn41FrvatbnyzjskHj5spTNQlIJl1UQhhNAC04EOQFWgvxCi6j1tKgHvA02llNWANwo6TqV40moErauUYPZLDZg1sD5SSl74dQ+v/LGfsOhE4PZVxatEaGBJ7Dk0FzfjP2M6Om9vLr8yjJTQ3N8BrihFjbWvKBoA56SUF6SUKcA8oNs9bYYA06WUNwGklJEFHKPyGGhZ2Zt1bzbn7acqs/l0JG0n/sv0zedINppo5NuIOl61+cndnZTNX6JzdaXUjz8iTSYuD/2fKkOuFHvWThR+wOV0z0PTtqUXAAQIIbYLIXYJIdpndiAhxFAhxD4hxL5r165l1kRRsmWj0/Jqq4psGNmC5gGefLPuNO2n/Md/Z6/zSu1hRGhgcUII/PMhNuXK4v/dVFIuXyb0tdeRqryHUoxZO1HkhA6oBLQE+gM/CSFc720kpZwppawnpazn5eVVsBEqxYq/mz0/PlfvTnfU87/uYfZGPdXca/GTd0mSd8+Af8fh0KABJb/8goTdu7n6yacUxxmEigKWD2FrCgNKpXvun7YtvVBgt5QyFQgWQpzBkjj2FkyIyuPK0h3lwU9bLzBt8zk0dg3Q+h1meWBr+mz5Cmyccek6jJRLl7k+bRr60qXwGjbM2mErSp6z9hXFXqCSEKKcEMIA9AOW39NmKZarCYQQnli6oi4UYIzKY8xGp2V460psGNmCpn6NMSX6Mc2YgKlKF1j3Phz4Hc9Xh+HSrSvXp35HzIoV1g5ZUfKcVROFlNIIDAfWASeB+VLK40KIz4UQXdOarQOihBAngM3A21LK3K9KoyiPwN/NnpnP16OMrhM3Uq+wrGZfqNAGVryGOLEUnzFjsK9fn6sffEjCXnWxqxQv6s5sRcmFMxHR9FzeDRcbN7Y9PRcxtxeE7oP+f2Hyqk9I/2cw3bhBmXl/YVOunLXDVZRcUXdmK0oeCCjhShOvXtyS5/nt6C545m/wDoS/n0Mbc4JSP/4AGg2X//cyxhs3rB2uouQJlSgUJZfGPzUIYXZk2oGfSNI6wnNLwLUUzO2LQReF/4zpGMPDifj6a2uHqih5QiUKRcklVzsHupTtQ6rNcT5ftxEcPOG5pWDnBr/3xN7PDudOnYj/d6uqNKsUCypRKMpDeLvxS2gwsOTCn5yJiAUXP3h+KWj18Ht3HGpVxBQTQ9LJk9YOVVEemUoUivIQXG1d6VGhJzrng7yzZCtmswSPCpZuqNREHIKnABC/Y4d1A1WUPKAShaI8pCG1X0QIOJmwivn70irRlKgGAxajk1HYeBlI2LnTukEqSh7I00QhhPDIy+MpSmHm5+hHh3LtsXXfw9h1B7kel2x5wb8utP4IB7cbJOzdizkpybqBKsojyusrCnWnkfJYean6S5hFMil22/hyVbrxiLoDsS/niEw1knjggPUCVJQ88NCJQghhFkKY0n6ahRAmoGy63xWl2KvsXpmmJZviVGIXSw6FsO3sdcsLelvse70GQhK/ep51g1SUR/QoVxTfA78DXlJKjZRSC1xM97uiPBYGVh9IkjmaEiWP89HSoySlWr4naZsMwq6EIP6/zVAMKyAoj4+HThRSyleBn4GlQogXbm/Ok6gUpQhp4NOAqh5VcfTeRkhUHDM2n7O8oDPg0LQ5SRGpmPYusG6QivIIHmmMQkq5DWiDpQLsesCQJ1EpShEihGBg9YFEJoXSpEY43/97nnORcQA4dB8MCOLnjQN1851SRD1wPQohxPfAYeAocERKGZv+9bQlTD8SQtQAmuZLlIpSyLUt3RZ/R3+M+k3Y6l/ggyVH+XtoI+xq10ZjZ0P8mWs4HV9MdPmuXItL5npsMtfikrkWm4yUMLBpWXRaNVtdKZxysnBRZ+B/pHUrCSEuAUfSHgeBf6WUUVLKo1iSiaI8dnQaHS9We5Evdn/BgJZmvl97gyFz9mEySzqXqEzV8CMkLPyYdskGTNw/hFfK3Y721X2tELmiPNgDE4WUslTa0qM173m8ATgARiHEEuAtKWVo/oWqKIVbt4rdmHF4BheNq2gb+CJHw2LwcrLhSoUaVA05Qqn4CH5rEExM5T54Otrg5WSDh4OBp6ZsZdGBMJUolEIrR0uhSimjga1pjzuEEJWBTliSxm4hRH0p5ZU8jlFRigRbnS39q/Rn+qHpLO76BpXc6gOQfM6XCxvnkpBcgeZXfoGer4Du7nBet9ol+W17CDfiU3B3UMN8SuHzqIPZp6WUk4DagBEYnQcxKUqR1a9yP+x0dsw6PuvONkOFCui8vIhPCoDoi3Dojwz79Azyx2iWrDyivmMphVOejJ5JKW8AU4EOeXE8RSmqXG1d6VmpJ6svrCY8PhywzIpyaNKY+GPBSL/6sHUCpN4t6xHo60ygrzOLDoRZK2xFyVZeTrO4CHjm4fEUpUh6rupzSCST90/GaDYCYN+4MaabN0ku+yLcCoP9szLs07OOH4cvR3P+WlzBB6woD/DARCGEWC2E+FoI0V8IUU0IkdVd1w2By3kbnqIUPX6OfgyuMZjVwasZsWkEsSmxODRuDED8pRQo+wT8NxFS4u/s0612STQClqirCqUQyskVhROW6bFzsUyJjRNCHBBCzBJCvCOEGCyEmAG8BvyZj7EqSpExvM5wPm70Mbuu7OK51c8Rbp+CoUIF4nfuhNYfQXwk7PnpTntvZ1ueqOTFkoNhlrUtFKUQeWCikFI+IaV0A8oC3YAxwBmgPvAlMBMYiKWcx9h8i1RRipi+lfvyw5M/cC3xGs+uepb42hVI2LcPs08QVGwL27+FpFt32vcM8iMsOpHdwTesGLWi3C/HYxRSyktSypVSyrFSyn5SymqAHeAHOEopX027S1tRlDQNfRsyt+NcXGxc+FazGZmUROLBQ9DqQ0i8Abt/uNO2XVUfHAxalhxUtyMphcujTo81SimvSilVWXFFyUJZl7L80fEPbOoFYRKwafFkzCVrQ5XOsOM7SLBcQdgZtHSs4cvqo+Ekpqj/pZTCQxWXUZQC4GLjwrddfuJmBU/Mew/xxuY3SHhiJCTHws5pd9r1CPIjLtnIPyfCrRitomSkEoWiFBC9Rk/ldn2oFC7Ye24Lz+//iquBHWHXD3BkAZhNNCrnQUkXW5YcVLOflMJDJQpFKUCOTZogzJJvXYYSFhdGf3mZw56lYPFgmNEYzYnF9Kjjy9Yz14iMVWttK4WDShSKUoDsatZE2NtT+vRN/uj4B3Z6RwY5mLjUZRIIDSx8idfODOQpsZvlalBbKSRUolCUAiQMBuzr1yN+504quFZgVvtZmKSJucZIeGUH9P4VG42Z7w3f0vrfXnByhVpGVbE6lSgUpYA5NG5MSnAwqVevUsKhBB3KdmDpuaXEGuOhei8YtoutNcYiU5Ph7wHw4xNwarVKGIrVqEShKAXMoUkTAOJ37ATg2arPkmBMYOm5pZYGGi3V2w+ho/Eblpf/BJLjYF5/mNkSLu+xTtDKY83qiUII0V4IcVoIcU4I8V427XoJIaQQol5Bxqcoec2mUiW0np6Wch5ANY9q1PGuw58n/8Rkttw/4e5g4InKvnxxuRamV/dCt+mQEAW/94CrR6wZvvIYsmqiSCswOB1LefKqQH8hRNVM2jkBrwO7CzZCRcl7QggcGjcmfudOZFp30jOBzxAaF8p/Yf/dadcryI/I2GS2X4iGOgNg0D9g6wJ/9oUYNdCtFBxrX1E0AM5JKS+klf+Yh6We1L3GAOMANV9QKRYcGjfGFBVF8pmzALQp3YYS9iX44+TdRY1aB3rjbKtj8YG0pOBcEp5dYKk6O7cPJMVYI3TlMWTtROFHxtLkoWnb7hBCBAGlpJSrsjuQEGKoEGKfEGLftWvX8j5SRclDDo0bARC/cwdguRmvX5V+7L66m7M3LcnDRqelc62SrDseQVyyZV0LSlSDvnPg+hmY/zwYVXk1Jf9ZO1FkSwihASYBbz2orZRyppSynpSynpeXV/4HpyiPQO/ri6FcuTvjFAC9K/XGVmvL3JNz72zrFeRHYqqJtcfSlfSo0Aq6fgcXtsCK19VsKCXfWTtRhAGl0j33T9t2mxNQHdgihAgBGgHL1YC2Uhw4NG5Mwt59yBTLVYGrrSudyndi5YWVRCdFAxBU2o0yHvZ3u59uq/0MtPwADv8JW74u4MiVx421E8VeoJIQopwQwgD0A5bfflFKGSOl9JRSlpVSlgV2AV2llPusE66i5B2HJo2RCQkkHrk7i+mZwGdINiWz8OxCwDLw3aOOHzsvRHElOjHjAVq8A7WfhX+/hoN/oCj5xaqJQkppBIYD64CTwHwp5XEhxOdCiK7WjE1R8pt9gwag0dy5nwIgwC2Ahj4NmXdq3p31tnvW8UdKWHronkKBQkCXb6F8K0sX1LmNBRm+8hix9hUFUsrVUsoAKWUFKeWXads+kVIuz6RtS3U1oRQXWmdn7GrUyDBOAfBs4LNEJESw8ZLlg7+0hz31yrix+EDYnem0dw+itwxue1WB+S9A+NGCCl95jFg9USjK48y+SWMSjxwhJfTu1UJz/+b4O/pnGNTuGeTPucg4joXdItlo4vKNBPYE32DZoTB+2H2NSV5fcNNkw/WZ3ej8xd/U+fwfzkXG5lmc5yJjeW/REWISUvPsmErRIe77hlIM1KtXT+7bpy48lMIv+exZQp7uh9Dr8f1qLE6tWwPw+4nfGb93PPM6z6OaRzViElKpP3YDGgFJqeb7juNko6Ox41W+TXifaIMPvZM/oaRPCf4e2hiNRjxSjEaTmZ7f7+BIaAzPNizNlz1qPNLxlMJLCLFfSnnfZCGVKBTFylJCQggdOZLkEydxe/45vEeNIp5k2i5oS5vSbRj7xFgA5u6+yLGwGHxd7PBxscU37VHC2RYnW73lYOc3wdw+hLvXo1noMMb0rEP/BqUfKb6f/7vAF6tOUquUK0dCo1n8ShPqlHZ71NNWCiGVKBSlEDOnpBA5/htu/vEHttWq4Td5EhOu/sGCMwtY33s9nnaeOT/YwbmwbBirHXrwXnx/NrzVAm8n24eK6/KNBNpN3krTih5Mfro2bSf9i4eDDcuHN0WnVT3XxU1WiUL9l1aUQkBjMODz0Yf4T/uOlMuXCe7Rk76X/TCajcw/PT93B6vzLDR8mY7xS2hj/I8xK08+VExSSj5YchStRjCme3WcbPV82qUaJ67eYs7Oiw91TKVoUolCUQoRp7ZtKb9kMTaVKpH60dd8stWbJcf+JsWUy1IdT46BUo0Yb/iJ00d2s/l0ZK5jWXIwjP/OXued9pXxdbEDoEN1H1oEeDFp/RnCYx6t9NraY1fpPn07UXHJj3QcJf+pRKEohYzez48yv8/BY8gQqm+/wqgfI9m0dU7uDqIzQN/Z6Oxd+NX2W75avJuEFGOOd4+KS2bMyhMElXZlQMMyd7YLIfi8WzVSTWbGrDyRu5jSOXDpJq/NO8Shy9EsP3zloY+jFAyVKBSlEBJ6Pd5vjcT/p5m4J2jxfW0SNxctuv8+iuw4+SD6zKYkkbydMJkp60/neNcxK08Ql2xkXK+a982aKuPhwPBWFVl19Cr/nsl9Ac7LNxIYMnsfvi62VPR2ZNkhlSgKO5UoFKUQc3riCS599wZnSkL4hx9xoWMnLr00iCvvvU/kpMnc+GMut/75h8TDh0m9ehWZes99DmUao3nqC57U7sdm5xSOhT24NPmW05EsPXSFYS0rUqmEU6ZthrYoT3kvBz5ZdoykVFOOzycmMZWBs/ZiNEt+fbE+fer6c+hyNBej4nN8DKXg6awdgKIo2etQ/xmefP4XBp3ypV20H6mRkSTvDsZ47RoY7+lOEgKtuzuGcmUpNWMGWmdnaPgyKRf38ubJBYz+uwaBb4xAm8W9FfHJRj5ccoyK3o4Ma1Uhy5hsdFq+6FadZ37ezYzN5xjZrvIDzyPVZGbY3P1cjIpnzksNqeDliF0tLV+tOcWyQ1d4rU2lXP1dlIKjEoWiFHL2ent6Vu7NVNPvdOw1A38HHwCk2Yzp5k2MERGkRkZijIzEGHmNlJAQbq1cSdx//+HSqRMIgaHHd9wKO8KbMeNYtDGIvk82y/S9Jv5zhrDoRBa+3BgbnTbbuJpU9KR77ZL88O8Futfxo7yXY5ZtpZR8vPQY289FMaFPLRpX8ACgpKsdDcq5s/RQGCNaV0SIR7s5UMkfqutJUYqAflX6IZFM2jeJhNQEAIRGg87DA9uqVXFq2RK3vn3xGv4qJcd9jcbFhfjtO+4ewOCA0/PzsNGYqbZtOFeu37jvPQ5djmbWjmCea1SGemXdcxTXh52qYqPX8PGyY9mOn/y49QLz9l5mROuK9K7rn+G1brVLcuFaPMev3MrReyoFTyUKRSkC/Bz9GFR9EGtC1tBjWQ+2hm7Nsq3Qai1rcm/fnuHDW3hWJL7T91QTwVyY9QrSfLcUSKrJzHuLjuDtZMs77R/cjXSbl5MN77SvwvZzUVnOXlpz9CpfrzlFl1olGflkwH2vd6zui14r1OynQkwlCkUpIl4Leo1Z7Wdhq7Pl1Y2vMurfUVxLyHzWkUPTJhgjIkg5dy7Ddq963TlQdjDN4tZyfOXUO9tnbr3AqfDYOzfW5cYzDUpTy9+FMStPEpOYcTD90OVo3vj7EHXLuPFN75qZdi25ORhoEeDF8kNXMJmLX6WI4kAlCkUpQuqWqMuCLgt4tfarbL60mW5LuzH/9HzMMmOhQMcmTQCI2779vmPUfPZr9umCCDgwhrgLu7lwLY5vN56lUw1fnqxaItcxaTWCL3vU4EZ8MhPW3Z2CG3ozgcGz91HC2ZaZz9XFVp/1mEfX2n6E30piT/D9XWKK9alEoShFjEFr4OVaL7Oo6yICPQIZs2sML6x5gbM3z95po/fzs6zJnX6cIo1Or8f26V+JkG6Y/hrAVwv+w1an4dOuVR86pup+LjzfuCx/7L7I4cvR3EpK5aVZe0kxmvj1xfp4ONpku3/bQG/sDVqWHw7Ltp1iHSpRKEoRVdalLD+3+5kvm31JyK0Q+q7oy7cHviXJaCmt4dCsGQl792JOvr9ERvVK5VhddTy2KTcZGzGUvyusxTvl0T6k32oXgJejDR8uPcqrcw9w4Vo8PwyoS0XvrGdD3WZv0PFUNR9WHw0n2Zjz+zKUgqEShaIUYUIIulboyvLuy+lUvhM/H/2ZHst6sCNsBw5NmyCTkkjcvz/TfZ/t3pU3bD/non11qlyYDd8FweyucHwJGHNZWwpwstXzSZeqHAu7xX9nrzO2Zw2aVMx51duutUsSk5jKv6dzf7e3kr9UmXFFKUb2XN3D57s+5+Kti9RxDOS9T0/g+twASr77XqbtE1KMGLQadPERcPAPODAbYi6DvaelCm3QC+CR9Y1395JS8tmKE/i72TH4ifK5ij3VZKbh2I00ruDB9GeCcrWvkjfUehSK8phINiWz5OwS/jr1FwNmnMUlWcPJyUN4uvLT+KTdrJclswnOb4b9v8HpNSBNUK4F1BsIlTtZig3mo4+XHmP+vsvs//hJHG3U/cAFTa1HoSiPCRutDf2q9GNpt6VUfKoXpSJMLNr1C08teoo3N7/J3vC9Wd8cp9FCpbbQby68eRxafwQ3gmHBizC1DsRG5Gvs3euUJNlo5p/j4fn6PkruqEShKMWUEIIq7fsBMMfjLQZWG8i+iH28tO4lei7vyfzT8+/c5Z0pZ19o/ja8fgiengu3Qi1dU9mQUrL47GL2R2Q+LvIgQaXd8HezY6mqKFuoqEShKMWYbWAgWjc39PuP80bdN1jfez1jmo5Br9EzZtcY2i5oy9yTc7M/iEYLgZ2hQmvYPwtMWa9r8eORH/l0x6d8seuLh4pXCEHXWiXZfu4612LVgkaFhUoUilKMCY0GhyZNiN+xE2k2Y6uzpXvF7vzd+W9+7/A7tbxr8fWer/n+8PcPXuui3iC4FQZn1mb68u8nfmf6oemUcirFuehznLl55qFi7lbbD5NZsvro1YfaX8l7KlEoSjHn0KwZpuvXST5z94NbCEFt79pMaz2NbhW6MePQDL47+F32ySKgPTj7wb5f7ntp8dnFjN87nral2zKnwxx0QseqC6seKt7KPk5U8XFi2SF1811hoRKFohRzDmnlPOK3bbvvNa1Gy+dNP6d3QG9+OvoTE/dNzDpZaHVQ90U4vwmizt/ZvDZ4LaN3jKapX1PGNR+Hp50njUs2Zk3wmvtKi+RUt9p+HLgUzaWotDEUUyoUwxmaRYVKFIpSzOlLeGNTqVKmdZ8ANELDJ40+oX+V/sw+MZuv93yddbIIeh40Otj3KwD/Xv6X9/97nzredZjccjIGrWX6bMfyHbkaf5WDkQcfKuYutXwBLCU9TKnwYwuY0w1Sshl8V/KNShSK8hhwaNqUxH37MScmZvq6EIL3G7zPC1Vf4M9Tf/L5rs8zvxpw8oEqneHgH+wJ/Y+RW0ZS2b0y09tMx05nd6dZ61KtsdPZsfrC6oeK19/Nnvpl3Vh66Ary0J8QeRyC/4W/nwVj1oPcl28k8PvOEMJjkh7qfZXMqUShKI8Bh2bNkKmpJGRzI6oQgrfqvcWQGkNYeGYhH2//GJM5k7pL9QdzWCYwfPPrlHYuzQ9tf8DRkLGek73enpalWrLu4jpSTan3HyMHutX241LkTVI3jwO/etBtuqXba8GLlquMdKSULNwfSodv/+PjZcdpPn4z7y8+Qsh1tRZ3XlCJQlEeA/b16iIMhkzHKdITQjCizgiG1R7G8vPL+WDbBxjNGafDnnb24hVfHzxNJmY+ORNXW9dMj9W5fGdikmPYceX+CrY50bGGL/10/2KIC4NWH0CdAdBxApxeDYsG35mmezM+hWFzDzBqwWGqlnRmwcuN6Vvfn0UHwmg9cQuvzzvIqXC1et6jUPfIK8pjQGNri329esRt386DVpwQQvBKrVfQa/R8e+BbUs2pjGs+Dr1GT3BMMEM3/A97vSM/BZ/GKzoU7L0yPU7jko1xtXFlVfAqWpRqkeuY3Q1m3rBZziFzIDXLtbJ8q20wxNL19M+HoLNla7XPGbXwKDcTUnivQxWGPFEerUZQv6w7r7WuxC/bgvlj10WWHbpC28ASDGtVgaDSbrmO5XFn9SsKIUR7IcRpIcQ5IcR9lcuEECOFECeEEEeEEBuFEGWsEaeiFHUOTZuScu48qeE5K48xuMZg3q73NusvrmfklpGExIQw5J8hAPz05Ez8hA3svX+q7G16jZ52Zdqx5fKW7O8Az8r+33A3Xefr5F7svXjz7vYmw0lt8QEcmUfo7y/jYqtjybCmvNyiAlrN3RX0vJ1teb9jINvfa82bbQPYd/EGPWfs4JmfdrH93PUH3zei3GHVRCGE0ALTgQ5AVaC/EOLe1VMOAvWklDWBhcD4go1SUYoHh2bNADJdzCgrz1d7ng8afsCWy1vosawHCcYEZj45k3LeNaBGHzi6EBKjs9y/U/lOJBoT2XR5U+6CTUmA/yZhKtOMw9oaLEu3nvbxKzF0PNCQ6cauPKPbxJoqq6le0jnLQ7naG3i9bSW2v9uajzoFci4yjmd/3k33GTs4eOlmlvspd1n7iqIBcE5KeUFKmQLMA7qlbyCl3CylvP11ZBfgX8AxKkqxYBNQCa2XJ/Hbsx+nuFf/Kv0Z3Xg0vo6+fN/2eyq7V7a8UH8QGBPh8Lws963tXRtfB9/cz37a+zPER6Jt8zHtqpVg9dGrJKWa+PHf83Sfvp2YJCPVn5sIjYah2/sjbPz8gfdZONjoGPxEef57txVf9azBtVtJvDRrL1eiM58Jlp8OXLpJ7+938Pq8g0zdeJZVR65yKvwWSamFc9Ema49R+AGX0z0PBRpm034QsCazF4QQQ4GhAKVLl86r+BSl2BBC4NikKXFbtiBNJoQ26zWs79UroBe9Anpl3OhbC/zrW+7Ubvg/EOK+/TRCQ8dyHZl1fBY3km7gbuv+4DdLjoXtU6BCGyjdiG4JESw7dIUO3/5H8PV42lfz4aueNXBzMEDAWEhNhG2TQG8PLd5+4OFtdFr6NyhNw3LudPluGyP+Osi8oY3Qawvme7PRZOb9RUcJv5VE+K0klqUrgCgElHKzp4KXAxW8HKng7Uh5N1uCynsWWHyZsfYVRY4JIQYA9YBvMntdSjlTSllPSlnPyyvzwTVFedw5NGuGKSaGpBMn8+aA9QbB9TMQ8l+WTTqW74hJmvgn5J+cHXP3j5AQBa0+BOCJSl54OBiIvJXE+N41+X5AkCVJgOWTtdMkqNUfNn8BO77LcejlvRz5uldN9l+8yTfrTud4v0f1+66LnI05Ra2qC1jYzY2DXT1ZXTOZOS4XmBa3gzf3/0W3P8dRd8xwyg3uiXOXFvz5wluYzdYbU7H2FUUYUCrdc/+0bRkIIdoCHwItpJSqpKSiPCSHJo0BiN++Hbsa1R/9gNV6wLr3LV1F5Zpn2iTALYCKrhVZdWEV/ar0y/54STGWD/uA9uBfFwC9VsPCV5pgp9fi42J7/z4aDXSdBsYk+Ocj0NlaZkflQJdaJdkbcoOZWy9Qr4wb7ao9YGGnR3Q9LplJ68/QLnUpIz4O5iaWbkABeAFeej06T090Xl7oKgWS7OJG5LFT1D6wgd/XHOSFTtZZ+c/aiWIvUEkIUQ5LgugHPJO+gRCiDvAj0F5KGVnwISpK8aHz8MCmaiDx27bh+fL/Hv2AelvL/Q27vodbVy1rWGSiU/lOfHvgW0JjQ/F3ymaYcecMSIq23DeRTjlPh+zj0Oqg50+WqbOrR4F3VSjbNEen8GGnQA5eiuatBYdZ5eNMaQ/7HO33MMavPUWyCKPuyRASHXT81srMEzW70r3RQHReXmhdXBCajB09/hcucKFjZy78NItjtStQ3c8l3+LLilW7nqSURmA4sA44CcyXUh4XQnwuhOia1uwbwBFYIIQ4JIRYbqVwFaVYcGzalIRDhzDF5dFdy3UHgtkIB+Zk2aRjuY4ArAnOdIjRIuEG7JoBgV0s4x+5pdVDr5/BpTSsfBOMKTnazUanZcazQQjg1T8PkGzMnwHlg5duMn9fKNUq7qfOBYln26dw7d6dsXIVBx2j0Lm53ZckAGzLl8e2ZUs6X9jOqDm7SEjJej2Q/GL1MQop5WopZYCUsoKU8su0bZ9IKZen/d5WSllCSlk77dE1+yMqipIdh6bNwGgkYc+evDmgRwXLwHM2ixqVdCxJkHcQqy6syvr+hZ3TLAPZLT/I/PWcMDhAp4lw/TTs+DbHu5Vyt2di39ocDYvhi5V5NH6Tjtks+XT5cbxckrAN3oZDEri3aceHDT+kvEt53vvvPSITsu4w8Rk6BMeUBAIObuHzFSfyPL4HsXqiUBSlYNkF1UHY2RGfRTXZh1J/EMRegTNZXzF0LNeR8zHnM1/QKP467PoBqveEEvfeSpVLAe2gajf495sM5dAf5MmqJRjavDy/77rI8sN5uxTr/H2XORIaQ4NaJ6lzxgh6PQ5NmmKvt2diy4kkGhN5Z+s795VLuc0+qA52derwQugO5u8OKfBFnVSiUJTHjMZgwL5B/QfWfcqVSk+Bs3+2d2q3K9vOsqBRcCYLGm2fYrkno8V9xRkeTvtxoDXAqrdytY7F209Vpm4ZN95fdITz1+LyJJSYhFTGrztN3bJ2HIpZwxMX7XBo0ACto2XcpYJrBT5u9DH7I/Yz/dD0LI/jMXgQ9jciGZB8nvcWHSHs3vs/Em5YJhVks1Ttw1KJQlEeQ45Nm5Jy8SIpoaF5c8Dbixpd2Jzlt3g3Wzea+DW5f0Gj2AjY8zPU6AteAXkTj7MvtPnEEs/RhTneTa/VMO2ZOhh0Gl6de4DElEcfr5i0/jTRCSk0qR2MfXgMrhHxOLZsmaFNlwpd6FWpFz8f/ZmtoVszPY5jq1YYypXj2ZCtmM2SN+cdwpR+yuzGz2H1O3Aj51dROaUShaI8hh6mnMcD3bOoUWY6lutIeHw4ByIO3N24bTKYUqDFO3kXC1i6w0oGWabvJua8VIevix1T+tXhdEQsnyw79kghnLhyi993XeTZhqVYH7aALuGWWWGOrVre1/a9Bu9R2a0yH2z7gKtx93ctCY0G95cGYj59iomVUtkTcoPpm89ZXgw7YBkjavg/8Kr8SDFnRiUKRXkMGcqVQ+frm7fjFE4lLDOWDv5huVs6E61KtbIsaBScVtIjJsySWGo/YxkUz0saLXSZYrl5b8PoO5tNMTEkHj9O7IYNGKOiMt21RYAXI1pVZMH+UObvu5xpmweRUjJ6+XFc7PQEVQ0lLC6MJy7aY1OpIgb/+6cI2+psmdhyIkazkVFbR2W6jodL165oPT0J/Hc53WuX5NuNZ9kfch1Wvw0OXtAyj7ru7qEShaI8hoQQODRtQvzOnUhjHvZp1x9suQ9iw2dwahVc2m3pikqKASmx19vTqlQr/rn4j+WD8L+JIM3Q/MGlN3LKnJJCcnAwcf/9x80tJ4gIb0bot8u40Lk9pxs05EzDRoT06k3o8BGE9H8GU2xspsd5vW0Ajct78MmyYw+1nsXyw1fYE3KDt5+qzMKzc6mk98P22AVLt1PSLcvf5h5lnMswuslojlw7wpQDU+57XWNjg/uAAcRv28bHVQ2UdLXln7mTIGwftBsDtvlzj4W1b7hTFMVKHJs2JWbhIhKPHsW+Tp28OWiZpuBbG3Z/b3mkpzWAvSednJxZbZPA9oX9aXl6CwQ9B26PvnqANJu5/sMPXJ/xPaRLfsJgQG9nh14bhn3HXuhLl0Hv7wdGI2Fvv8PVDz7Ab+pUxD21qrQawbf9a9Np6jZe/n0/Y7pXp1lFz/vaZSYu2cjY1Sep4edCQJkovjp1hPGp3cB4EUe3cJhSw5JQu0233LCYTvuy7TkQcYA5J+YQVCKINqXbZHjdrX8/rs+cSdIfs5k2ZBj+f8zhvH0Nytfoy4MjezgqUSjKY8qhcWMQgus//IBTq9bofEqg9/FBV6IEWlfXHH0g3kcIGLwR4sItU17jr0PCdYi/dud54/hruKaeYvWt07R08YMnRj3yuZgTErjy/gfErluHU/v2OLZsgaFUKfT+/ui8vBBn1sC8Z6CtCzQbyM4rO/nu4Hf079+Yin9s4Mbs2Xi8+OJ9x/V2smXGs0G88scBnvtlDwElHHmpaTm61/HDVp91UcXvNp0l4lYy3w+oy5yTn+FscKbalsMk2UjsQmZC5acs94ysegt8aoJvzQz7j6o3iiPXjvDxto8J6BJAKae7lY60Li649enNjbl/UrVmPFoRz4DoZxh88Aq96uZPcW1RHBfvqFevntyXzdrAiqJYhI0cya2168BszrBd2NhYEod3CXQ+Puh9SmAoWw6Xrl0Qev0jv+8Xu75g2bllbHl6Cw76B5TneIDUK1e4/Opwkk+fxnvUKNwHvph5kpv3LJdCNvNNnY5sCd+NXqNHms3M2xGEefteysyZg31QJldW0ZdIPf8vW6JLMPmInhMRCbg7GHi2YWmea1QGb+eM9afOX4uj/ZStdKvtx4g29nRZ0ZPBsUbazzThEOCG3+RpUKo+xF2DH5uDzgBDt4BdxpX3QmND6buiL6WcSzGr/SySjcnEpMQQnRxN7KVg3J/7gMhqyRzrU4P5saW5kRhN3fI2fNd2Is6GrNfnyI4QYr+Ust5921WiUJTHmzQaMUZFYQwPJzU8AmNE2s/wcFIj0n5GRkJqKnZ16uA3eRJ6n0crnncw8iDPr3mesc3G0qVCl4c+TsL+/YSOeA2ZkoLfpIk4Ns+8MGFsSiwz907kj7MLMQgNQ4Nep1P5TnRf1p2mTrV5ZeoFZEoK5ZYsRueeVgo96rylfPnheZYSJYA0OBHtWZstCRX4O9Kf46IiT9Yqx0tNy1HdzwUpJc//uodjl66z9ckwvj3+LYttBKsuexI9NwK/yZNw7tDhbmCX98BvHaDik9DvT0uBw3Q2XdrE65tfz/ScRiw3Ue+sZNirOoSTM3EJemw0jizr/TN+zg/33yerRKG6nhTlMSd0OvQlSqAvUQK7LEosSbOZW2vWEP7xJwT36EnJ8eNxfKLZQ79nLa9alHQoyargVQ+dKKIXLuTqZ59jKFkS/+9nYFO+/H1tTGYTS88tZerBqdxMukl3l0BeO7YBz3ru4ODDkBpDmHJgCk9//BHOw8dyZdTblPpyJGLHFDi20DKuUm8Q1HkWrp9FXNyB26Vd9IjeRg+DxCh0HDteju1HAljtVQ/H8o3wubCEGU4rMW0KZ1lpfzr7NEIbVQV0s3Boek+hwlIN4KmxsOYd2D4Znngrw8utS7dmYouJnIs+h4uNi+VhcMHl4m6cyn5D0nF3VjAcr2f+x5qjV/l0+XGSkh0f6u+ZHXVFoShKjiVfCCbsjTdIPnsWj5f/h9fw4blaACm9bw98y2/HfmNVz1X4OfrleD9pNBIxbjw3f/8dh6ZN8Zs0Ea3L/bN99oXvY9zecZy6cYo63nV4t8G7VHOtDD+1grgIGL6XZL0t3ZZ2w0HvwE8RTYmY+BOe1WLxCjJb7sNoPNwy7fdeiTctVwOXdmIM3oG4cgCtvDud1exbh58qNWDa5TUs6boEzfMj0bq7U2b2rExOSMKiwXB8MTy3BMq3zP4PkBQD39UD19Jc2lWWpDNnqLhxIxqDgfhkIw42D//9P6srCqSUxe5Rt25dqShK/jAlJMiwDz6QJypXkSHPvyBTIyMf6jhXYq/I+n/Ul0PWDZFmszlH+xhv3pQXBw6UJypXkeFffS3Nqan3tQmNDZVvbn5TVp9VXbZd0FauubAm4/FD90n5qYuUK9+SUkq5dt80WX1WdbngG18Z1rGCPFGliozdsDp3J5OSKFMvbJNnFn8hI/YulkmpibL5vOby5fUvy+TLofJE5Sry+m+/Zb1/UqyU0xpIOa68lNGh2b/Xmvcs8YcdkHHbt8sTlavImwsW5C7eLAD7ZCafqeqKQlGUhxK9eAnhn3+OxskRvwkTcWjYIEf7mW7dInbTJmLXriP64B4u2iXgEVCDijWaYyhXDkPZshjKlr1TC+m25HPnuDzsVYxXr+IzejSuvXred+ytoVt5c/ObaISGl2q8xIvVXsROZ3d/EKvfgT0zoXQj5KWdvOjnR4itAyvbzePa4DcxXr9OuSWLH3osZvHZxXy641N+avcTARvPE/HFF1RYuwZD2bJZ73TtjOVqxzsQXlxtGeS+V8Rx+OEJy13wXaYgpSS4Vy9kYhLlV63MtEx5bqjBbEVR8lzS6TOEvfEGKRcv4vXaa3gMHZLph5Xp1i1iN24idu1a4nbsgNRUdL6+2DdswOFT/+IYHoNXDBkK+Om8vO4kDZ2XFzdmz0bY2eE/dWqms5Mux17m6ZVP4+/oz9TWU/FxyOZDPukWfN/EstBR09c4VrYh/dcPYnCNwbzs2pWQ3r2xCQigzO9zcj3LyyzN9FjWA4PWwPzO87k8eAipYWFUWJvNWhy3HV8KC16ABv+DjuMzviYlzOoMkcdhxAGwtwy6x6xcxZVRo/CfMR2n1q1zFeu91GC2oih5zrZyAGUXLCD8k0+4NmUKCfv3U3L8OHRubphiYojduIlb69YSv2OnJTmU9MX92Wdxbv8UtjVrWpJKXBg9lvWgrmtNJld6h5SQEFJCLlp+BgcTu2EDpps3sa1aFf/p09D73r+KXrIpmbe2WAaCJ7WclH2SALB1hmG7LIsd6WyoDnQu35k5x+fQu0dvfL/8grA3RxI5YSIl3s9dWYxtYdu4EHOBsc3GYo5PIGHPHtwGDHjwjgDVukPocMvaHKUaQI3ed187tgguboPOU+4kCQDn9k9xbdIkon759ZETRVZUolAU5ZFoHR0oOXEC9g3qE/HlWIJ79MSmckDG5DBgwN3kcM89Dn6OfrxV9y2+2P0FKyodpk+7Pve9hyk2Fo2DQ5ZdK+P3jOfkjZNMbTU1+6VW07PJODvo9aDX2XBxA1P2T+GbDt+QsP8AN2bPxi4oCOen2uXsmMDs47Pxtvemfbn2xG/cgkxNva9abLbajoaw/bB8BJSoZumKSo61rAfuW9vS7ZSO0Olwf/FFIsaOJeHgwby7yz4dVetJUZRHJoTArV8/ysz7C42DAylnz+E+YABl5/9NxY0bKfHuO9jVqpXl3d59KvehoU9DJuydwJW4+xcN0jo5ZZkkVl5Yyfwz8xlYfSCtSrd66HPwcfDhxeovsjZkLYciD1HinbexrVWTqx9+SEpISI6OcSLqBHvC9zAgcAB6jZ64zVvQODllfiNfVrR66DMLDI7w9wBLN9m/4yH2qmX1Ps39s8xce/VE4+LCjV+zrtz7KFSiUBQlz9hVq0aFVSupuCktOWRyBZEZjdDwWdPPkEhG7xid9XKp9zgffZ7Pd35OkHcQr9V57VHDZ2C1gXjbeTN+73ikXof/5MkIrZbQ19/AnJh5RdzbTGYTM4/MxEHvQO+A3kizmbh//8XxiSdyfze7k48lWdwItpQe2TUD6jwH/vfPXAXQODjg9kx/YjdsJDk4OHfvlQMqUSiKUijc7oLaeXUni84uemD7hNQERm4ZiZ3Ojm9afINO8+g96fZ6e14Leo2j14+yOng1+pIlKTnhG5LPnOHqhx9mmcAu3rrIC2tfYOOljbxQ9QWcDE4kHTuGKSoq07UncqRsU0s3VMh/lrXA247Otrn7gAGUHD8+0xLmj0olCkVRCo07XVD7JmS6eM9tUko+2/kZIbdCGN98PN723nkWQ5cKXQh0D2TK/ikkGhNxfOIJvN8aya3Va4j6cWaGtmZp5q9Tf9FnRR+CY4IZ98Q4Xq71MgCxmzeDRnNnkaiH0mQEtHgXeswEB89sm+o8PHDp0jlPanHdSyUKRVEKjdtdUGZp5tMdn2b5DX7BmQWsDl7Nq7VfpaFvwzyP4Z367xCREMGc43MAcB80COcuXbg2ZQqxGzcCEB4fzv/W/4+xu8cSVCKIJd2W0LF8xztdbXFb/sUuqA46N7cs3+uBhIBWH0Dl9o98Xo9CJQpFUQqVB3VBHb9+nK/3fE0zv2YMrjE4X2Ko51OPtqXb8suxX4hMiEQIge+Yz7GtUYMrb7/Dmg0/0mNZDw5fO8wnjT/h+zbfZ7iqSQ0PJ/nkSZxyM9upEFOJQlGUQierLqiY5Bje+vctPOw8+KrZV2hE/n2Ejaw7EqPZyHcHvwNAY2uLw4TPidWbsPtoCjUN5VjUdRF9AvrcN2Aft2ULQO6mxRZiKlEoilLoZNYFZZZmPtr2EREJEUxoMQFXW9d8jaGUcymeDXyWZeeWcSLqBBsubqD3zqGM7wFe8VreX6HH3zbzG/viNm9BX6oUhgp5vA64lahEoShKoXRvF9Ss47PYErqFUfVGUcsri3roeWxozaG42rjyv/X/480tb+Lj4MNXLy/C/4svSdy9h4ivvr5vH3NiIvG7duHYsuXDrRJYCKk7sxVFKbT6VO7D+ovrGb93PCmmFNqVacczVZ4psPd3MjjxZt03+Xzn57xS6xWG1ByCXqOHbhVIOn2GG7/+ik1AAG79nr6zT/zOXcjkZJxatSywOPObShSKohRat7ugeizrgb+TP581+azAv6X3qNSDjuU7YqO1ybDd+62RJJ89S/gXX2BToTz29esDlvEJjb099vUyvzmuKFJdT4qiFGp+jn7M7zyfOR3m4GjI+9XbcuLeJAEgtFr8Jk7AUKoUoa+9TkpoGFJK4rZswaFZM4QhkzLhRZRKFIqiFHplXcribuv+4IYFTOvsjP+M6UijkdBXXyVx3z6MkZE4tnr4mlOFkUoUiqIoj8CmXDn8Jk0i+exZLr86HITAsfkT1g4rT1k9UQgh2gshTgshzgkh7iv8LoSwEUL8nfb6biFEWSuEqSiKkiXHJ5rh/c7bmG/dwq5mTXQeHtYOKU9ZdTBbCKEFpgNPAqHAXiHEcinliXTNBgE3pZQVhRD9gHHA0/cfTVEUxXrcX3gBjEZsAgOtHUqes/aspwbAOSnlBQAhxDygG5A+UXQDRqf9vhCYJoQQsjiu4aooSpElhMBjcP6UFLE2a3c9+QGX0z0PTduWaRsppRGIAe67rhNCDBVC7BNC7Lt27Vo+hasoivL4sXaiyDNSyplSynpSynpeXl7WDkdRFKXYsHaiCANKpXvun7Yt0zZCCB3gAkQVSHSKoiiK1RPFXqCSEKKcEMIA9AOW39NmOfBC2u+9gU1qfEJRFKXgWHUwW0ppFEIMB9YBWuBXKeVxIcTnwD4p5XLgF+B3IcQ54AaWZKIoiqIUEGvPekJKuRpYfc+2T9L9ngT0Kei4FEVRFAtrdz0piqIohZxKFIqiKEq2RHEcFxZCXAMuWjuOXPIErls7iAKmzvnxoM656Cgjpbzv/oJimSiKIiHEPill8SlgnwPqnB8P6pyLPtX1pCiKomRLJQpFURQlWypRFB4zrR2AFahzfjyocy7i1BiFoiiKki11RaEoiqJkSyUKRVEUJVsqURSwBy39mtamrxDihBDiuBDiz4KOMa/lYLnb0kKIzUKIg0KII0KIjtaIM68IIX4VQkQKIY5l8boQQkxN+3scEUIEFXSMeS0H5/xs2rkeFULsEELUKugY89qDzjldu/pCCKMQondBxZbnpJTqUUAPLIUPzwPlAQNwGKh6T5tKwEHALe25t7XjLoBzngm8kvZ7VSDE2nE/4jk3B4KAY1m83hFYAwigEbDb2jEXwDk3SfdvusPjcM5pbbTAJiz17HpbO+aHfagrioJ1Z+lXKWUKcHvp1/SGANOllDcBpJSRBRxjXsvJOUvAOe13F+BKAcaX56SUW7FUOs5KN2COtNgFuAohfAsmuvzxoHOWUu64/W8a2IVl7ZkiLQf/nQFGAIuAIv3/sUoUBSsnS78GAAFCiO1CiF1CiPYFFl3+yMk5jwYGCCFCsXzzGlEwoVlNTv4mxdkgLFdUxZoQwg/oAXxv7VgelUoUhY8OS/dTS6A/8JMQwtWaARWA/sAsKaU/lm6Z34UQ6t9mMSSEaIUlUbxr7VgKwBTgXSml2dqBPCqrr0fxmMnJ0q+hWPpvU4FgIcQZLIljb8GEmOdycs6DgPYAUsqdQghbLEXVivTlejZy8jcpdoQQNYGfgQ5SysdhOeN6wDwhBFj+PXcUQhillEutGtVDUN/aClZOln5diuVqAiGEJ5auqAsFGGNey8k5XwLaAAghAgFb4FqBRlmwlgPPp81+agTESCmvWjuo/CSEKA0sBp6TUp6xdjwFQUpZTkpZVkpZFlgIDCuKSQLUFUWBkjlb+nUd0E4IcQIwAW8X5W9fOTznt7B0sb2JZWD7RZk2ZaQoEkL8hSXZe6aNu3wK6AGklD9gGYfpCJwDEoCB1ok07+TgnD8BPIAZad+wjbKIV1fNwTkXG6qEh6IoipIt1fWkKIqiZEslCkVRFCVbKlEoiqIo2VKJQlEURcmWShSKoihKtlSiUIqktAq7L1o7DmsQQrQTQrxh7TiUx4dKFEpR1Rd40dpBWEk74A1rB6E8PlSiUJRCQAhh9zi+t1I0qEShFDlCiFlAL6CFEEKmPUanvdZNCLFPCJEkhAgXQowXQujT7TtaCHFdCNEwrV2iEGJbWokRbyHEUiFEnBDipBCi9T3vGyKEmCCE+Djt2HFCiLlCCJd72rkLIWYKISLS4tghhGh4TxsphBgphJgihLgGHE3b3kkIsT5tQZxbaRWE26WPH8ud7GXSnfustNe2CCEW3vM+LdPaVE97Xjbt+bNCiDlCiGhgRU7jVh5PqoSHUhSNAUoDrsCwtG2hQoi+wF/Aj8AHQAXgKyxfiEal298ey2JJ44F4YCrwO5CMpfz1DOAdYIEQopSUMiHdvv2xlN4YAvimHeNnoA+AEMIG2JAW29tYChu+AmwQQlSSUoanO9bbwFbgOe5+aSuH5YN7AmDGssjPGiFEcynl9rT3qgS0xlLCGh6uLtYELLWX+gCmXMatPG6svXKSeqjHwzywFFnbku65AC4Cv93T7iUgEfBIez4aSz2pFunaDEvb9km6bVXTtnVIty0Ey0I1jum2PYvlAz0w7fkgIAWolK6NDssqf9+k2yaBAw84R03avuuw1Mi6vX0CmawCCGwBFt6zrWXae1VPe1427fmSe9rlKG71eDwfqutJKS4CsFxlzBdC6G4/sCxDaQtUT9c2Bfgv3fNzaT83ZbLt3gWF1ksp49I9X4IlSdVPe94W2I+lRPztGAD+xVJ2Or3V956EEMJfCDFbCBEGGIFULIPXAZmc86NYdc/z3MStPGZU15NSXHim/bzvwzdN+vUfYmXGxWRS0n5G394gpUxJq3Jqe89xMqyRIaVMEELEYemGuh1HIywf8Pc6f8/ziPRP0hZrWg44Yam2eg5L19jngHdmJ/UIIu55npu4lceMShRKcXF77eKhwMFMXg/Oo/fJ8IEthLAHHIHb60ncAPZh6d+/V/I9z+8t3VwRqIOlu2ttuvfI6aykJMBwzza3LNre+965iVt5zKhEoRRVKWT8tn8ayypxZaWUP+Xj+z4phHBM1/3UA8uH7r605xuxdBVdklLmdoW+2wnhzgezEKIM0BQ4kq7dved+WyjQ/J5t7TJpl5lHiVsp5lSiUIqqU0A3IUR3LB+QV7BMG/1dCOGMZfZSClAe6A70lhlnLz2sRGCVEOIbLN1N32AZGD6R9voc4GVgixBiApbVCT2ABkC4lHLyA84pFJgohPgYSxfUZ9y/TOopoETanenHgOtSyhAs4yWDhBCTsYxBtCJtidkceJS4lWJOJQqlqJqBpZvmVyzdK59JKUcLIW5hmRr7EpYVAi8AK7k7DvGo5gGxwC9YupyWk667RkqZJIRohWVc4TOgBJZxjT3cvwRsBlLKZCFET2A6llldocCXWGYupR+Mn48lCYwHvIDZWFYFXCWE+ADLLK7BwDLg9bSf2XqUuJXiT61wpyg5JIQIwTL9dNSD2ipKcaKmxyqKoijZUolCURRFyZbqelIURVGypa4oFEVRlGypRKEoiqJkSyUKRVEUJVsqUSiKoijZUolCURRFydb/AYIHBNt70k6bAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# visualize results\n", "for n,u4_beta in zip(n_list,u4_list_n):\n", " plt.plot(temp_list, np.array(u4_beta), label='n={}'.format(n))\n", "\n", "plt.legend()\n", "plt.ylabel('$U_4$', fontsize=15)\n", "plt.xlabel('temperature', fontsize=15)\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "There is variation in the data due to insufficient statistics. However, we can see that the four system sizes are intersected at a single point at a temperature near 1.0, which is roughly the phase transition point. The estimation of the phase transition point with the Binder cumulant is a popular method used at the forefront of numerical analysis.\n", "\n", "> Of course, academic studies need to be done diligently, not only obtain sufficient statistics, but also to evaluate these errors (computation of error bars). As this calculation is limited to and overview, accurate error evaluation and more is omitted.\n", "\n", "## Conclusion\n", "\n", "We introduced a method of MonteCarlo sampling with annealing. We show an example of a phase transition in statistical mechanics. OpenJij can be applied in a variety of ways, depending on our ideas." ] } ], "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.9.3" } }, "nbformat": 4, "nbformat_minor": 4 }