{ "cells": [ { "cell_type": "markdown", "id": "311704d6-fa8c-42f2-bb05-4cd0bf202d4f", "metadata": {}, "source": [ "# Symmetrization of correlations on the lattice" ] }, { "cell_type": "markdown", "id": "e5a41c27-92c1-4048-8466-7e8f86ef20f3", "metadata": {}, "source": [ "The pyALF analysis offers an option to symmetrize correlation functions, by averaging over a list of symmetry operations on the Bravais lattice. This feature is meant to be used as an improved estimator, meaning to explicitly restore a symmetry of the model lost due to imperfect sampling, to increase the quality of the data.\n", "\n", "For this feature, the user has to supply a list of functions $f_i$, taking as arguments an instance of {class}`py_alf.Lattice` and an integer corresponding to a $\\boldsymbol{k}$-point of the Bravais lattice and returning an integer corresponding to the transformed $\\boldsymbol{k}$-point of the Bravais lattice. The analysis then averages the correlation over all transformations:" ] }, { "cell_type": "markdown", "id": "d027e2c6-d2b4-47d1-a619-2beb7ecdfc37", "metadata": {}, "source": [ "$$\n", "\\tilde{C}(n_{\\boldsymbol{k}}) = \\frac{1}{N} \\sum_{i=1}^{N} C \\left( f_i(latt, n_{\\boldsymbol{k}}) \\right)\n", "$$" ] }, { "cell_type": "markdown", "id": "4fcae61c-4261-4156-a5fa-db766d6e8cf6", "metadata": {}, "source": [ "```{note}\n", "This symmetrization feature does not affect custom observables, but only the default analysis. Improved estimators would have to be included directly in the definition of custom observables.\n", "```" ] }, { "cell_type": "markdown", "id": "35fcd7b6-01a7-4f17-ac9d-fd7351615635", "metadata": {}, "source": [ "This demonstration begins, as usual, with some imports:" ] }, { "cell_type": "code", "execution_count": 1, "id": "ef359177-7cbb-4e43-ab52-28a3ffc1a0c6", "metadata": {}, "outputs": [], "source": [ "# Enable Matplotlib Jupyter Widget Backend\n", "%matplotlib widget\n", "\n", "import numpy as np # Numerical libary\n", "import matplotlib.pyplot as plt # Plotting library\n", "from py_alf.analysis import analysis # Analysis function\n", "from py_alf.ana import load_res # Function for loading analysis results\n", "from py_alf import Lattice # Defines Bravais lattice object\n", "from custom_obs import custom_obs # Custom observable specifications\n", " # from local file custom_obs.py " ] }, { "cell_type": "markdown", "id": "61db120d-49cc-420a-a9fe-8364a20e8f11", "metadata": {}, "source": [ "The Hubbard model on a square lattice possesses a fourfold rotation symmetry ($=C_4$ symmetry). To restore this symmetry, a list of all possible realizations of it has to be handed to the analysis. These are: rotation by 0 or $2\\pi$ ($=\\text{identity}$), rotation by $\\pi/2$, rotation by $\\pi$ and rotation by $3\\pi/2$." ] }, { "cell_type": "code", "execution_count": 2, "id": "39300200-101c-4a52-827c-bba402200be4", "metadata": {}, "outputs": [], "source": [ "# Define list of transformations (Lattice, i) -> new_i\n", "# Default analysis will average over all listed elements\n", "def sym_c4_0(latt, i): return i\n", "def sym_c4_1(latt, i): return latt.rotate(i, np.pi*0.5)\n", "def sym_c4_2(latt, i): return latt.rotate(i, np.pi)\n", "def sym_c4_3(latt, i): return latt.rotate(i, np.pi*1.5)\n", "\n", "sym_c4 = [sym_c4_0, sym_c4_1, sym_c4_2, sym_c4_3]" ] }, { "cell_type": "markdown", "id": "2ca46eac-4745-4ac5-86c6-3393a6e0fce7", "metadata": {}, "source": [ "Set directory to be analyzed." ] }, { "cell_type": "code", "execution_count": 3, "id": "3751b4b0-c4eb-49e6-bc61-41a5077253da", "metadata": {}, "outputs": [], "source": [ "directory = './ALF_data/Hubbard'" ] }, { "cell_type": "markdown", "id": "0655443e-3974-4e17-9f49-02b60b2f6b8e", "metadata": {}, "source": [ "Analyzed without symmetrization and load results." ] }, { "cell_type": "code", "execution_count": 4, "id": "d9c8d92f-8cdc-4a3b-91af-b2b05968615d", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "### Analyzing ./ALF_data/Hubbard ###\n", "/scratch/pyalf-docu/doc/source/usage\n", "Custom observables:\n", "custom E_squared ['Ener_scal']\n", "custom E_pot_kin ['Pot_scal', 'Kin_scal']\n", "custom R_Ferro ['SpinT_eq']\n", "custom R_AFM ['SpinT_eq']\n", "custom SpinZ_pipi ['SpinZ_eq']\n", "custom SpinXY_pipi ['SpinXY_eq']\n", "custom SpinXYZ_pipi ['SpinT_eq']\n", "Scalar observables:\n", "Ener_scal\n", "Kin_scal\n", "Part_scal\n", "Pot_scal\n", "Histogram observables:\n", "Equal time observables:\n", "Den_eq\n", "Green_eq\n", "SpinT_eq\n", "SpinXY_eq\n", "SpinZ_eq\n", "Time displaced observables:\n", "Den_tau\n", "Green_tau\n", "SpinT_tau\n", "SpinXY_tau\n", "SpinZ_tau\n", "./ALF_data/Hubbard\n" ] } ], "source": [ "analysis(directory, symmetry=None, custom_obs=custom_obs, always=True)\n", "res_nosym = load_res([directory]).iloc[0]" ] }, { "cell_type": "markdown", "id": "56fa65c4-9d1c-4445-9930-9d9ae2e66e51", "metadata": {}, "source": [ "Analyze with symmetrization and load results." ] }, { "cell_type": "code", "execution_count": 5, "id": "7262eb4e-4144-48be-b4ab-df773e9044fd", "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "### Analyzing ./ALF_data/Hubbard ###\n", "/scratch/pyalf-docu/doc/source/usage\n", "Custom observables:\n", "custom E_squared ['Ener_scal']\n", "custom E_pot_kin ['Pot_scal', 'Kin_scal']\n", "custom R_Ferro ['SpinT_eq']\n", "custom R_AFM ['SpinT_eq']\n", "custom SpinZ_pipi ['SpinZ_eq']\n", "custom SpinXY_pipi ['SpinXY_eq']\n", "custom SpinXYZ_pipi ['SpinT_eq']\n", "Scalar observables:\n", "Ener_scal\n", "Kin_scal\n", "Part_scal\n", "Pot_scal\n", "Histogram observables:\n", "Equal time observables:\n", "Den_eq\n", "Green_eq\n", "SpinT_eq\n", "SpinXY_eq\n", "SpinZ_eq\n", "Time displaced observables:\n", "Den_tau\n", "Green_tau\n", "SpinT_tau\n", "SpinXY_tau\n", "SpinZ_tau\n", "./ALF_data/Hubbard\n" ] } ], "source": [ "analysis(directory, symmetry=sym_c4, custom_obs=custom_obs, always=True)\n", "res_sym = load_res([directory]).iloc[0]" ] }, { "cell_type": "markdown", "id": "6ffca4fd-7a3c-49ae-a722-e146ec45763d", "metadata": {}, "source": [ "We now compare results for the points $(\\pi, \\pi) + \\boldsymbol{b}_1$ and $(\\pi, \\pi) + \\boldsymbol{b}_2$, where $\\boldsymbol{b}_1= (2\\pi/L, 0)$ and $\\boldsymbol{b}_2= (0, 2\\pi/L)$ are the primitive vectors of the Bravais lattice in k-space, with and without symmetrization." ] }, { "cell_type": "code", "execution_count": 6, "id": "3077c5b4-eb8f-41d3-81f3-3745bd45808e", "metadata": {}, "outputs": [], "source": [ "latt = Lattice(res_nosym['SpinT_eq_lattice'])\n", "n = latt.k_to_n((np.pi, np.pi))\n", "n1 = latt.nnlistk[n, -1, 0]\n", "n2 = latt.nnlistk[n, 0, -1]" ] }, { "cell_type": "code", "execution_count": 7, "id": "da66366e-c5fd-4c5a-9e47-049bd1333471", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1.1329733101215775 0.011981792090571623\n", "0.9815172263634666 0.0945777371428819\n" ] } ], "source": [ "print(res_nosym.SpinT_eqK_sum[n1], res_nosym.SpinT_eqK_sum_err[n1])\n", "print(res_nosym.SpinT_eqK_sum[n2], res_nosym.SpinT_eqK_sum_err[n2])" ] }, { "cell_type": "code", "execution_count": 8, "id": "c2b7aa9c-0212-4ebf-8914-61b72824d77a", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1.057245268242522 0.04896785403524968\n", "1.057245268242522 0.04896785403524977\n" ] } ], "source": [ "print(res_sym.SpinT_eqK_sum[n1], res_sym.SpinT_eqK_sum_err[n1])\n", "print(res_sym.SpinT_eqK_sum[n2], res_sym.SpinT_eqK_sum_err[n2])" ] }, { "cell_type": "code", "execution_count": 9, "id": "6c66620e-a2fa-4290-a06d-1d985e5f8f4f", "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "7b09b669e7f249d0ba34ea464f67f8d4", "version_major": 2, "version_minor": 0 }, "image/png": "iVBORw0KGgoAAAANSUhEUgAAAoAAAAHgCAYAAAA10dzkAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dD1DU953/8fcXkD9BQUURHDFqbPFfNR6aiFWj0WjRs5qaTHuTMcbE3JgzicaaGuxMNL1a0juuJUaDkhqtpWrawT9pNZ7ORSA540UUGy+JTpPYQghqMRWQUxDYH5/PL8uA7rKA8P0u+3k+Zj6ju3y/u5/dBfbF+/NnLVcDAQAAgDGCnO4AAAAA7EUABAAAMEyI0x1A11NfXy9ffvml9OjRQyzLcro7AAA0UjPbKisrpX///hIURJ3LGwIg2kyFv4SEBJ45AIDfKi4ulgEDBjjdDb9FAESbqcqf+4crKiqKZxAA4DcqKip0kcL9XgXPCIBoM/ewrwp/BEAAgD9iilLLGBwHAAAwDAEQAADAMARAAAAAwxAAAQAADEMANFxaWpqeKLtixQqnuwIAAGxCADTYiRMnJCsrS0aPHu10VwAAgI0IgIa6evWqPPLII/L6669Lr169nO4OAACwEQHQUMuWLZM5c+bIjBkzfB5bXV2tN9Zs2gAAQNfFRtAG2r17t5w6dUoPAbd2nuBLL73Uyb0CAAB2oQJoGPXxbcuXL5fs7GwJDw9v1TmpqalSXl7e2NRtAACArstyNXC6E7DPvn375MEHH5Tg4ODG6+rq6vRK4KCgID3c2/Rrnqgh4OjoaB0G+Sg4AIA/4T2qdRgCNsz06dPlzJkzza5bvHixDBs2TFavXu0z/AEAgK6PAGiYHj16yKhRo5pdFxkZKTExMbdcDwAAAhNzAAEAAAxDBRCSm5vLswAAgEGoAAIAABiGAAgAAGAYAiAAAIBhCIAAAACGIQACAAAYhgAIAABgGAIgAACAYQiAAAAAhiEAAgAAGIYACAAAYBgCIAAAgGEIgAAAAIYhAAIAABgmxOkOAG31xRdfGPWk9e7d2+ku2M7lcjndBdtZluV0F2z11VdfOd0F2yUkJDjdBduZ+LPcVVABBAAAMAwBEAAAwDAEQAAAAMMQAAEAAAxDAAQAADAMARAAAMAwBEAAAADDEAABAAAMQwAEAAAwDAEQAADAMARAAAAAwxAAAQAADEMABAAAMAwBEAAAwDAEQAAAAMMQAAEAAAxDAAQAADAMARAAAMAwBEAAAADDEAABAAAMQwAEAAAwDAEQAADAMARAAAAAwxAAAQAADEMABAAAMAwBEAAAwDAEQAAAAMMQAAEAAAxDAAQAADAMARAAAMAwBEAAAADDEAABAAAMQwAEAAAwDAHQQJmZmTJ69GiJiorSLTk5Wd5++22nuwUAAGxCADTQgAED5OWXX5aCggLd7r//fpk3b5589NFHTncNAADYIMSG+4CfmTt3brPL69ev11XB48ePy8iRIx3qFQAAsAsB0HB1dXXy+9//XqqqqvRQsCfV1dW6uVVUVNjVPQAA0AkYAjbUmTNnpHv37hIWFiZLly6VvXv3yogRIzwem5aWJtHR0Y0tISHB5t4CAICORAA0VGJiopw+fVoP+z711FOyaNEi+fjjjz0em5qaKuXl5Y2tuLjY5t4CAICOxBCwoUJDQ2Xo0KH6/+PGjZMTJ07IK6+8Ilu2bLnlWFUlVA0AAAQGKoDQXC5Xs3l+AAAgcFEBNNCaNWskJSVFz+WrrKyU3bt3S25urhw6dMjprgEAABsQAA108eJFWbhwoZSWlupFHWpTaBX+HnjgAae7BgAAbEAANNDWrVud7gIAAHAQcwABAAAMQwAEAAAwDAEQAADAMARAAAAAwxAAAQAADEMABAAAMAwBEAAAwDAEQAAAAMOwETQAAEALrl+/LjU1NS0c4VtoaKiEh4ff1m10JAIgAABAC+EvIiLCy1dbLy4uTs6fP+83IZAhYAAAAC9ut/LnduHChQ67rY5ABRAAAKAVLMtqxVG3crlc7TqvMxEAAQAAWhH+2hsA/TEEEgABAAB8IAACAAAYJigo6LaGgOvr6zu4R7eHCiC6HH9ZQWUXfxs2sENIiHm/mmpra53ugq0iIyOd7oLtLl686HQX4GAF0N+Y91sWAACgjQiAAAAAhrECrALIPoAAAACtDIDtbe2Vlpamz1+xYoXXY3Jzcz3e59mzZ72ewxAwAACAH1YAT5w4IVlZWTJ69OhWHX/u3DmJiopqvNy3b1+vx1IBBAAAaMUq4NtpbXX16lV55JFH5PXXX5devXq16pzY2Fj9kXPuFhwc7P3xtLlHAAAAhrE6YAi4oqKiWauurvZ6f8uWLZM5c+bIjBkzWt3HsWPHSnx8vEyfPl2OHj3a4rEEQAAAABsCYEJCgkRHRzc2Nb/Pk927d8upU6e8fv1mKvSpoeKcnBzZs2ePJCYm6hCYn5/v9RzmAAIAALQyAN6O4uLiZnP0wsLCPB6zfPlyOXz4cKv3vVWBTzW35ORkfTvp6ekyZcoUj+dQAQQAALChAqjCX9PmKQCePHlSLl26JElJSXpTfNXy8vJkw4YN+v91dXU+evr/TZgwQf785z97/ToVQAAAgFYGwM6mhm7PnDnT7LrFixfLsGHDZPXq1S0u7GiqsLBQDw17QwAEAADwQYW/9qzmVdryOcA9evSQUaNG3fLRiTExMY3Xp6amSklJiezYsUNfzsjIkEGDBsnIkSOlpqZGsrOz9XxA1bwhAAIAAHRiBbCjK4elpaVSVFTUeFmFvlWrVulQGBERoYPggQMHZPbs2d775DLxk+ZxW9TSdbV6qby8vNlkVruUlZXZfp9OUj/MplHzXExTW1vrdBdspd6wTHPjxg2nu2A7tS9dV3+Pqvj69lQF7nYqgJcvX3bsfdMT837LAgAAdOEKYEcgAAIAAPhAAAQAADCMRQUQAADALEHt/Exff8UQMAAAgA9UAAEAAAxjMQQMAABgFosACAAAYBaLAAgAAGAWiwAIAABgliBWAQMAAJjFogIIAABgFosACAAAYBaLAAgAAGBmCAwUfBIIAABAJy4Ccblc7TqvMxEAAQAAfGAIGAAAwDBWgM0BbF8tEwAAwMAh4KB2tvZKS0vTAXLFihUtHpeXlydJSUkSHh4uQ4YMkc2bN7f8eNrdI3RZ6ptp/Pjx0qNHD4mNjZX58+fLuXPnnO4WAAB+XwG02tna48SJE5KVlSWjR49u8bjz58/L7NmzZfLkyVJYWChr1qyRZ599VnJycryeQwA0kPorYdmyZXL8+HE5cuSI1NbWysyZM6WqqsrprgEA4JeCbK4AXr16VR555BF5/fXXpVevXi0eq6p9AwcOlIyMDBk+fLgsWbJEHn/8cUlPT/f+eNrcI3R5hw4dkscee0xGjhwpY8aMkW3btklRUZGcPHnS6a4BAOCXLJsrgKpQM2fOHJkxY4bPY99//31dyGlq1qxZUlBQIDdu3PB4DquAIeXl5fpZ6N27t8dno7q6Wje3iooKnjUAgFGsDlgEcvP7Z1hYmG432717t5w6dUoPAbfGhQsXpF+/fs2uU5fVCF9ZWZnEx8ffcg4VQMOpvYlWrlwpkyZNklGjRnmdMxgdHd3YEhISbO4lAADOCg4Ovq2mqPfPpu+n6v31ZsXFxbJ8+XLJzs7WCzpa6+Zw6t570FtopQJouKefflo+/PBDee+997wek5qaqkOim/oLhhAIADBJUAdsBK3CXVRUVOP1nqp/ajrWpUuX9Ipet7q6OsnPz5eNGzfqETl3oHSLi4vTVcCm1G2EhIRITEyMxz4RAA32zDPPyFtvvaW/qQYMGOD1OG8lagAATGF1wBCwCn9NA6An06dPlzNnzjS7bvHixTJs2DBZvXr1LeFPSU5Olj/84Q/Nrjt8+LCMGzdOunXr5vF+CIAGUn+JqPC3d+9eyc3NlcGDBzvdJQAA/JrVEOLaWwGsr69v9bFqi7abp2RFRkbqSp77ejUyV1JSIjt27NCXly5dqquDarTuySef1ItCtm7dKrt27fJ6PwRAA6mVRTt37pT9+/frbzR32VjNR4iIiHC4dwAA+GcAtPzkk0BKS0v17h1uqpBz8OBBee6552TTpk3Sv39/2bBhgyxYsMB7n1z++AnF6FTevhHVdjBqexhf1BxAFRbV6mFfpezOoFY0mcTEUK7mrZhGrdYzSU1NjdNdsJ237TgCmfqwAbt19HtUxde3d88997T7d5P6+f7ggw8ce9/0xLzfsmicjAoAADp/EUh7z+tMBEAAAIAuNATcEQiAAAAAPlABBAAAMIxFBRAAAMAsFgEQAADALEEsAgEAADCLRQUQAADALEFUAAEAAMxiUQEEAAAwi0UABAAAMC8ABrXzEz3Uuf6GjaABAAB8oAIIAABgmCAWgQAAAJjFYg4gAACAWYKoAAIAAJjFogIIAABgFosACAAAYBYrwAJg+za0AQAAMDAAWu1sbZGZmSmjR4+WqKgo3ZKTk+Xtt9/2enxubq7H+zx79qzXc9gHEF3OHXfc4XQXbBUSYt6PaWhoqNNdsJ1pr3NdXZ3TXbBdt27dnO4CboNlYwVwwIAB8vLLL8vQoUP15V//+tcyb948KSwslJEjR3o979y5czowuvXt29frsWb9xgEAAPDzADh37txml9evX6+rgsePH28xAMbGxkrPnj1bdR8MAQMAALQyALa3KRUVFc1adXV1q6rlu3fvlqqqKj0U3JKxY8dKfHy8TJ8+XY4ePdrisQRAAAAAGwJgQkKCREdHN7a0tDSv93fmzBnp3r27hIWFydKlS2Xv3r0yYsQIj8eq0JeVlSU5OTmyZ88eSUxM1CEwPz/f6+0zBAwAAGDDRtDFxcXN5uipcOeNCnGnT5+WK1eu6GC3aNEiycvL8xgC1bGqualKobqv9PR0mTJlisfbJwACAADYMAfQvaq3tYvh3ItAxo0bJydOnJBXXnlFtmzZ0qrzJ0yYINnZ2V6/TgAEAABohfYGwI7gcrlaNWfQTa0YVkPD3hAAAQAA/GgV8Jo1ayQlJUXPGaysrNSLQNRef4cOHdJfT01NlZKSEtmxY4e+nJGRIYMGDdIrhGtqanTlTw0bq+YNARAAAMCPAuDFixdl4cKFUlpaqheLqE2hVfh74IEH9NfV9UVFRY3Hq9C3atUqHQojIiJ0EDxw4IDMnj3be58aSoqudj0aGEstXVffkOXl5a2ey9CR/u///s/2+3SSaRsEm7oRtGkbI6stLUzT3gUEXZlaxdrV36Mqvr69f/qnf2r37yYV0Hbt2uXY+6Yn5r2zAAAA+HEF0A4EQAAAAB8IgAAAAIaxqAACAACYxSIAAgAAmMUiAAIAAJjFIgACAACYxSIAAgAAmMUiAAIAAJjFIgACAACYxSIAAgAAmMUiAAIAAJjFIgACAACYxSIAAgAAmMUiAAIAAJgZAgNFiNMdAAAA8HcWFUAAAACzWAEWAIOc7gAAAEBXCYBWO1tbZGZmyujRoyUqKkq35ORkefvtt1s8Jy8vT5KSkiQ8PFyGDBkimzdvbvF4AiAAAIAPdgbAAQMGyMsvvywFBQW63X///TJv3jz56KOPPB5//vx5mT17tkyePFkKCwtlzZo18uyzz0pOTo7X+yAAGig/P1/mzp0r/fv319+U+/btc7pLAAD4NcvGAKjeo1Wg++Y3v6nb+vXrpXv37nL8+HGPx6tq38CBAyUjI0OGDx8uS5Yskccff1zS09O93gcB0EBVVVUyZswY2bhxo9NdAQCgSwgKCrqtplRUVDRr1dXVPu+3rq5Odu/erd+71VCwJ++//77MnDmz2XWzZs3S1cMbN254PIdVwAZKSUnRDQAAtE57KnlNz1USEhKaXb927VpZt26dx3POnDmjA9/169d19W/v3r0yYsQIj8deuHBB+vXr1+w6dbm2tlbKysokPj7+lnMIgPBJ/YXS9K8U9VcLAAAmsTogABYXF+tFHW5hYWFez0lMTJTTp0/LlStX9Fy+RYsW6YUe3kLgzX1zuVwer3cjAMKntLQ0eemll3imAADGsjogALpX9bZGaGioDB06VP9/3LhxcuLECXnllVdky5YttxwbFxenq4BNXbp0SUJCQiQmJsbj7TMHED6lpqZKeXl5Y1N/wQAAYGIAtGxYBOKJquh5mzOohoqPHDnS7LrDhw/r4NitWzeP51ABhE+qRN1SmRoAgEBn2bgRtNrGRc3VV3MGKysr9SKQ3NxcOXToUGNhpqSkRHbs2KEvL126VC/sXLlypTz55JN6UcjWrVtl165dXu+DAAgAAOBHAfDixYuycOFCKS0tlejoaL0ptAp/DzzwgP66ur6oqKjx+MGDB8vBgwflueeek02bNult3jZs2CALFizweh8EQANdvXpVPv3002YbSKqJpr1799b7CAEAAOcCoKretWT79u23XHfffffJqVOnWn0fBEADqX2Bpk2b1nhZlYwVtcLI0zcVAACmswLss4AJgAaaOnVq4/JwAADgGwEQAADAMBYVQAAAAPMCYNDXH+nWnnP9DUPAAAAAPlABBAAAMIzFEDAAAIBZLAIgAACAWSwCIAAAgFksAiAAAIBZLAIgAACAWSwCIAAAgFksAiAAAIBZLAIgAACAWYKCgtr9SSDtPa8z8UkgAAAAPlABBAAAMIzFEDAAAIBZrAALgP43KG2IZ555RsrKypzuBgAAaEMAbG9ri7S0NBk/frz06NFDYmNjZf78+XLu3LkWz8nNzfV4v2fPnvV4PAHQIVOnTpVZs2bJT3/6U7l27ZpT3QAAAH4WAPPy8mTZsmVy/PhxOXLkiNTW1srMmTOlqqrK57kqKJaWlja2b3zjGx6PIwA6ZMGCBfI///M/Eh0dLRMmTJBf/epX4nK5nOoOAADwwY7wpxw6dEgee+wxGTlypIwZM0a2bdsmRUVFcvLkSZ/nqophXFxcYwsODvZ4HAHQQSEhIfK9731PfvjDH8oLL7wgo0aNkoMHDzrZJQAA4HAF8Gbl5eX63969e/s8duzYsRIfHy/Tp0+Xo0ePej2ObWAckpKSIp988okMGDBA7rnnHnn11Vflm9/8pmzatEmXe3/5y1861TUAANAJi0AqKiqaXR8WFqZbS9To4MqVK2XSpEm6UOSNCn1ZWVmSlJQk1dXV8pvf/EaHQDU3cMqUKbccTwC0ycWLF6Vfv36Nl3/2s5/J6NGjbynNvvHGGzJs2DACIAAAARYAExISml2/du1aWbduXYvnPv300/Lhhx/Ke++91+JxiYmJurklJydLcXGxpKenEwCdnvOnUrga9nWXaN3U5E739QrDwC1r+lyZwLTHq9TV1al/jOKP20R0pvDwcKe7YLv6+nqnuwCHPwmkuCGQRUVFNV7vq/qndgx56623JD8/X48YtpVaY5Cdne25T22+NbRLr1699At5s8uXL8uMGTOaXTdkyBCeZQAAAmwOYFRD+GvavAVANeyrKn979uyRd955RwYPHtyuPhcWFuqhYU/MKy04RI3Fq7l+arXvkiVL9HVqDuA//uM/6lU+AADAf1k2bgSttoDZuXOn7N+/X+8FeOHCBX292jkkIiJC/z81NVVKSkpkx44d+nJGRoYMGjRIZ4qamhpd+cvJydHNEwKgTXr27KlfhPvuu0++9a1vyd///nf5wQ9+IP/8z/8sP//5z+3qBgAA8PMAmJmZ2bhncFNqOxi1PYyi9vhTW8O4qdC3atUqHQpVSFRB8MCBAzJ79myP90EA7ETz5s2Tu+++W8/3U/+q4KdW+c6ZM0euX7+u/79o0aLO7AIAAPCTOYCt1Zp9gbdv397s8o9+9CPdWosA2InU7tv//d//La+99pqe66eqgGpDR/XCPvLIIzoU3rhxQ7p169aZ3QAAAF2oAmgHAmAnUkuv3b744gs5ffq0bjExMfJf//Vfej6gWuGptn3505/+1JldAQAAt8EiAKI91PJt1dSiD7erV6/qFTpqfx8AAOC/LAIgOkr37t1l8uTJugEAAP9lEQABAADMEmTjIhA7MAcQAADAByqAAAAAhobAQEEFEAAAwAcqgAAAAIYJYg4gAACAWSxWAQMAAJjFIgACAACYxSIAAgAAmMUiAAIAAJgliEUgAAAAZrGoAAIAAJjFIgACAACYxSIAAgAAmCUowOYA+l+PAAAA/LQCaLWztUVaWpqMHz9eevToIbGxsTJ//nw5d+6cz/Py8vIkKSlJwsPDZciQIbJ582avxxIAAQAAfLAzAKogt2zZMjl+/LgcOXJEamtrZebMmVJVVeX1nPPnz8vs2bNl8uTJUlhYKGvWrJFnn31WcnJyPB4f0qYeAQAAGMiycQ7goUOHml3etm2brgSePHlSpkyZ4vEcVe0bOHCgZGRk6MvDhw+XgoICSU9PlwULFtxyPBVAg7322msyePBgXSpWJeN3333X6S4BAODXAdC6jQpgRUVFs1ZdXd2q+y4vL9f/9u7d2+sx77//vq4SNjVr1iwdAm/cuHHL8QRAQ7355puyYsUK+fGPf6xLxapknJKSIkVFRU53DQAAv2M1hDj3QpC2NnWukpCQINHR0Y1NzfXzxeVyycqVK2XSpEkyatQor8dduHBB+vXr1+w6dVkNH5eVld1yPEPAhvrFL34hTzzxhCxZskRfViXj//zP/5TMzMxWfUMCAGASqwOGgIuLiyUqKqrx+rCwMJ/nPv300/Lhhx/Ke++91+r7aRoePV2vEAANVFNTo+cRvPDCC82uV6XjY8eO3XK8KlE3LVOrsjUAACaxOiAAqvDXNAD68swzz8hbb70l+fn5MmDAgBaPjYuL01XApi5duiQhISESExNzy/EMARtIlYLr6uo8lopv/uZRVEWwaclalbABADAxAFo2rAJWlTtV+duzZ4+88847er6+L8nJyXrFcFOHDx+WcePGSbdu3W45ngBosJu/IdU3nKdv0tTUVD0B1d1UCRsAAJMEtXP+n7u1hdoCJjs7W3bu3Kn3AlTFGdWuXbvW7L350Ucfbby8dOlS+etf/6rnC37yySfyxhtvyNatW2XVqlUe74MhYAP16dNHgoODPZaKb64KuucotGaeAgAAgcqycRsYNR9fmTp16i3bwTz22GP6/6Wlpc0Wbqoq4cGDB+W5556TTZs2Sf/+/WXDhg0et4BRCIAGCg0N1du+qFLxgw8+2Hi9ujxv3jwHewYAAFxfL95oyfbt22+57r777pNTp0616gkkABpKlYgXLlyo5waoeQNZWVn6LwlVQgYAAM5VAO1AADTU97//fbl8+bL85Cc/0WVktbeQKh3feeedTncNAAC/YxEAESj+5V/+RTcAANAyAiAAAIBhLCqAAAAAZrEIgAAAAGaxCIAAAABmsQIsAPJJIAAAAIZhGxgAAADDKoAEQAAAAB8IgAAAAIaxqAACAACYGQIDBUPAAAAAPlABBAAAMIzFEDAAAIBZrAALgOwDCAAAYBjmAAIAAPhABRAAAMDQAGi1s7VFfn6+zJ07V/r376/P3bdvX4vH5+bmerzPs2fPej2HCiAAAIAP7QlyTc9ti6qqKhkzZowsXrxYFixY0Orzzp07J1FRUY2X+/bt6/VYAiAAAIAfBcCUlBTd2io2NlZ69uzZqmNZBAIAAOBDe4d+byc4ttXYsWMlPj5epk+fLkePHm3xWCqAAAAANlQAKyoqml0fFham2+1SoS8rK0uSkpKkurpafvOb3+gQqOYGTpkyxeM5BEAAAAAbAmBCQkKz69euXSvr1q1r1202lZiYqJtbcnKyFBcXS3p6OgEQAADAyQBY3BDKmi7S6IjqnzcTJkyQ7Oxsr1+nAggAAGADFf6aBsDOVFhYqIeGvSEAAgAA+NEq4KtXr8qnn37aePn8+fNy+vRp6d27twwcOFBSU1OlpKREduzYob+ekZEhgwYNkpEjR0pNTY2u/OXk5OjmDQEQXc7169ed7oKtIiMjne6C7YKDg53ugu3UL22T1NbWOt0F29XX1zvdBXSRAFhQUCDTpk1rvLxy5Ur976JFi2T79u1SWloqRUVFzX5/rFq1SofCiIgIHQQPHDggs2fP9t4nV4M2Pg4YTq1iio6OlvLycttK2Tffv0kIgGYgAAY+EwNg9+7du/x7VMXXt3fixIl2Px5V0Rs/frxj75ueUAEEAADwowqgHQiAAAAAPhAAAQAADGMFWAWQj4IDAAAwDEPAAAAAXbSS114EQAAAAB8YAgYAAECXRgUQAADAsAogARAAAMAHAiAAAIBhrACrALINDAAAgGEYAgYAADCsAkgABAAAMCwAMgQMAABgGCqAAAAAhlUACYAAAAA+EAABAAAMYwVYBZA5gAAAAH4kPz9f5s6dK/3799fhcd++fT7PycvLk6SkJAkPD5chQ4bI5s2bWzyeAAgAANDKCmB7W1tUVVXJmDFjZOPGja06/vz58zJ79myZPHmyFBYWypo1a+TZZ5+VnJwcr+cwBxAAAMCPhoBTUlJ0ay1V7Rs4cKBkZGToy8OHD5eCggJJT0+XBQsWeDyHCiAAAIAfVQDb6v3335eZM2c2u27WrFk6BN64ccPjOQRAA61fv14mTpwod9xxh/Ts2dPp7gAAYISKiopmrbq6ukNu98KFC9KvX79m16nLtbW1UlZW5vEcAqCBampq5OGHH5annnrK6a4AAGBMBTAhIUGio6MbW1paWof2rymXy+XxejfmABropZde0v9u377d4Z4AAGDOHMDi4mKJiopqvD4sLKxD+hYXF6ergE1dunRJQkJCJCYmxuM5BEAAAAAbqPDXNAB2lOTkZPnDH/7Q7LrDhw/LuHHjpFu3bh7PYQgYPqk5CjfPWwAAwCSWjYtArl69KqdPn9bNvc2L+n9RUZG+nJqaKo8++mjj8UuXLpW//vWvsnLlSvnkk0/kjTfekK1bt8qqVau83gcBMECsW7fO5zefWg3UHmqOQtM5C2oOAwAAJrFsDIDq/Xrs2LG6KSrYqf+/+OKL+nJpaWljGFQGDx4sBw8elNzcXLn77rvlX//1X2XDhg1et4DRj8flniWILk2t8vG20sdt0KBBeodwNzUHcMWKFXLlyhWfFcCmK5VUBVCFwPLy8k4pZftiWgUyMkcvEcEAABEHSURBVDLS6S7YLjg42OkuOLI4yyRqdaJp6uvrne6C7bp37+7Ie4QqVnTUe1TF17enqnA9evRo121UVlbqkObU+6YnzAEMEH369NGtM6hJqh01URUAADiPAGggVTb+6quv9L91dXWNcwyGDh3qyF9rAAB0BVYnb+hsJwKggdQcgl//+teNl91zDI4ePSpTp051qlsAAPh1+LNs+ig4O7AIxEBq7p+a+nlzI/wBAGAGKoAAAACGVQAJgAAAAIYFQIaAAQAADEMFEAAAwAcqgAAAAOjSqAACAAD4QAUQAAAAXRoVQAAAAB+oAAIAAKBLowIIAADgAxVAAAAAdGlUAAEAAHygAggAAIAujQogAACAD1QAAQAA0Olee+01GTx4sISHh0tSUpK8++67Xo/Nzc1tDKlN29mzZz0eTwUQAADAzyqAb775pqxYsUKHwG9/+9uyZcsWSUlJkY8//lgGDhzo9bxz585JVFRU4+W+fft6PC6ozT0CAABAp/rFL34hTzzxhCxZskSGDx8uGRkZkpCQIJmZmS2eFxsbK3FxcY0tODjY43EEQAAAABtUVFQ0a9XV1R6Pq6mpkZMnT8rMmTObXa8uHzt2rMX7GDt2rMTHx8v06dPl6NGjXo8jAAIAAPjgaX5dW5qiKnjR0dGNLS0tzeN9lZWVSV1dnfTr16/Z9eryhQsXPJ6jQl9WVpbk5OTInj17JDExUYfA/Px8j8czBxAAAMAGxcXFzebnhYWFtWnuoMvl8jqfUAU+1dySk5P1/aWnp8uUKVNuOZ4KIAAAgA0VQBX+mjZvAbBPnz567t7N1b5Lly7dUhVsyYQJE+TPf/6zx68RAAEAAPxIaGio3vblyJEjza5XlydOnNjq2yksLNRDw54wBAwAAOBn28CsXLlSFi5cKOPGjdPDuWp+X1FRkSxdulR/PTU1VUpKSmTHjh36slolPGjQIBk5cqReRJKdna3nA6rmCQEQAADAz3z/+9+Xy5cvy09+8hMpLS2VUaNGycGDB+XOO+/UX1fXqUDopkLfqlWrdCiMiIjQQfDAgQMye/Zsj7dvNUwodNnySBAw1NJ1tXqpvLy82WRWO+/fJJGRkU53wXbe9q0KZOqXt0lqa2ud7oLt6uvrne6C7bp3797l36Mqvr49Fcbae3vqNmJiYhx73/SEOYAAAACGYQgYAADAz+YAdjYCILoctTmmSa5fv+50F2wXEmLerybThkSvXbvmdBdsxxAw/Il5v2UBAAAMrwAyBxAAAMAwVAABAAB8oAIIAACALo0KIAAAgA+VlZXtnsunzvU3BEAAAIAWPpc3Li5OEhISvBzROuo21G35CwIgAACAF+Hh4XL+/Pnb/rQeFf7UbfkLAiAAAEALVHDzp/DWEdgGBgAAwDAEQAAAAMMQAAEAAAxDAAQAADAMARAAAMAwBEAAAADDEAABAAAMQwAEAAAwDAEQAADAMARAAAAAwxAAAQAADEMANMxf/vIXeeKJJ2Tw4MESEREhd911l6xdu/a2P+QaAAB0HSFOdwD2Onv2rNTX18uWLVtk6NCh8r//+7/y5JNPSlVVlaSnp/NyAABgAMvVwOlOwFn//u//LpmZmfL555+36viKigqJjo6W8vJyiYqK6uTe3ervf/+77ffppNDQUKe7YLuQEPP+Nq2trXW6C7a6du2a012wnfrj2zSxsbG236fT71FdBUPA0D8kvXv35pkAAMAQ5v2ZjWY+++wzefXVV+U//uM/vD4z1dXVujX96woAAHRdVAADxLp168SyrBZbQUFBs3O+/PJL+c53viMPP/ywLFmyxOttp6Wl6XK6uyUkJHT2wwEAAJ2IOYABoqysTLeWDBo0SMLDwxvD37Rp0+Tee++V7du3S1BQUJsqgCoEMgfQHswBNANzAAMfcwDtwRzA1mEIOED06dNHt9YoKSnR4S8pKUm2bdvWYvhTwsLCdAMAAIGBAGgYVfmbOnWqDBw4UG/78re//a3xa3FxcQ72DAAA2IUAaJjDhw/Lp59+qtuAAQOafY0dgQAAMAOLQAzz2GOP6aDnqQEAADMQAAEAAAxDAAQAADBMkNMdAAAAgL0IgAAAAIYhAAIAABiGAAgAAGAYAiAAAIBhCIAAAACGIQACAAAYhgAIAABgGAIgAACAYQiAAAAAhiEAAgAAGIYACAAAYBgCIAAAgGEIgAAAAIYhAAIAABiGAAgAAGAYAiAAAIBhCIAAAACGIQACAAAYhgAIAABgGAIgAACAYQiAAAAAhiEAAgAAGIYACAAAYBgCIAAAgGEIgAAAAIYhAAIAABiGAAgAAGCYEKc7ALRVr169eNIQcMLCwpzugq0iIyOd7gJgNCqAAAAAhiEAAgAAGIYACAAAYBgCIAAAgGEIgAAAAIYhAAIAABiGAAgAAGAYAiAAAIBhCIAAAACGIQACAAAYhgAIAABgGAIgAACAYQiAAAAAhiEAAgAAGIYACAAAYBgCIAAAgGEIgAb67ne/KwMHDpTw8HCJj4+XhQsXypdfful0twAAgE0IgAaaNm2a/O53v5Nz585JTk6OfPbZZ/LQQw853S0AAGATy9XApvuCn3rrrbdk/vz5Ul1dLd26dfN5fEVFhURHR0t5eblERUXZ0EMAAFqH96jWoQJouK+++kp++9vfysSJE1sV/gAAQNdHADTU6tWrJTIyUmJiYqSoqEj279/v9VhVGVR/UTVtAACg6yIABoh169aJZVkttoKCgsbjn3/+eSksLJTDhw9LcHCwPProo+JtNkBaWpoe8nW3hIQEux4WAADoBMwBDBBlZWW6tWTQoEF65e/NvvjiCx3qjh07JsnJyR4rgKq5qQqgOp45gAAAf8McwNYJad1h8Hd9+vTRrT3clb+mIa+psLAw3QAAQGAgABrmgw8+0G3SpEnSq1cv+fzzz+XFF1+Uu+66y2P1DwAABB7mABomIiJC9uzZI9OnT5fExER5/PHHZdSoUZKXl0eVDwAAQ1ABNMy3vvUteeedd5zuBgAAcBAVQAAAAMNQAUS7F42wHyAAwN+435v4oLOWEQDRZpWVlfpf9gMEAPjze5XauxaesQ8g2qy+vl6+/PJL6dGjh95g+mbufQKLi4sD+rOCeZyBhdczsPB6mvt6qsqfCn/9+/eXoCBmunlDBRBtpn6gBgwY4PM49UMayAHQjccZWHg9Awuvp5mvJ5U/34jGAAAAhiEAAgAAGCZ4XQOnO4HAExwcLFOnTpWQkMCeZcDjDCy8noGF1zOwmPJ62oVFIAAAAIZhCBgAAMAwBEAAAADDEAABAAAMQwAEAAAwDAEQne673/2uDBw4UMLDwyU+Pl4WLlyoP0kkkPzlL3+RJ554QgYPHiwRERFy1113ydq1a6WmpsbprnW49evXy8SJE+WOO+6Qnj17Ot2dDvPaa6/p1099nyYlJcm7777rdJc6VH5+vsydO1d/OoL6BJ99+/Y53aVOkZaWJuPHj9efVBQbGyvz58+Xc+fOOd2tDpeZmSmjR49u3Bg5OTlZ3n77bae7Zcvrq75/V6xY4XRXujwCIDrdtGnT5He/+53+JZyTkyOfffaZPPTQQwH1zJ89e1Z/RN6WLVvko48+kl/+8peyefNmWbNmjdNd63Aq1D788MPy1FNPOd2VDvPmm2/qN5Qf//jHUlhYKJMnT5aUlBQpKipyumsdpqqqSsaMGSMbN250uiudKi8vT5YtWybHjx+XI0eOSG1trcycOVM//kCiPo3p5ZdfloKCAt3uv/9+mTdvnv79E6hOnDghWVlZOviiA7gAm+3fv9/V8BecqyFIBPRz/2//9m+uhoqS093oNNu2bXNFR0c73Y0Occ8997iWLl3a7Lphw4a5XnjhBYd61LnUr/69e/c63Q1bXLp0ST/ehmDodFc6Xa9evVy/+tWvnO5Gp6isrHR94xvfcDWEetd9993nWr58udNd6vKoAMJWX331lfz2t7/VQ4jdunUL6Ge/vLxcevfu7XQ30IqK5smTJ3WVqCl1+dixYzx/AfBzqATyz2JdXZ3s3r1bVznVUHAgUlXdOXPmyIwZM5zuSsAgAMIWq1evlsjISImJidHDag1VwIB+5tUw96uvvioNVSWnuwIfysrK9Btov379ml2vLl+4cIHnrwtrKHLIypUrZdKkSTJq1Cinu9Phzpw5I927d5ewsDD9u6ahqisjRoxwulsdToXbU6dO6fl/6DgEQLSL+gRBNRG3pabmpbg9//zzem7V4cOH9cf5PProo/qXc6A9TkUtcPnOd76j58ktWbLEoZ53/uMMNOoxNqW+P2++Dl3L008/LR9++KHs2rXL6a50isTERDl9+rSe76jm5C5atEg+/vhjp7vVoYqLi6VhuFeys7P1Ai10HD4KDu2umqjWkkGDBnn8gf3iiy8kISFBD6/5+3BFWx+nCn9q0cu9994r27dvl6CgoIB9PdXjUwsnrly50tnd6/QhYLWi+fe//708+OCDjderNx315qoWFQQaFWxVtUitkA1UzzzzjF7prFY/q9XdJlDDo2oHArUYLVCo11D9XKrCgZuq2KvvYfX7tbq6utnX0Hp8ojLapU+fPrq1h7vyp35wA+lxlpSU6PCnthDZtm1blwl/t/t6dnWhoaH6NVMrRpsGQHVZrapE16J+v6jwpwJubm6uMeHP/di7wu/Vtpg+fboe6m5q8eLFMmzYMD21iPDXfgRAdKoPPvhANzUHp1evXvL555/Liy++qP9K9ffqX1uoyt/UqVP1fofp6enyt7/9rfFrcXFxDvas46k5nGoxj/pX/SWuqmTK0KFD9XykrkjNE1P7U44bN05/X6qtJtTjC6Q5nFevXpVPP/208fL58+f1a6cWR6jv20BaLLBz5049z1jtBeiexxkdHa336AwUaosptVWRGk2prKzU8+RU4D106JDTXetQ6jW8ef6mez55IM7rtJWTS5AR+Brm37gaqmKuhjcZV8NEZVfDMKLebqNhGNjprnX4lijqx8lTCzQN84w8Ps6jR4863bXbsmnTJtedd97paqgIuv7hH/4h4LYNUa+Pp9dNvZ6BxNvPofoZDSSPP/544/dr3759XQ2VMlfDHGunu2ULtoHpGMwBBAAAMEzXmaQEAACADkEABAAAMAwBEAAAwDAEQAAAAMMQAAEAAAxDAAQAADAMARAAAMAwBEAAAADDEAABAAAMQwAEAJv98Ic/lLlz5/K8A3AMARAAbHb69Gm5++67ed4BOIYACAA2+9Of/iRjx47leQfgGAIgANiouLhYLl++3FgBvHLlih4OnjhxopSWlvJaALAFARAAbB7+jY6OlsGDB8uZM2dk/PjxEh8fL7m5ufpfALADARAAbA6AY8aMkV27dsmUKVNk1apVkpWVJaGhobwOAGxjuRrYdm8AYLgFCxbI0aNH9f//+Mc/6qFfALAbFUAAsLkCqELg9evX9fy/m917771SUFCg/79o0SLJzMzk9QHQ4agAAoBNKisr9fy/kydP6pXAy5cvl2PHjsnIkSMbj1FVwTfeeEO+/e1vy+effy6bNm3i9QHQ4QiAAGCTd999V+6//365evWqhIWFyfPPPy85OTnywQcfSJ8+fRqPUyuEY2Nj5eDBgxISEsLrA6DDMQQMADZRVb9hw4bp8Kf8/Oc/lxEjRsj3vvc9qamp0depMKiGhnv27En4A9BpqAACgJ8oKSmRlJQU2b9/vw6FO3fulOHDhzvdLQABiAogAPiBa9euyUMPPSQbN27UewT+6Ec/kp/+9KdOdwtAgKICCAAAYBgqgAAAAIYhAAIAABiGAAgAAGAYAiAAAIBhCIAAAACGIQACAAAY5v8B5totgmk7fGYAAAAASUVORK5CYII=", "text/html": [ "\n", "