{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Quick Start"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The pairinteraction software allows for calculation properties of Rydberg atoms and Rydberg pair potentials, optionally considering electric and magnetic fields in arbitrary directions.\n",
    "The software is written so that it works for atoms with a single valence electron as well as for atoms which require multi-channel quantum defect theory.\n",
    "To achieve a high performance, the backend of the software is written in C++ and wrapped with nanobind, providing a Python interface and graphical user interface.\n",
    "The software can be installed via `pip` with the following command on Linux, MacOS 13+, and Windows:\n",
    "\n",
    "```bash\n",
    "pip install pairinteraction\n",
    "```\n",
    "\n",
    "The construction of Hamiltonians is accelerated by using pre-calculated matrix elements, which are stored in database tables.\n",
    "These tables are automatically downloaded from GitHub when needed \\[[1](https://github.com/pairinteraction/database-sqdt/releases),[2](https://github.com/pairinteraction/database-mqdt/releases)\\].\n",
    "Once the tables are downloaded, they are cached locally and the software can be used without an internet connection.\n",
    "To download tables manually, for example, the tables for Rubdium and Yb171 (for a full list of supported species, see the \"Identifier\" column in the \"quantum defect references\" table of pairinteraction's [README](https://www.pairinteraction.org/pairinteraction/sphinx/html/index.html#quantum-defects)), run:\n",
    "\n",
    "```bash\n",
    "pairinteraction download Rb Yb171_mqdt\n",
    "```\n",
    "\n",
    "To test the software, run:\n",
    "\n",
    "```bash\n",
    "pairinteraction --log-level INFO test\n",
    "```\n",
    "\n",
    "To start the graphical user interface, run:\n",
    "\n",
    "```bash\n",
    "pairinteraction --log-level INFO gui\n",
    "```\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Using the Python Interface"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Import the necessary libraries and initialize pairinteraction's database that contains pre-calculated atomic states and matrix elements. By setting `download_missing=True`, database tables are automatically downloaded from GitHub when needed."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%pip install matplotlib numpy pairinteraction\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import pairinteraction.real as pi  # Use \"import pairinteraction.complex as pi\" if complex\n",
    "from pairinteraction.visualization.colormaps import alphamagma\n",
    "\n",
    "if pi.Database.get_global_database() is None:\n",
    "    pi.Database.initialize_global_database(download_missing=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The library is structured around Python classes that can be used to model systems of Rydberg atoms. In the following, we describe a system of a single Rydberg atom and calculate the energy shift due to an electric and magnetic field."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Energy shift: -0.009 GHz\n"
     ]
    }
   ],
   "source": [
    "# Construct and diagonalize the Hamiltonian of a single Rydberg atom in an electric and\n",
    "# magnetic field\n",
    "ket = pi.KetAtom(\"Rb\", n=60, l=0, m=0.5)\n",
    "basis = pi.BasisAtom(\"Rb\", n=(ket.n - 3, ket.n + 3), l=(ket.l - 3, ket.l + 3))\n",
    "system = (\n",
    "    pi.SystemAtom(basis)\n",
    "    .set_electric_field([0, 0, 0.5], unit=\"V/cm\")\n",
    "    .set_magnetic_field([0, 0, 10], unit=\"gauss\")\n",
    ")\n",
    "pi.diagonalize(\n",
    "    [system], diagonalizer=\"lapacke_evd\", float_type=\"float64\"\n",
    ")  # Use \"float32\" for faster calculations\n",
    "\n",
    "# Obtain the eigenenergies and the index of the eigenstate having the largest overlap with\n",
    "# the ket we are interested in\n",
    "eigenenergies = system.get_eigenenergies(unit=\"GHz\")\n",
    "overlaps = system.get_eigenbasis().get_overlaps(ket)\n",
    "index = np.argmax(overlaps)\n",
    "\n",
    "# Calculate the energy shift\n",
    "energy_shift = eigenenergies[index] - ket.get_energy(unit=\"GHz\")\n",
    "print(f\"Energy shift: {energy_shift:.3f} GHz\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Two single-atom systems can be combined to model a pair of atoms. In the following, we show how to calculate the pair potential for the ket state that we have specified above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of states in the pair basis: 343\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjMAAAG2CAYAAACKxwc0AAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjEsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvc2/+5QAAAAlwSFlzAAAPYQAAD2EBqD+naQAAbcJJREFUeJzt3Xd8Tff/B/DXOTcLGUIkkYgEURFEQogdq2LTZRWh+H3tkWo1tWqGFlW1laqS0mXUJkLtHVuMhBhJ0EgiO7n3/P5Ic+tKkHtvkntv7uvpcR511ue8j7bOO58pSJIkgYiIiMhAiboOgIiIiEgbTGaIiIjIoDGZISIiIoPGZIaIiIgMGpMZIiIiMmhMZoiIiMigMZkhIiIig8ZkhoiIiAwakxkiIiIyaExmiIiIyKAxmSEiIjJSf//9N7p16wYnJycIgoBt27a99Z7Dhw+jQYMGMDc3h7u7O9avX1/scb4NkxkiIiIjlZqaivr162PZsmWFuj46OhpdunRBmzZtEBERgfHjx2Po0KHYt29fMUf6ZgIXmiQiIiJBELB161b07NnztddMmjQJu3btwtWrV5XH+vTpg8TEROzdu7cEoiyYic6erCMKhQKPHz+GlZUVBEHQdThERKTHJEnCixcv4OTkBFEsvsaMjIwMZGVlaV2OJEn5vm3m5uYwNzfXumwAOHnyJNq3b69yLCAgAOPHjy+S8jVldMnM48eP4eLiouswiIjIgDx48ABVqlQplrIzMjJQpkyZIinL0tISKSkpKsemT5+Or776qkjKj4uLg4ODg8oxBwcHJCcnIz09vcjeQ11Gl8xYWVkByP0P09raWsfREBGRPktOToaLi4vy21Ec/quRkWldVkpKSr7vW1HVyugzo0tm8qrfrK2tmcwQEVGhlEy3BAECNH+OhNwusMX5fXN0dER8fLzKsfj4eFhbW+usVgYwwmSGiIhIPwmANklTCQznadq0KXbv3q1y7MCBA2jatGnxP/wNODSbiIhIL4hFsKknJSUFERERiIiIAJA79DoiIgIxMTEAgODgYAwcOFB5/fDhwxEVFYXPP/8cN2/exPLly/Hrr79iwoQJGr1xUdGLZGbZsmVwc3ODhYUF/Pz8cObMmddeu379egiCoLJZWFiUYLRERESlw7lz5+Dj4wMfHx8AQFBQEHx8fDBt2jQAQGxsrDKxAYBq1aph165dOHDgAOrXr4+FCxfihx9+QEBAgE7iz6PzZqYtW7YgKCgIK1euhJ+fHxYvXoyAgABERkbC3t6+wHusra0RGRmp3OcQayIiMnSCln1mAPVbmlq3bo03TTdX0Oy+rVu3xsWLF9V8UvHSec3MokWLMGzYMAwePBienp5YuXIlypYti3Xr1r32HkEQ4OjoqNxeHSZGRERkaARRpvVmrHSazGRlZeH8+fMqE/CIooj27dvj5MmTr70vJSUFrq6ucHFxQY8ePXDt2rXXXpuZmYnk5GSVjYiIiEoPnSYzz549g1wuL3ACnri4uALvqVWrFtatW4ft27dj48aNUCgUaNasGR4+fFjg9SEhIbCxsVFunDCPiIj0U8l3AC4tDO7NmzZtioEDB8Lb2xv+/v74888/UalSJaxatarA64ODg5GUlKTcHjx4UMIRExERvZ0giFpvxkqnHYDt7Owgk8kKnIDH0dGxUGWYmprCx8cHd+7cKfB8Ua5JQURERPpHp2mcmZkZGjZsiLCwMOUxhUKBsLCwQk/AI5fLceXKFVSuXLm4wiQiIioBInKXNNB0Y82MzgQFBSEwMBC+vr5o3LgxFi9ejNTUVAwePBgAMHDgQDg7OyMkJAQAMHPmTDRp0gTu7u5ITEzEN998g/v372Po0KG6fA0iIiKtGHtTkTZ0nsz07t0bT58+xbRp0xAXFwdvb2/s3btX2Sk4JiZGZdn158+fY9iwYYiLi4OtrS0aNmyIEydOwNPTU1evQERERDokSG+aLacUSk5Oho2NDZKSkrjQJBERvVFJfDPynmFqUkmrmhlJUiA756lRft90XjNDREREbGbSBpMZIiIivWDcc8Vog39qREREZNBYM0NERKQH2MykOSYzREREekAUZVp2ANZuxW1DxhSQiIiIDBprZoiIiPSAABEC6xg0wmSGiIhID7DPjOb4p0ZEREQGjTUzREREekAQBC1rZoxqQn8VTGaIiIj0gPZ9ZpjMEOmN6GtxuHb6HsrZlEHTTrVhUdZM1yEREZEeYzJDeiX8j0tYMPJXyOUSIEmoXq8yvtnxfyhnbaHr0IiIipUgyCAIMl2HYZDYAZj0RnZWDr6b8CfkcgVEmQBRJiL6ehy2rT6u69CIiIpd3mgmbTZjxZoZ0hvJCWnITM/+txOcAAgAFEB8zHNdh0ZEVOxELfvMSEbcZ8Z40zjSO+XtysHatiwgAZIkQSFXQJIkuHo46Do0IiLSY0xmSG/ITGT4cl1flLEyh0KuACSgUfta6D60qa5DIyIqfto2MbGZiUg/1G9RA2vPfIo7lx7BsnwZ1GpQBaJovP+DEpHxEEUZRG06ABtvKxOTGdI/tpUs0ah9LV2HQUREBoLJDBERkR7InTRP85oZwYirZpjMEBER6QFth1cb84rbxvvmREREVCqwZoaIiEgPaLs2kzHXzDCZISIi0gPaLmdgzH1mjDeNIyIiolKBNTNERER6QPz3lzYlGCsmM0RERHpAwL/r0ml6v6T5vYaOyQwREZEeECDjPDMaMt46KSIiIioVWDNDRESkB2SCdmszCYLx1swwmSEiItIDuTMAa5PMKIowGsPCZiYiIiIyaKyZISIi0gOcAVhzTGaIiIj0AOeZ0ZzxvjkRERGVCqyZISIi0gvazTMDGG8HYCYzREREeoB9ZjRnvG9OREREpQJrZoiIiPSAqOWkeZIRzzPDZIaIiEgPsJlJc0xmiIiI9IBMkEGmRc0MjLhmxnjTOCIiIioVWDNDRESkBwRJgCBp0cwkCUUYjWFhMkNERKQH2GdGc8b75kRERFQqsGaGiIhID4iQQdRiBmBJq9mDDRuTGSIiIj3AZibNGe+bExERUanAmhkiIiI9IP77S1OSEddPMJkhIiLSA0xmNGe8b05ERESlAmtmSK9lZmRjb+hZxD1IgHtdZ7R93xuCYLwTQxFR6SVCBhlHM2mEyQzprazMHHzacyWunbsPURQgz1Hg4rE7+HTRh0xoiKjU4WgmzRnvm5PeO7z9Eq6euafcFwRg14bTiL4Rp7ugiIiKSW4qo80v4/0hj8kM6a3EpykQRAF5lTCCKECSJDx/mqLbwIiISpFly5bBzc0NFhYW8PPzw5kzZ954/eLFi1GrVi2UKVMGLi4umDBhAjIyMkoo2oKxmYn0lkcDF0ACFAoJgghICgnmZcxQrbajrkMjIipy4r8NTdrcr64tW7YgKCgIK1euhJ+fHxYvXoyAgABERkbC3t4+3/WhoaH44osvsG7dOjRr1gy3bt3CoEGDIAgCFi1apHHs2mLNDOktr6bVMWJWN4gyEZIElClnjulr+6OCvZWuQyMiKnK5q2Zrt6lr0aJFGDZsGAYPHgxPT0+sXLkSZcuWxbp16wq8/sSJE2jevDn69esHNzc3dOjQAX379n1rbc6rbt26hcWLF2Pbtm1IT09XO+5XsWaG9Fqvkf7o0KshnsUmobJrRZSzttB1SEREei05OVll39zcHObm5vmuy8rKwvnz5xEcHKw8Jooi2rdvj5MnTxZYdrNmzbBx40acOXMGjRs3RlRUFHbv3o0BAwa8MSZJknDixAls374dO3bswKNHj+Dv74/o6Gj069cP7dq1Q/fu3dGtWzc4Oqpf+64XNTPqttfl2bx5MwRBQM+ePYs3QNKp8naWcK/nzESGiEo17Tr//jfhnouLC2xsbJRbSEhIgc979uwZ5HI5HBwcVI47ODggLq7ggRb9+vXDzJkz0aJFC5iamqJGjRpo3bo1vvzyy9e+1++//w4HBwd89NFHeP78Ob755hs8ffoUO3fuxLVr13D58mW0bdsWmzZtQtWqVeHn56f2n53Oa2bUba/Lc+/ePUycOBEtW7YswWiJiIiKR1H1mXnw4AGsra2VxwuqldHU4cOHMXfuXCxfvhx+fn64c+cOxo0bh1mzZmHq1KkF3uPm5oYdO3bAz8+vwGk13N3dMWHCBEyYMAEJCQnYvXu32nHpvGZG3fY6AJDL5fj4448xY8YMVK9evQSjJSIi0m/W1tYq2+uSGTs7O8hkMsTHx6scj4+Pf21Tz9SpUzFgwAAMHToU9erVw3vvvYe5c+ciJCQECoWiwHt8fX3RpEmTQs0PVqFCBfTv3/+t171Kp8lMXntd+/btlcfe1l4HADNnzoS9vT2GDBny1mdkZmYiOTlZZSMiItI3onLaPM03dZiZmaFhw4YICwtTHlMoFAgLC0PTpk0LvCctLQ2iqJo6yGS5Mw9LkqTmGxcdnSYzmrTXHTt2DGvXrsWaNWsK9YyQkBCVtkMXFxet4yYiIipqgiDTelNXUFAQ1qxZg59++gk3btzAiBEjkJqaisGDBwMABg4cqNJBuFu3blixYgU2b96M6OhoHDhwAFOnTkW3bt2USU1BsrOz8fnnn8Pd3R2NGzfO1/oSHx//xvvfRud9ZtTx4sULDBgwAGvWrIGdnV2h7gkODkZQUJByPzk5mQkNERERgN69e+Pp06eYNm0a4uLi4O3tjb179yorGWJiYlRqYqZMmQJBEDBlyhQ8evQIlSpVQrdu3TBnzpw3PmfOnDnYsGEDJk6ciMTERAQFBeH06dNYtWqV8hptanYESYf1QllZWShbtix+//13lRFJgYGBSExMxPbt21Wuj4iIgI+Pj0r2ltdGJ4oiIiMjUaNGjTc+Mzk5GTY2NkhKSlLpIEVERPSqkvhm5D2jnV0QTETNO+vmKDIR9myRXn7fatasiW+//RZdu3YFANy5cwedOnVCixYtsG7dOjx58gROTk6Qy+Uala/TZiZ12+s8PDxw5coVREREKLfu3bujTZs2iIiIYI0LEREZLEHL/jKCHq/N9OjRI9StW1e57+7ujsOHD+PEiRMYMGCAxklMHp03MwUFBSEwMBC+vr5o3LgxFi9enK+9ztnZGSEhIbCwsFD5wwCA8uXLA0C+40RERIZEEIRCjfh50/36ytHREXfv3oWbm5vymLOzM8LDw9GmTRsMGjRIq/J1nsyo215HREREhqVt27YIDQ1Fu3btVI47OTnh0KFDaN26tVbl67TPjC6wzwwRERVWSfaZCaj0GUy16DOTrcjEvqff6OX37f79+7h58yYCAgIKPP/48WMcOHAAgYGBGpXPKg8iIiI9IBbBpq/Wrl2LihUrvva8k5OTxokMoN/vTkRERKXAw4cP0blzZ1SpUgUjRozAnj17kJWVVWTlM5khIiLSA0IR/NJX69atQ1xcHH755RdYWVlh/PjxsLOzwwcffIANGzYgISFBq/KZzBAREekBmShovekzURTRsmVLfP3114iMjMTp06fh5+eHVatWwcnJCa1atcKCBQvw6NEj9csuhniJiIiI3qh27dr4/PPPcfz4cTx48ACBgYE4evQofvnlF7XL0vnQbCIiIvpvoUlt7jdUlSpVwpAhQwq1gHRBWDNDRESkBwRB+81Q3b17F23bttX4fiYzREREpFMpKSk4cuSIxvezmYmIiEgP5M4Vo00zk/5asmTJG89r0un3ZUxmiIiI9IC2TUX63Mw0fvx4VK5cGWZmZgWe13bOGSYzREREekDbuWL0eZ4ZV1dXzJ8/H7169SrwfEREBBo2bKhx+fpcK0VERESlQMOGDXH+/PnXnhcEAdosFcmaGSIiIj0gCrmbNvfrq5kzZyItLe215z09PREdHa1x+UxmiIiI9IAgihBFzRtMBD1ubPH09HzjeVNTU7i6umpcvv6+OREREVEhMJkhIiLSA2IRbIasbdu2mDVr1hubo17H0N+djJAkSVp1FCMi0kfGPAMwAFStWhVhYWHw8PBQ+172mSGDIZcrsHzWDvz549+QyxXo+GEjfDqvF8wtTHUdGhERaWn9+vUAgOTkZLXvZTJDBmPDd/uxcelBQJIAQcCOTSdhamaCz77urevQiIi0pm1TUWlparG2tlb7ntLy7mQEwndGQJIkyExkkMlESAoJh/6K0HVYRERForQ3M6Wnp+PYsWO4fv16vnMZGRnYsGGDxmUzmSGDYWFhBgFQ9peRADYxEREZgFu3bqF27dpo1aoV6tWrB39/f8TGxirPJyUlYfDgwRqXz2SGDEbfkW0hiAIUCgVycuQQBAEDxrTXdVhEREWiNI9mmjRpEurWrYsnT54gMjISVlZWaN68OWJiYoqkfPaZIYPRpqs3Fmz8H7ZuOA55thwBHzZCp16NdR0WEVGRKM1rM504cQIHDx6EnZ0d7Ozs8Ndff2HkyJFo2bIlwsPDUa5cOa3KZzJDBqVFQD20CKin6zCIiIqcKAoQtViTQNTjZCY9PR0mJv+lHIIgYMWKFRg9ejT8/f0RGhqqVflMZoiIiKhYeXh44Ny5c6hdu7bK8aVLlwIAunfvrlX5+tzERkREZDTyFprUZtNX7733Hn755ZcCzy1duhR9+/bVajJUJjNERER6QCiCTV8FBwdj9+7drz2/fPlyKBQKjctnMkNERETFqmrVqhgzZgz279+PnJycIi+fyQwREZEeKM3NTD///DPMzMwwatQo2NnZoXfv3ti0aRMSExOLpHwmM0RERHpBgCBovulzQ5O/vz8WLlyI27dv4/jx4/D29sb3338PR0dHtG3bFosXL0ZUVJTG5TOZISIiohJTp04dBAcH49SpU7h37x769u2LsLAw1K1bF3Xr1sWuXbvULpNDs4mIiPSAMS406ejoiGHDhmHYsGFIS0vDvn37YG5urnY5TGaIiIj0gCDmbhrfr/nIZp3LycnBs2fP8N5772l0vyEmckRERFSKXLt2DdWqVdP4ftbMEBER6QFjbGYqKkxmiIiI9IAg5G7a3K+vGjRo8Mbz6enpWpXPZIaIiEgPaDtXjD7PM3P9+nX06dPntU1JsbGxuHXrlsblq5XMKBQKHDlyBEePHsX9+/eRlpaGSpUqwcfHB+3bt4eLi4vGgRAREVHpVLduXfj5+WHEiBEFno+IiMCaNWs0Lr9QTWzp6emYPXs2XFxc0LlzZ+zZsweJiYmQyWS4c+cOpk+fjmrVqqFz5844deqUxsEQEREZK7EINn3VvHlzREZGvva8lZUVWrVqpXH5haqZeeedd9C0aVOsWbMG7777LkxNTfNdc//+fYSGhqJPnz6YPHkyhg0bpnFQRERExqY095n57rvv3ni+Ro0aCA8P17j8QiVy+/fvx6+//orOnTsXmMgAgKurK4KDg3H79m20bdtW44CIiIjIuI0cORLPnj0r9PWFSmZq165d6AJNTU1Ro0aNQl9PREREgCgKWm+lxcaNG5GcnFzo69VuYqtevToGDx6MzMxMlePPnj1D9erV1S2OiIiIkLtMpLZbaSFJ6k1nrHYyc+/ePRw/fhwtW7ZEXFyc8rhcLsf9+/fVLY6IiIhIK2onM4IgYO/evahSpQoaNmyIs2fPFkdcRERERiVvnhltNmOldjIjSRIsLS3x559/YuDAgfD398fGjRuLIzYiIiKjIRTBL2Ol9gzAwktjv0JCQlCnTh0MGzYMffv2LdLAiAorMyMbZuYmKv9tEhGR8dCoZuZl/fv3x6FDh7B79+4iC4qoMGKi4tG79Vdo4jICLaqNxtaNR3UdEhGRxgQtm5hK089z/fv3h7W1daGvV7tmRqFQ5DvWtGlTXLp0CTdv3lS3OCKNZGfnYGSvxXgc8wwSJKSkZGDGhJ/gVNUOfq0KP5UAEZG+0HZEkj7nMs+ePUNqaipcXV2Vx65du4YFCxYgNTUVPXv2RL9+/ZTnVqxYoVb5RTb7sYODA/z9/YuqOKI3enjvKR7efwqFJEEURYiiAJko4GT4VV2HRkSkEZmg/aavxowZgyVLlij3nzx5gpYtW+Ls2bPIzMzEoEGD8PPPP2tcfqFrZnx8fArVJ+HChQsaB0NUWGXLWUAAIEFC3s8jkgSUKWuu07iIiCi/U6dOYf369cr9DRs2oEKFCoiIiICJiQkWLFiAZcuWYcCAARqVX+hkpmfPnsrfS5KEkJAQDB8+HBUqVNDowUTacHCyRecPm2DX76egkCsgiiKsbMqiR78Wug6NiEgjpXltpri4OLi5uSn3Dx06hPfffx8mJrlpSPfu3RESEqJx+YVOZqZPn66yv3DhQowbN46z/pLOzPh+MGp6VsGls3dQoZI1PhnXGY7OTK6JyDAJkCBAvZlvX71fX1lbWyMxMVHZZ+bMmTMYMmSI8rwgCPlWFlCH2h2AifSFiYkMg8Z01HUYRET0Fk2aNMGSJUuwZs0a/Pnnn3jx4oXKotS3bt2Ci4uLxuUzmSEiItID2s7iq88zAM+aNQvt2rXDxo0bkZOTgy+//BK2trbK85s3b9ZqEBGTGSIiIj1Qmodme3l54caNGzh+/DgcHR3h5+encr5Pnz7w9PTUuPxCJzMvD6kCgJycHKxfvx52dnYqx8eOHatxMERERFQ62dnZoUePHgWe69Kli1ZlFzqZ+fbbb1X2HR0d840JFwSByQwREZEGSnMzU2Fs2LABzZs3R40aNdS+t9CT5kVHR791i4qKUjsAAFi2bBnc3NxgYWEBPz8/nDlz5rXX/vnnn/D19UX58uVRrlw5eHt7azXRDhERkT4QBUHrzZANGjQInp6eGDNmjNr3FtkMwJrasmULgoKCMH36dFy4cAH169dHQEAAnjx5UuD1FSpUwOTJk3Hy5ElcvnwZgwcPxuDBg7Fv374SjpyIiKjo5M0zo82mCXUqFAAgMTERo0aNQuXKlWFubo533nmnSNZnVCgUuHnzJmrXVn9JGkF6deXI10hPT0dYWBi6du0KAAgODlYZEy6TyTBr1ixYWFioFYCfnx8aNWqEpUuXAsh9GRcXF4wZMwZffPFFocpo0KABunTpglmzZr312uTkZNjY2CApKUmtRayIiMj4lMQ3I+8ZIT5fwUKm3jf0ZRnyDARf/EqtWLds2YKBAwdi5cqV8PPzw+LFi/Hbb78hMjIS9vb2+a7PyspC8+bNYW9vjy+//BLOzs64f/8+ypcvj/r162scu7YK3Wfmp59+wq5du5TJzNKlS1GnTh2UKVMGAHDz5k04OTlhwoQJhX54VlYWzp8/j+DgYOUxURTRvn17nDx58q33S5KEQ4cOITIyEvPnzy/wmszMTJWkKzk5udDxERERlRRd9JlZtGgRhg0bhsGDBwMAVq5ciV27dmHdunUFViisW7cOCQkJOHHiBExNTQFAZWbft4mLi8Pp06cRFxcHAMqRTY6OjuoH/5JCJzObNm3C559/rnIsNDRUOQPwxo0bsWzZMrWSmWfPnkEul8PBwUHluIODwxtX4E5KSoKzszMyMzMhk8mwfPlyvPvuuwVeGxISghkzZhQ6JiIiIl0pil4vr/7Qbm5uDnPz/OvWaVKhsGPHDjRt2hSjRo3C9u3bUalSJfTr1w+TJk2CTCZ7bUypqan43//+h82bN0MQBOVSSAkJCZAkCX379sWqVatQtmxZTV658H1m7ty5g3r16in3LSwsIIr/3d64cWNcv35doyDUZWVlhYiICJw9exZz5sxBUFAQDh8+XOC1wcHBSEpKUm4PHjwokRiJiIh0wcXFBTY2NsrtdWsevalCIa/m5FVRUVH4/fffIZfLsXv3bkydOhULFy7E7Nmz3xjTuHHjcObMGezatQsZGRmIj49HfHw8MjIysHv3bpw5cwbjxo3T7IWhRs1MYmKiSnPN06dPVc4rFAq111Wws7ODTCZDfHy8yvH4+Pg3VjmJogh3d3cAgLe3N27cuIGQkBC0bt0637Wvy0iJiIj0SVE1Mz148EClz0xRfgMVCgXs7e2xevVqyGQyNGzYEI8ePcI333yTbw3Hl/3xxx/YtWsXmjVrpnJcJpOhQ4cOWLduHbp27Yo1a9ZoFFeha2aqVKmCq1evvvb85cuXUaVKFbUebmZmhoYNGyIsLEx5TKFQICwsDE2bNi10OZokUkRERPpExH8JjUbbv+VYW1urbK9LZjSpUKhcuTLeeecdlSal2rVrIy4uDllZWa99N4VCATMzs9eeNzMzg0KheO35tyl0MtO5c2dMmzYNGRkZ+c6lp6djxowZGs3gFxQUhDVr1uCnn37CjRs3MGLECKSmpio7Iw0cOFClPS8kJAQHDhxAVFQUbty4gYULF+Lnn39G//791X42ERGRsdKkQqF58+a4c+eOSuJx69YtVK5c+Y3JSteuXfF///d/uHjxYr5zFy9exIgRI9CtWzeN36XQzUxffvklfv31V9SqVQujR4/GO++8AwCIjIzE0qVLlQtHqat37954+vQppk2bhri4OHh7e2Pv3r3KNryYmBiVvjmpqakYOXIkHj58iDJlysDDwwMbN25E79691X42ERGRvtBmrpi8+9UVFBSEwMBA+Pr6onHjxli8eHG+CgVnZ2dlv5sRI0Zg6dKlGDduHMaMGYPbt29j7ty5b539f+nSpejXrx8aNmwIW1tb5bDvJ0+eIDExEQEBAcopWjRR6HlmgNxZgEeMGIEDBw4g7zZBEPDuu+9i+fLlypFN+ozzzBARUWGV5Dwzi32no4yJ5vPMpOdkYPy5GWrHunTpUnzzzTfKCoUlS5YoF4Js3bo13NzcsH79euX1J0+exIQJExAREQFnZ2cMGTLkraOZ8ty8eRMnT55UGZrdtGlTeHh4qPeyr1ArmcmTkJCAO3fuAADc3d2VQ6wMAZMZIiIqLGNIZkqDQjczvaxChQpo3LhxUcdCRERktHTRzFRSMjMzIYqicqK9u3fvYt26dYiJiYGrqyuGDBmCatWqaVx+oToADx8+HA8fPixUgVu2bMGmTZs0DoiIiMgYiUWw6auAgABs374dAHD8+HHUqVMHO3fuRHZ2Nnbv3o26desWaub/1ylUzUylSpVQp04dNG/eHN26dYOvry+cnJxgYWGB58+f4/r16zh27Bg2b94MJycnrF69WuOAiLQhSRKyMrNhbvH6XvVERPpIF8sZlJSLFy8q126aPHkyRo4ciUWLFinPT506FZ999hmOHTumUfmFSuRmzZqFW7duoXnz5li+fDmaNGmCqlWrwt7eHrVq1cLAgQMRFRWF1atX49SpU/Dy8tIoGCJt7Pj1KJq4D0ODqoPRveXniLr1SNchERERALlcDrlcDiC3E3BgYKDK+UGDBuHSpUsal1/oPjMODg6YPHkyJk+ejOfPnyMmJgbp6emws7NDjRo1IOhzYx2VehdORyJ4zEpICgmCAETdeoxhveZj96kFrKUhIoMg/PtLm/v1lZ+fH/766y94eHigRo0auHTpksoq2xEREVoNJtKoA7CtrS1sbW01fihRUTt++ApEQYAk5k4XoFBIiH30D6LvxMKjrquuwyMieitRkCAKag8wVrlfX82ePRudOnVCamoq+vbti08//RS3b99G7dq1ERkZiSVLlqhMkKsujZIZIn1Ttqw5FNLLK85KgCCgbDnNhzkSEVHRaNq0Kfbs2YOgoCCcPn0aADBnzhwAgJOTE7766quSWWiSSJ9179USP63cg4RnyZDLFYAgoEPXRnBxs9d1aEREhSIAWjUU6W8jU66mTZvi5MmTePr0KaKioqBQKFC5cmW4ublpXTaTGSoVKjmUx5b9s7Bm8XY8iU+EV4MaGDyqC/tyEZHBELQczWQof91VqlQJlSpVKtIymcxQqVHZuSKmffOJrsMgIqI3iImJgampKSpXrqw8Fhsbi+zsbFStWlWjMtWeY2f69Om4f/++Rg8jIiKigpXmSfNe5ubmhnbt2qkca9u2bfHPAPyy7du3o0aNGmjXrh1CQ0ORmZmp8cOJiIgolyBqvxmC8PBw/PTTTyrHNmzYgEOHDmlcptqvHhERgbNnz6JOnToYN24cHB0dMWLECJw9e1bjIIiIiMg4+Pv7o1GjRirHGjVqBH9/f43L1CiP8/HxwZIlS/D48WOsXbsWDx8+RPPmzeHl5YXvvvsOSUlJGgdE/5H+iYZ0929I0SchZSTrOhwiIipGxtLMBABPnjzB0aNHcfToUTx58kTr8rR6d0mSkJ2djaysLEiSBFtbWyxduhQuLi7YsmWL1sEZM+nxZeDSn0D0CeDOEeB8KKT0RF2HRURExSRvbSZtNn334sULDBgwAM7OzvD394e/vz+cnZ3Rv39/rSpCNEpmzp8/j9GjR6Ny5cqYMGECfHx8cOPGDRw5cgS3b9/GnDlzMHbsWI2DMnaSpADuHgUUOYC5FWBhCaQ9Bx5G6Do0IiIqJgIkrTd9N3ToUJw+fRo7d+5EYmIiEhMTsXPnTpw7dw7/+9//NC5X7aHZ9erVw82bN9GhQwesXbsW3bp1g0wmU7mmb9++Ws3kZ/QUciAnC5CZ/jtxwL/pdk6GTsMiIiLSxs6dO7Fv3z60aNFCeSwgIABr1qxBx44dNS5X7WSmV69e+OSTT+Ds7Pzaa+zs7KBQKDQOytgJMlNINk7AP1G5yYxCDggywLry228mIiKDpG1TkSE0M1WsWBE2Njb5jtvY2Gi15qPazUxTp059YyJDRcSzE2BbFZAUuYmMa2PAqZ6uoyIiomIiitpv+m7KlCkICgpCXFyc8lhcXBw+++wzTJ06VeNy1a6ZCQoKKvC4IAiwsLCAu7s7evToodVS3gQIFtaQGvYFMlMBmSkEUy6YSEREhsfHx0dlaZnbt2+jatWqytl+Y2JiYG5ujqdPn2rcb0btZObixYu4cOEC5HI5atWqBQC4desWZDIZPDw8sHz5cnz66ac4duwYPD09NQqKcgmCCFhY6ToMIiIqAdp24tXXDsA9e/Ys9meonczk1br8+OOPsLa2BgAkJSVh6NChaNGiBYYNG4Z+/fphwoQJ2LdvX5EHTEREVBqV1j4z06dPL/ZnqN3C9s0332DWrFnKRAbI7bjz1Vdf4euvv0bZsmUxbdo0nD9/vkgDJSIiotIpNjYWo0eP1vh+tZOZpKSkAmfre/r0KZKTc2epLV++PLKysjQOioiIyPhIRbDpr2vXrmHp0qVYvXo1EhMTAQDPnj3D+PHjUb16dYSHh2tcttrJTI8ePfDJJ59g69atePjwIR4+fIitW7diyJAhynaxM2fO4J133tE4KCIiImNTmpcz2LFjB3x8fDB27FgMHz4cvr6+CA8PR+3atXHz5k1s3boV165d07h8tfvMrFq1ChMmTECfPn2Qk5OTW4iJCQIDA/Htt98CADw8PPDDDz9oHBRRUYh99AxHDl6AIAho06Eh7B05wo6ISBdmz56NUaNGYdasWfjhhx8QFBSEsWPHYvfu3fkWndSEIElSoeul5HI5jh8/jnr16sHU1BRRUVEAgOrVq8PS0lLrYEpCcnIybGxskJSUpNLvh0qX61ei0b/HVLxITgMA2NhaYfPO2XCv5aLjyIjIkJTENyPvGb+3nopyJppPw5Gak4EPD8/Sy++bjY0Nzp8/D3d3d8jlcpibm2Pv3r1o3759kZSvVq2UTCZDhw4dkJiYCEtLS3h5ecHLy8tgEhkyHjMmrUFKSgYEQYAgCHiRlIrZk3/UdVhERK8lCJLWm7568eKFMsGSyWQoU6YMqlevXmTlq93MVLduXURFRaFatWpFFgRRUXt4/wkUcgVkstx8XS5XIOZe3FvuIiLSHW37vehznxkA2Ldvn3IpA4VCgbCwMFy9elXlmu7du2tUttrJzOzZszFx4kTMmjULDRs2RLly5VTO61vVFhknj7puePY0EXmtqDKZCM96TMCJiHQlMDBQZf/V2X4FQYBcLteobLWTmc6dOwPIzZ5enp5YkiStAiEqSjMX/B/695iOBzHxAIBq1Rwxde4QHUdFRPR6oiBB1KKpSJt7i1txLz6tdjKjzThwopLi7GKPXUe/RcT5WxAEAT6+78CijLmuwyIiei1ByN20ud9YqZ3M+Pv7F0ccREWubDkLNGvlpeswiIiomKmdzADA0aNHsWrVKkRFReG3336Ds7Mzfv75Z1SrVg0tWrQo6hiJiIhKvdK6NlO1atVUuqUU1vjx4zF27NhCXat2MvPHH39gwIAB+Pjjj3HhwgVkZmYCyF3mYO7cudi9e7e6RRIRERk9EVomM0UWSdFav369Rve5ubkV+lqNRjOtXLkSAwcOxObNm5XHmzdvjtmzZ6tbHBEREZViJdE9Re1ELjIyEq1atcp33MbGRrlwFBEREalHgKT1pu+ePn362nNXrlzRuFy1kxlHR0fcuXMn3/Fjx44V6Wx+RERExiSvz4w2m76rV68edu3ale/4ggUL0LhxY43LVTuZGTZsGMaNG4fTp09DEAQ8fvwYmzZtwsSJEzFixAiNAyEiIqLSLSgoCB988AFGjBiB9PR0PHr0CO3atcPXX3+N0NBQjctVu8/MF198AYVCgXbt2iEtLQ2tWrWCubk5Jk6ciDFjxmgcCBERkTHTtqnIEJqZPv/8c7z77rsYMGAAvLy8kJCQAD8/P1y+fBmOjo4al6t2MiMIAiZPnozPPvsMd+7cQUpKCjw9PbnYJBERkRZK8wzAL3N3d0fdunXxxx9/AAB69+6tVSIDaDjPDACYmZnB09NTq4cTERFRLkHLfi+GMAPw8ePH0b9/f1SoUAGXL1/G8ePHMWbMGOzevRsrV66Era2tRuWqncykpqZi3rx5CAsLw5MnT/KttxAVFaVRIERERFS6tW3bFhMmTMCsWbNgamqK2rVro02bNujfvz/q1auHhw8falSu2snM0KFDceTIEQwYMACVK1fWaFY/IiIiKohhNBVpav/+/fnmnalRowaOHz+OOXPmaFyu2snMnj17sGvXLjRv3lzjhxIREZEqmZi7aXO/vnvdBHqiKGLq1Kkal6v2q9va2qJChQoaP5BIF1JT0jF76hp82GUixv3va9y/F6vrkIiIjEp6ejqOHTuG69ev5zuXkZGBDRs2aFy22snMrFmzMG3aNKSlpWn8UKKSJJfL8Um/6Vi7ahvOnrqGv7YeQc8OExAf94+uQyMi+o9QBJueunXrFmrXro1WrVqhXr168Pf3R2zsfz9UJiUlYfDgwRqXr3Yys3DhQuzbtw8ODg6oV68eGjRooLIR6Zub1+/h9IkrkBQSZDIRkiQh8Xkydm0/quvQiIiUREhab/pq0qRJqFu3Lp48eYLIyEhYWVmhefPmiImJKZLy1e4z07NnzyJ5MFFJycrMBqA6bFEQBGT+e5yIiIrXiRMncPDgQdjZ2cHOzg5//fUXRo4ciZYtWyI8PBzlypXTqny1k5np06dr9UCikla7TjU4VbFH3ONnuVMJCAJkMhFt2zfSdWhEREqledK89PR0mJj8l3IIgoAVK1Zg9OjR8Pf312opA0CDZiYASExMxA8//IDg4GAkJCQAAC5cuIBHjx5pFQxRcbAoY47QP0Pg5fMOTM1M4excCT/8PB21PN10HRoRkZIgCFpv+srDwwPnzp3Ld3zp0qXo0aMHunfvrlX5atfMXL58Ge3bt4eNjQ3u3buHYcOGoUKFCvjzzz8RExOjVW9kouLiVt0J2/Z9q+swiIiM0nvvvYdffvkFAwYMyHdu6dKlUCgUWLlypcblq10zExQUhEGDBuH27duwsLBQHu/cuTP+/vtvjQMhIiIyZqW5A3BwcDB279792vPLly/Pt6KAOtSumTl79ixWrVqV77izszPi4uI0DoSIiMiYCYIEQYt+L9rca+jUTmbMzc2RnJyc7/itW7dQqVKlIgmKiIjI2IhaLjSpzb2GTu1mpu7du2PmzJnIzs4b7iogJiYGkyZNwgcffFDkARIRERG9iUaT5qWkpMDe3h7p6enw9/eHu7s7rKystFokioiIyJgJovabsVL71W1sbHDgwAH89ddfWLJkCUaPHo3du3fjyJEjGk96s2zZMri5ucHCwgJ+fn44c+bMa69ds2YNWrZsCVtbW9ja2qJ9+/ZvvJ6IiMgQ5M0zo81mrNTuM5OnRYsWaNGihdYBbNmyBUFBQVi5ciX8/PywePFiBAQEIDIyEvb29vmuP3z4MPr27YtmzZrBwsIC8+fPR4cOHXDt2jU4OztrHQ8REREVn9u3byM8PBxPnjzJN4Jp2rRpGpUpSJKkdioXFhaGsLCwAgNZt26dWmX5+fmhUaNGWLp0KQBAoVDAxcUFY8aMwRdffPHW++VyOWxtbbF06VIMHDjwrdcnJyfDxsYGSUlJsLa2VitWIiIyLiXxzch7xvn3P4WlqbnG5aRkZ6Lhnwv1+vu2Zs0ajBgxAnZ2dnB0dFSZ6E8QBFy4cEGjctWumZkxYwZmzpwJX19fVK5cWasZB7OysnD+/HkEBwcrj4miiPbt2+PkyZOFKiMtLQ3Z2dmoUKFCgeczMzORmZmp3C9oJBYREZGuCdBwWv6X7td3s2fPxpw5czBp0qQiLVftZGblypVYv359gbP4qevZs2eQy+VwcHBQOe7g4ICbN28WqoxJkybByckJ7du3L/B8SEgIZsyYoXWsREREpJ3nz5/jo48+KvJy1U4Cs7Ky0KxZsyIPRBPz5s3D5s2bsXXrVpXZiF8WHByMpKQk5fbgwYMSjpL0UVzsM8yfvQ6TJnyLrb+FQYPWViKiIpU3aZ42m7776KOPsH///iIvV+2amaFDhyI0NBRTp07V+uF2dnaQyWSIj49XOR4fHw9HR8c33rtgwQLMmzcPBw8ehJeX12uvMzc3h7m55m2QVPrEx/2DTm1GIOGfJEAQsOmnXbh5IxrB04bqOjQiMmK6mjRv2bJl+OabbxAXF4f69evj+++/R+PGjd963+bNm9G3b1/06NED27Zte+11S5YsUf7e3d0dU6dOxalTp1CvXj2YmpqqXDt27FiN3kHtZCYjIwOrV69WJhGvBrJo0aJCl2VmZoaGDRsiLCwMPXv2BJDbATgsLAyjR49+7X1ff/015syZg3379sHX11fdVyAjt2HdDiT8kwSFQqHs87ViyRaMHNcHNjaWOo6OiKjkqDuiOM+9e/cwceJEtGzZ8q3P+PZb1UV+LS0tceTIERw5ckTluCAIJZfMXL58Gd7e3gCAq1ev5gtEXUFBQQgMDISvry8aN26MxYsXIzU1FYMHDwYADBw4EM7OzggJCQEAzJ8/H9OmTUNoaCjc3NyU60FZWlrC0pIfInq7xOcvIAiC8r9XQQAkSUJyUgqTGSLSGV2szbRo0SIMGzZM+c1duXIldu3ahXXr1r12RLFcLsfHH3+MGTNm4OjRo0hMTHzjM6Kjo9WOS11qJzPh4eFFGkDv3r3x9OlTTJs2DXFxcfD29sbevXuVnYJjYmIgiv917VmxYgWysrLw4YcfqpQzffp0fPXVV0UaG5VOfs288NPaHQCk3KRGFOHkXAmVnbi2GBHpjrYT3+Xd++qo3dd1t9B0RPHMmTNhb2+PIUOG4OjRo2rFOHPmTEycOBFly5ZVOZ6eno5vvvlG43lmNB4FdufOHezbtw/p6ekAoFUHytGjR+P+/fvIzMzE6dOn4efnpzx3+PBhrF+/Xrl/7949SJKUb2MiQ4XVrac/Jnw+AKIoAwA4O9vjp1/mwMREpuPIiMiYCYKWyxn82zji4uICGxsb5ZbXsvGqN40ozmv1eNWxY8ewdu1arFmzRqN3nDFjBlJSUvIdT0tL02rksdo1M//88w969eqF8PBwCIKA27dvo3r16hgyZAhsbW2xcOFCjYMhKgmCIODTLwLxv1EfISkpBY6VK0ImYyJDRKXDgwcPVCbNK6pBMC9evMCAAQOwZs0a2NnZaVSGJEkFdkm5dOnSa+eLKwy1k5kJEybA1NQUMTExqF27tvJ47969ERQUxGSGDIalVVlYWpV9+4VERCWgqPrMWFtbF2oGYHVHFN+9exf37t1Dt27dlMfyVgEwMTFBZGQkatSoUeCzbG1tlX0V33nnHZWERi6XIyUlBcOHD3/7S76G2snM/v37sW/fPlSpUkXleM2aNXH//n2NAyEiIjJmJT00W90RxR4eHrhy5YrKsSlTpuDFixf47rvv4OLi8tpnLV68GJIk4ZNPPsGMGTNgY2OjEoebmxuaNm2q3gu8RO1kJjU1NV/HHQBISEjgfC5EREQGRJ0RxRYWFqhbt67K/eXLlweAfMdfFRgYCACoVq0amjVrlm9aF22pncy0bNkSGzZswKxZswDk9j9QKBT4+uuv0aZNmyINjoiIyGgIUu6mzf1qUndEsSaSk5OVzV4+Pj5IT09XDh56laYLZKq9avbVq1fRrl07NGjQAIcOHUL37t1x7do1JCQk4Pjx469tL9MXXDWbiIgKqyRXzY7sPwZWZpq3cLzIykStjd/r3fdNJpMhNjYW9vb2EEWxwA7AeR2D5XK5Rs9Qu2ambt26uHXrFpYuXQorKyukpKTg/fffx6hRo1C5cmWNgiAiIqLS6dChQ8qRSocOHdJogt23UTuZAQAbGxtMnjy5qGMhIiIyWoLw31wxmt6vj/z9/ZW/b926dbE8Q6PlDAoiCAIsLCxQtWpVdgQmIiJSk64WmixJrVq1QuvWreHv74/mzZvDwsKiSMpVO5nx9vZWVhHldbd5ucrI1NQUvXv3xqpVq4osSKLidujgKez66wjMzEzRr39X1Kv/jq5DIiIqdTp06IC///4bixYtQk5ODnx9fVWSm4JGSxeG2l2Ut27dipo1a2L16tW4dOkSLl26hNWrV6NWrVoIDQ3F2rVrcejQIUyZMkWjgIhK2pbQPfi412fYvGk3Nvy4HZ3f/R/Onrn69huJiIqQIEpab/puypQp2L9/PxITExEeHo6uXbvi3Llz6NKlS8nOADxnzhx89913CAgIUB6rV68eqlSpgqlTp+LMmTMoV64cPv30UyxYsEDjwIhKyrw5qyFJgEKR+xeBJMnx3cIN2Ljlax1HRkTGpLT2mSlIVFQUrly5gkuXLuHy5cuwsrJCq1atNC5P7WTmypUrcHV1zXfc1dVVOTOgt7c3YmNjNQ6KqCQlJ6UC+K+5VKFQ4Pnz5DfdQkRU5LStXTGEmpl+/frhyJEjyMzMRKtWreDv748vvvgCXl5eWo1yUruZycPDA/PmzUNWVpbyWHZ2NubNmwcPDw8AwKNHj/Ktwkmkr1r6N4RMJkKhUEChUEAQBLRp21jXYRERlTqbN29GdnY2hg4diuHDh2PYsGGoX7++1sO11a6ZWbZsGbp3744qVarAy8sLQG5tjVwux86dOwHkVh+NHDlSq8CISsqi77/AkIFTcOpEBAABH/UJwNigAboOi4iMjZbNTDCAZqZ//vkHR48exeHDhxEcHIwbN27A29sbrVu3RuvWrdGhQweNylV7BmAgdxnwTZs24datWwCAWrVqoV+/frCystIoiJLEGYCpIJIkITHxBUxNTLiSNhEpleQMwPeGjYK1FjMAJ2dlwm3NMoP6vt25cwezZ8/Gpk2boFAoSm4GYACwsrLSaqluIn0jCAJsbQ3jf34iIkP1zz//4MiRIzh8+DAOHz6M69evo3z58ujWrZvK5HrqKlQys2PHDnTq1AmmpqbYsWPHG6/t3r27xsEQEREZK2PoAGxvbw87Ozu0bNkSw4YNQ+vWrVGvXj2tyy1UMtOzZ0/ExcXB3t4ePXv2fO112iwSRUREZNy0XDUb+p/MXL58GXXq1CnycguVzCgUigJ/T0RERFRYxZHIABr2mSEiIqKiJYgCBC0WWNLmXkNX6HlmOnfujKSkJOX+vHnzkJiYqNz/559/4OnpWaTBERERGQtB1H4zVoV+9X379iEzM1O5P3fuXCQkJCj3c3JyEBkZWbTREemAQqHAb1v2YPaM5djw41ZkZWXrOiQiInqDQjczvTodjQbT0xDpPUmSMHr4DPz+6z6YmJpAniPHtq1h+PXPxTAxYassERUfQZAgaNEBWJt7DZ0RV0oR5Rdx4QZ+/20fAEAhz+3sfuzvc9i355guwyIiI5DXZ0abzRCEhYWha9euqFGjBmrUqIGuXbvi4MGDWpVZ6GRGEIR8aydou5YCkb558uSf3NGRr/y3/ST+Hx1FRERGQ5S03/Tc8uXL0bFjR1hZWWHcuHEYN24crK2t0blzZyxbtkzjctVqZho0aBDMzXOnWs7IyMDw4cNRrlw5AFDpT0NkqOrUrQlTUxNkZ+cAyP3vXhAAnwbs3E5EpK25c+fi22+/xejRo5XHxo4di+bNm2Pu3LkYNWqURuUWumYmMDAQ9vb2sLGxgY2NDfr37w8nJyflvr29PQYOHKhREET6ooqLI5avngELCzNIkgSZTMSc+Z/Cu0FtXYdGRKWcIGi/6bvExER07Ngx3/EOHTqojJhWV6FrZn788UeNH0JFL7cDdjZyZ3w0hWDMY/KKWI/32sG/dSPExMTCydkedna2ug6JiIyAtsOrDeEz0L17d2zduhWfffaZyvHt27eja9euGpfL4RkGSJIUkKRngPTi3yNlANEegsB/nUWlvK01ynPhSSKiIuXp6Yk5c+bg8OHDaNq0KQDg1KlTOH78OD799FMsWbJEee3YsWMLXS6/fgbpBSA9hwRzAMK/+yYQBHtdB0ZERJoSod0YYwOomVm7di1sbW1x/fp1XL9+XXm8fPnyWLt2rXJfEAQmM6WdJGVCggBBMP133xSQMnQcFRERacMYmpmio6OLpVwDeHV6lQATANK/zU0SgBwwLyUiImPFL6AhEqwhSOkAUv89YA5BrKDLiEq9nJwcZGZmo1y5MroOhYhKqYLmc1P3fkPw8OFD7NixAzExMcjKylI5t2jRIo3KZDJjgATBBBAdAaQjdzSThbLJiYqWJElYMH8NFi1ci5xsORr61sW6DV/DyYn9k4ioiAnQrr3EAHKZsLAwdO/eHdWrV8fNmzdRt25d3Lt3D5IkoUGDBhqXy2YmAyUIMgiCJQTBiolMMdocuhPzQ1YhOysHkiTh4oVrGNR/oq7DIiIySMHBwZg4cSKuXLkCCwsL/PHHH3jw4AH8/f3x0UcfaVwukxmiNzh86BREUYAoihBFEXK5AhcvXMeLF6lvv5mISA3GsDbTjRs3lBPsmpiYID09HZaWlpg5cybmz5+vcblMZojewNK6HARBUFklXiYTYWFhpsOoiKhUEotg03PlypVT9pOpXLky7t69qzz37Nkzjcs1gFcn0p1h/+sDMzMziKKA3P5JwLgJg2FqyqY9IipiQhFseq5JkyY4duwYAKBz58749NNPMWfOHHzyySdo0qSJxuWyAzDRG3h4VMeB8A1YvfIXvEhOhX8bP3w8oIeuwyIiMkiLFi1CSkoKAGDGjBlISUnBli1bULNmTY1HMgFMZojeyqN2DSz6boquwyCi0k7276bN/XquevXqyt+XK1cOK1euLJJymcwQERHpAWOZZ6Y4MJkhIiKiYmNra1voRCshIUGjZzCZISIi0gfaduLV04qZxYsXF/szmMwQqSEp6QXGjJyOA/uOokwZC0yYOASjxwYadfUuERURmZC7aXO/HgoMDASQuyxMaGgoAgIC4ODgUKTP4NBsIjWM+L/J2L0zHJmZWXj+PAnTp3yLLb/8peuwiIj0nomJCYYPH46MjIwiL5vJDFEh5eTkYP/ev6FQKCAIubMCA8BfO8J0HBkRlQZ5HYC12fRd48aNcfHixSIvl81MRIUkiiJMTEyQmZmlnBVYFAVYWJjrOjQiKg20ncXXAKonRo4ciU8//RQPHz5Ew4YNUa5cOZXzXl5eGpXLZIaokERRxIjRA/DdonVQKBT/1swIGPp/fXQdGhGRQejTJ/fvy7FjxyqP5f1wKAgC5HK5RuUymSFSw9TpY+BgXxG7dx+GpWVZjBo9EE2bab5sPRGRkiAA2iwWaQDNTNHR0cVSLpMZIjWIoojho/pj+Kj+ug6FiEoZQczdtLlf37m6uhZLuQbw6kREREZAFLXfDMDPP/+M5s2bw8nJCffv3weQOxfN9u3bNS7TMN6cSI9lZGQiJuYRMjOzdB0KEZFeW7FiBYKCgtC5c2ckJiYq+8iUL19eq8n1mMwQaWH7tv2oVrUZ6nl2QLWqzbF71yFdh0REhkosgk3Pff/991izZg0mT54Mmey/lTF9fX1x5coVjcs1gFcn0k+3b9/DJ4M/Q1pa7gRQqalpGNg/CNHRD3QcGREZIkEUtN70XXR0NHx8fPIdNzc3R2pqqsblMpkh0tDZM5eQnZUN4L/VarOzs3Dh/FVdhkVEpLeqVauGiIiIfMf37t2L2rVra1wuRzMRaahCxfIFzLgpoEIFG53EQ0QGThC0G15tAEOzg4KCMGrUKGRkZECSJJw5cwa//PILQkJC8MMPP2hcLpMZIg21a9cMTZo2wJnTFyGKIuRyBVq08EXLVo11HRoRGSJRy3lmDKCZaejQoShTpgymTJmCtLQ09OvXD05OTvjuu++UE+ppQufNTMuWLYObmxssLCzg5+eHM2fOvPbaa9eu4YMPPoCbmxsEQSiRZcWJXsfU1BTbdqzBVzMm4OP+72H2nIn4Y9tqmJjwZwQiotf5+OOPcfv2baSkpCAuLg4PHz7EkCFDtCpTp3/rbtmyBUFBQVi5ciX8/PywePFiBAQEIDIyEvb29vmuT0tLQ/Xq1fHRRx9hwoQJOoiYSFWZMhYYN0G7/wmJiAAYRc3M7Nmz8fHHH6NatWooW7YsypYtWyTl6rRmZtGiRRg2bBgGDx4MT09PrFy5EmXLlsW6desKvL5Ro0b45ptv0KdPH5ibc3E/0j9ZWVk4cvgk9u45jGfPEnQdDhEZEiOYNO+3336Du7s7mjVrhuXLl+PZs2dFUq7O3jwrKwvnz59H+/bt/wtGFNG+fXucPHmyyJ6TmZmJ5ORklY2oOCQmJqNdmz7o1mUwen04HA3qd8TFCxzZRESU59KlS7h8+TJat26NBQsWwMnJCV26dEFoaCjS0tI0LldnycyzZ88gl8vh4OCgctzBwQFxcXFF9pyQkBDY2NgoNxcXlyIrm+hlIXO+x5XLN5X7yckpGDbkMx1GRESGRBAErTdDUKdOHcydOxdRUVEIDw+Hm5sbxo8fD0dHR43L1P86KS0FBwcjKSlJuT14wAnNqHhcu3oLcrlc+ZeKQqHAnTv3oFAodB0aERkCAdrN/msYuYyKcuXKoUyZMjAzM0N2drbG5egsmbGzs4NMJkN8fLzK8fj4eK2ys1eZm5vD2tpaZSMqDtVruMLERAZJkiBJEgRBgLNzZYgG0I5NRHpAFLTsM2MY2Ux0dDTmzJmDOnXqwNfXFxcvXsSMGTO0apXR2d+yZmZmaNiwIcLCwpTHFAoFwsLC0LRpU12FRaSxyVPGoIqLEwAJAGBuYYaVq0N0GxQRkR5p0qQJ3N3d8fvvv2Pw4MG4f/8+wsLCMGTIENjYaD7hqE6HZgcFBSEwMBC+vr5o3LgxFi9ejNTUVAwePBgAMHDgQDg7OyMkJPeDkJWVhevXryt//+jRI0RERMDS0hLu7u46ew8iAHBwrITjJ7dh755wpKVlwL91E1Srxj5aRFRIRjA0u127dli3bh08PT2LtFydJjO9e/fG06dPMW3aNMTFxcHb2xt79+5VdgqOiYlRqaJ//PixygJVCxYswIIFC+Dv74/Dhw+XdPhE+VhbW6JX727K/SuXb+D48XMoX94a3bq/i3LlimZOBSIqhYwgmZkzZw4AKIdk29nZFUm5Op+qdPTo0Rg9enSB515NUNzc3CBJUglERaS933/fhSGDP/23D40CHh41cfDQFtjYWOk6NCKiEpeYmIjJkydjy5YteP78OQDA1tYWffr0wezZs1G+fHmNy2bPRKJikJ2djVEjvoRCofg3ARdw61YUlny3VtehEZG+yltoUptNA+osK7RmzRq0bNkStra2sLW1Rfv27d94fZ6EhAT4+fnhp59+wgcffICFCxdi4cKFeP/997F+/Xo0bdpUmeBogskMUTFISEhCWlo6gJfnjgBiYh7pODIi0leCKGi9qStvWaHp06fjwoULqF+/PgICAvDkyZMCrz98+DD69u2L8PBwnDx5Ei4uLujQoQMePXrz320zZ86EmZkZ7t69i1WrVmH8+PEYP348Vq9ejTt37sDU1BQzZ85UO/48TGaIikGlShVgZ1cBgiAoh2orFBLq1HlH16ERESmpu6zQpk2bMHLkSHh7e8PDwwM//PCDciTym2zbtg0LFizIN1EuADg6OuLrr7/G1q1bNX4PJjNExUAURWzY+J1Kh9927Vtg5KhAHUZFRHqtiNZmenUJn8zMzAIfVxTLCqWlpSE7OxsVKlR443WxsbGoU6fOa8/XrVtXq3lmdN4BmKi0atnSD5evHsSliOuwtrFEw4ZeOHsmAjExj+BZpxbq1q2l6xCJSJ8I0G4W33/vfXXZnunTp+Orr77Kd/mblhW6efNmvusLMmnSJDg5OakkRAWxs7PDvXv3UKVKlQLPR0dHvzUhehMmM0TFqFKlimj/bktIkoSxo6fgxx+3AJAgCCJC5n2J0WMG6zpEIiplHjx4oDLbvbm5ebE8Z968edi8eTMOHz4MCwuLN14bEBCAyZMn48CBAzAzM1M5l5mZialTp6Jjx44ax8JkhqgE7Nt3+N9EBgBy120K/mIuOndpi+rVXXUaGxHpiZeaijS+Hyj00j3aLCu0YMECzJs3DwcPHoSXl9dbnzVz5kz4+vqiZs2aGDVqFDw8PCBJEm7cuIHly5cjMzMTP//881vLeR32mSEqAbdvRUMUVVe3lSQJd+7c03VoRKQv8ibN02ZTg6bLCn399deYNWsW9u7dC19f30I9q0qVKjh58iQ8PT0RHByMnj174r333sPkyZPh6emJ48eP52seUwdrZohKgLu7GxSK/yZ8zFuIsnr1qjqMioj0iqBlzYyg/r3qLis0f/58TJs2DaGhoXBzc1N22rW0tISlpeUbn1WtWjXs2bMHz58/x+3btwEA7u7uWvWVycNkhqgEdOzUBh9//B42bdqqTGRmzvoM7u7VdB0aERkxdZcVWrFiBbKysvDhhx+qlPO6TsYFsbW1RePGjYvsHQBAkIxsfYDk5GTY2NggKSmpUG2KREVFkiQcPXoaMTGPULt2TezZfQgrlq+HXC5Hr9498M2CqcXWUY+INFMS34y8ZyT+OQfW5d7ckfaN5aRmoPz7k43y+8aaGaISIggCWrVqAgD4bvEahMxdolxrbN26XyCTifh2seYzYBKRgdNiSQLl/UaKHYCJdOD333ZCkiSIoghRFKGQK/DH77t0HRYRkUFizQyRDphbmCtHNAn//jT16twLRGRkimhotjEy3jcn0qHRo/+bLE+hUEAQBIwYGYirV28iNjb+DXcSUamlo1WzSwMmM0Q60PO9TgjdvAJt2jRDy5Z+mPDp/2HRtyvRuFFH1Kjuhy+/nAsj65tPRKQxNjMR6UiPHgHo0SMAmZmZcHdviqTEJOW5xd+uho9PPXz0UTcdRkhEJUqE2hPf5bvfSBnxqxPph3v3HuCfZwlQKCTl7MAmJiY4ffqCrkMjopJURKtmGyPjfXMiPWFnV1HZGRjInY9GoVDAwd5Ox5ERERkGJjNEOlaxoi2+/HIcBEGATCaDKIqwt6+IH3/8BfaVaqNTxz548OCxrsMkouJWwmszlSbsM0OkByZPGQ+v+p44deo8srOzsXzZj1AoFFAoFDh69BS6dvkY587vh6mpqa5DJaLiwqHZGjPeNyfSM926dcCcOcEoX94agKQcsp2TI0dk5B1cv35L1yESUXESRO03I2W8b06kp0xNTfHqqGxBAJ49S0BMzCMO2SYiegWTGSI906tXD5QtWwYymQySJEEmE2FlZYUunfvinZp+CAjohaSkZF2HSURFTSiCzUgxmSHSM66uVXD4yFZ07NQW3t51Ub26G168eAFJkiBJEo4dPY2Jn07XdZhEVNQELYdlG3EzEzsAE+mhOnU88Mcf6wAA7jUaKfvPAIBcLsfff5/UZXhERHrFeNM4IgNRqZIdRFFU1syIooi0tHQ0a9YJvXsNwY0b7BhMVCpw0jyNGe+bExmI2XOCIYoiTExMYGJiAkmS8OzZM5w/F4G//toH/1bd8ODBI12HSUTaYjKjMeN9cyID0a5dKxw99hfGjB2Kvn3fhyTJlc1OcrkcKSkp+PXX7boOk4hIZ5jMEBmABg28MG/eVIwYGVjg+Z/W/wL/Vt0wdcpcZGRklHB0RFQkOAOwxtgBmMiA1K1bG9WruyEm5gFycuT/1s4ocPv2XSgUEs6cuYDbd6KwefMPug6ViNQmQrs6BuOtnzDeNycyQObm5tiz91c0a9YY5cqVhY2NNURRVK64rVAosG3rbty6dRdpaWm6DpeIqEQwmSEyMG5uLjhw8E8kPL+LAQN7QXyl058kSahXtzlsy1fHsGHjkZ2draNIiUgtbGbSGJMZIgPWvXtHyOW5zU2SJCnnoskbxv3zhi2YM2eRjqMkokIRBS1HMzGZISID1KpVM2zYsBzOzpVRpowFRFGAJOWOdMprdlr/Yyj69h2G4OBZSEh4ruuQieh1ODRbY8b75kSlRK/ePXE36jwSk6Lh6GivrKXJW5Dy8eM4bP1zJ75bvAKtWnZFSkqqjiMmIipaTGaISpGvvpoEAJDJZJDJ/vvfW5Ik5OTIcevWHYwa9RnmzFmI48dP6ypMIiqIIGi/GSkOzSYqRQYG9oFdpYr484+diIuPx/594fmu+eWXPyCTiZgxYz5WrfoWgwd/rINIiSgfbZuKjLiZickMUSnTufO76Nz5XTx79g88avkhJSUVCoVCeV6SFMjJyd0fPfpz3LhxC1ZWlhg4sA9cXV10FTYRkcaMN40rRSRJgpSTCEXWQ0jZsZAkDsUlwM6uInbsCIXDv/1oTE1NIIhQdg4GgKysLHz//SrMmbMQDRu2xq1bd3UcNZERYzOTxpjMlAY58VBk3ICUeQ+KzLtQpN9gQkMAgGbNG+PevQg8eXoLM2d9qewYnNc5GABycuT/rvGUis6dP4K3d0v06zcMDx8+1mHkREZI0HIkk2C8n3Q2Mxk4SZJDkfUQACDIrCBJCkCRDOQ8BUyddBwd6QNBEGBjY41x4/6Hc+cu4I8//lJJZvJqaeRyOe7ffwAAiIy8jdOnzyE09AcIAlCvnifKlCmjk/iJiN7GeNO4UkMBQA4IuXmp8G9mLklyHcZE+sjExAShoT/g8uXjOHZsL3x8vGBiIstXUyMIAnJy5IiJeYgWLQLQrFkH1K7th5s3b+kweiIjIIjab0bKeN+81DABREtAkQFJyoKkSAMggyBa6jow0kOCIMDDoyb8/Bpi8+a1qFbNTVkzU5C8JCc2Ng5t23ZHo0Zt8O67PREefrSkQiYyHpw0T2NsZjJwgiBANK8BBQAoUgDBBKKZMyArr+PISN9Vr+6Gy5eP4c6daGRkZCAg4AMkJSVDLs+t1Xt5eQS5XI4nT57iyZOnEEURR4+ewMSJYyGTyVC3bm18+GGPNyZFRETFiclMKSCIFhAtPAEpGxBkEASZrkMiA2FiYgIPj5oAgAMHtmL48CDcunUHpqamSEhIUCY2efKWSJAkYN68b2Fqaors7Gzs2XMAPXt2QU5ODlq1ag47u4q6eB0iw6btYpFGvDaTIL3cWG4EkpOTYWNjg6SkJFhbW+s6HCK9dP/+AzRp0h4JCQnKPjQAlEsl5Hl1HwAqVqyAX375ARkZmbCzq4hGjRqw1oYMVkl8M/KekXj1V1hbldW8nBdpKF+3l1F+31gzQ0T5uLq64OLFvxEa+jvS09Oxe/cBnDt3MV9NTUGeP09EQMD7kCQ5BAF4//3uGD16OJ48eYp69ergnXfcS+ANiAyQtnPFGPEPDayZIaK3SkxMwqhRn2LfvkOwsLBAXFx8gdf9V1MjIXeknSqZTIa5c6fDzMwMkiShR48ucHWtWrzBE2mhRGtmrv2mfc1MnY+M8vvGZIaI1LZw4VJ88cVXkCQJoigCkKBQ5P1VkvtPQcj9Z0F/xeTNQly2bFkEBvZDfPwTVK7sgM8+Gw8np8ol9BZEb1eSyczzG3/A2qqc5uW8SIVt7Q+M8vvGZIaINPL4cSwePnyMMmUs0LPnx8oJ93Llb44qqH+NIACSBIiiCFEUYWVlibJlyyIxMRHNmzfBd999jfv3H8LU1ARNmjSCubl5Mb8VkaoSTWYi/9Q+man1vlF+39hnhog04uRUWVmLcvHi39i/PxzZ2dl4/DgWX3wxTTnyKU/+ROa/5EaSJOTk5OD580Q8f54IAAgLOwIvrybIzs4BAHh6eiAgoB3u3o2Gm1tVBAdPhCgKSEpKRtWqLjAx4V9nRMaK//cTkdasra3x4Yc9lPv+/i1w6tQZVKhgi19++R179uxXnhNF8d/h3a+vFBYEAXK5HC/3N75+/SauX78JQRAgk4lYv34TkpOTAQBublUxfvxoXLx4CWXKWOB//xuCqlWrICbmIVxcnFG+fPkif2eiIifIcjdt7jdSbGYiomKlUCiwd+9BxMbGwcXFGZ9++iVu3IgEkNshWJIklRoc4M1DwPN+/19H4/8SJJlMhCCIyqHgOTnZMDExxYQJo3Ht2g2kpKSge/fO6N69C06ePANra2t06NAW5ubmyMnJYe0O5VOizUx3dmjfzOTe3Si/b0xmiKhEZWZm4vz5CEiSBHNzM3z88VBERUXDwsIcOTlyKBSKNw4BV010Xt/J+E33yWQy5TPc3WvgxYsXePLkKapXr4bPPhuPY8dOICsrC++/3wN163ri1KmzsLOriI4d20Mul+OffxLg6OgAmcx4fxI2FkxmDAOTGSLSueTkZFhaWuLy5auYOXMe4uKeoGbNGti6dQfS0zOU1+Xvh5N/+HdBtTgFnStIXg2P+O8aN7m/F5QjtWrUqI6YmAfIycmBvX0l9OnzIQ4eDIckSQgM/BiWlpY4f/4iHB3tMWrU/3DmzHk8efIUjRs3hKenByIjb8PExAQ1a9aAKIqQy+VMiPRciSYzd3dqn8zU6GqU3zcmM0Skt6Kj72HfvjCYmJjAxaUKRo0Kwv37MbCxsYZcLkdaWioEQUR2drYyEXnZ25KZNyU+b/O6601MTCBJEkxMTJCZmQkgN0mqXNkRjx/HAgAaNPBGcvIL3LlzF46ODhg+fCj27TuAf/5JQEBAezRo4I29ew/AwsIC//vfJ7h58xZu3bqD6tWr4cMPe2D//kNISUmFv39zVKhgi4iIy7C2tkb9+vWQkZGB2Ng4ODs7wcLCAgqFAgqFgk1oGirRZCZqt/bJTPXORvl904tkZtmyZfjmm28QFxeH+vXr4/vvv0fjxo1fe/1vv/2GqVOn4t69e6hZsybmz5+Pzp07F+pZTGaIDFt6ejosLCzw+HEsVq9eh8TEJPj4eOOHH9bj1KkzEATAxaUK7t9/AFNTU+Tk5KjdDPXqfm4fHAl5p4sjKXr1n6IovrQWlqR8l3LlyiIlJRUAYG5uBjMzU7x48QKCIKBevTqIjLyFrKxsWFqWQ7t2bbB79z7k5OQgIKA9TExMcPToCVSsWAHDhw/FgQNhuH//ARo1aoB27drg99+3QaFQoG/fXkhI+AenT5+Do6MDhgwJxNatO/DkyVP4+TWCt7cXduzYBVEU8dFH7yMm5gEuXbqMKlWqoFOnDjh06AgSEp6jSZNGqFixAk6dOgtzczO0aNEMT548xe3bd1CtmhtcXasiKioaL16kwMPjHZiamuL+/RhYWpaDvb09srOz8fTpM9jbV4KJiQmysrKQk5ODsmVzJ5bLycmBTCYr1uUymMwYBp0nM1u2bMHAgQOxcuVK+Pn5YfHixfjtt98QGRkJe3v7fNefOHECrVq1QkhICLp27YrQ0FDMnz8fFy5cQN26dd/6PCYzRKWTJEl4/vw5LCwsYGZmhnXrfsbt23dQo0Y1JCe/wIIF3yEtLR2+vj44c+a8stbEzMzs35FT8kInIuokM9pc+6Znqsrf3FZwGbnz+rycMOU1q+XVar3atyiPJEmQyWTIzs5WOS6TyZCVlaUsw9bWFs+fPwcAmJqaoEyZMkhOfgEAqFLFGY8fx0KhUPybfHniypXrACQ4ODjA3NwcMTG58xW1bt0Kp0+fRXp6OqytreHv3wK7du2FQqFAixbNkJmZiXPnzqNcuXLo378vdu3ai8ePY1GvXh00adIYv/22FXK5HAMG9MPXX8+Gqalpof6MXlWiyUz0Xlhba5HMJKfCtlpHo/y+6TyZ8fPzQ6NGjbB06VIAuW3ULi4uGDNmDL744ot81/fu3RupqanYuXOn8liTJk3g7e2NlStXvvV5TGaI6P79GJw4cRrW1laoX78e5s9fhKioe6hTpzZevHiBdet+hkIhR4MG3rh06Sqys7MBAGXKlEFmZiZEUUROTk6+crVtqlI/8VG8tI9CJUnaJFDqvl9RKcxzX03Q8oiiiKCgsZg/f7ZGzy7RZOb+fu2TGdcORvl902kjalZWFs6fP4/g4GDlMVEU0b59e5w8ebLAe06ePImgoCCVYwEBAdi2bVtxhkpEpYira1WVNaG+/36hyvnvv18IuVwOc3NzREffw9mzF1CxYgV4edXB0qWr8OhRLLy9vZCamooffvgJkiTh3Xfb4q+/diMuLh4mJjLUqvUOrl69DgAwMZFBoci/9EPucPL/RlblNjNB2Zz18gf81Y957n5u89fLTWDaUidhKa4aqcJe++rEi692EM/b37Rps8bJTIkSxNxNm/uNlE6TmWfPnkEul8PBwUHluIODA27evFngPXFxcQVeHxcXV+D1mZmZyupkAMpJtoiIXsfExETZYbZaNTdUq+amPDdjxhSVaydN+u+Hq6VLFyI2Ng4VKtjCwsICx46dxNOnz9CgQX0kJSVj69YdkMlk6NGjKzZv/h0REZfh6loVvr4N8N13y5CQ8BytWjVDfPxT/P33cYiiiHbt/HHo0N/K/iFOTpXx4MFDAICFhQWysjIgCLm1EaampsjOzoaJiYmyr5Aoivk+/q82L71cm/H2BKrgfXWufZU21xaGmZlmTUxkOEp99/aQkBDMmDFD12EQkRGQyWSoUsVZud+qVXOV897eXsrf169fT+Xc0KGBKvtpaWkwMTGBmZkZ/vnnH0RH34era1VYW1vh8OGjSElJRbNmfnj8OBbh4UdgaWmJ7t27IDR0C6KiouHp6QFPz9pYu/YnZGZmokuXjrhxIxIHDx5CpUp2eO+9Hvjhhx/x8OEj+PjUh7OzMzZv/g2SpEDnzh1x8WIE7t2LgZmZGdq1a429e/dDkoAKFWxhZWWJe/diIAiAq6srHjx4ALk8tx+Mo6MDYmNzf7g0MTGBKOaONpMkCeXKlUVqahpkMhlycnJgbm6G7OwcyGQicnLkyo7OkqSAJEGZeMnl8n/77IgABJV5iPJqtvL+WVBiNnbsKO3/5ZYEUQBELWpXxOLrCK3vdJrM2NnZQSaTIT4+XuV4fHw8HB0dC7zH0dFRreuDg4NVmqWSk5Ph4uKiZeRERMUrb8QOAFSsWBEVK1ZU7gcEtFf+vnJlRzRs6KPcnzhxvEo5bdu2fu0zhg8fqrK/Zs0yAP/V4Dx//hzW1tYwMTFBQkICnj37B66uVSEIAq5cuQZRFOHlVRcxMQ9w/fpNVKniDE9PD4SFHcbz58/h59cIAHDgwCFYWFigW7dO+Pvv47hx4yaqV6+Gli2bY+PGX5CSkoJ3322HlJRU7N9/EFZWlujT5yP8+ed2REXdg4fHO2jWrAnWr9+IjIwMdOvWGY8ePUZ4+N+oWLEC+vbthV9++RX3799HgwY+8PSsjY0bQ5GTI0f//n3xf//3idp//roh/rtpc79x0osOwI0bN8b3338PILcNuWrVqhg9evRrOwCnpaXhr7/+Uh5r1qwZvLy82AGYiIiKVIl2AH4YDmtrSy3KSYFtlTZG+X3TeTNTUFAQAgMD4evri8aNG2Px4sVITU3F4MGDAQADBw6Es7MzQkJCAADjxo2Dv78/Fi5ciC5dumDz5s04d+4cVq9ercvXICIi0o4oatnMZLw1MzpPZnr37o2nT59i2rRpiIuLg7e3N/bu3avs5BsTE6OcWhzIrYUJDQ3FlClT8OWXX6JmzZrYtm1boeaYISIi0ltcNVtjOm9mKmlsZiIiosIq0Wam2GPaNzNVbmGU3zed18wQERER2MykBSYzRERE+oDNTBoz3jSOiIiISgXWzBAREekDUZa7aXO/kWIyQ0REpBc4aZ6mmMwQERHpA1EGiFp8lo24ZsZ40zgiIiLCsmXL4ObmBgsLC/j5+eHMmTNvvP63336Dh4cHLCwsUK9ePezevbuEIn09JjNERET6QDDTflPTli1bEBQUhOnTp+PChQuoX78+AgIC8OTJkwKvP3HiBPr27YshQ4bg4sWL6NmzJ3r27ImrV69q+/Za4aR5REREr1GSk+Zp+wxNyvHz80OjRo2wdOlSALnrI7q4uGDMmDGvXR8xNTUVO3fuVB5r0qQJvL29C7U+YnExuj4zeblbcnKyjiMhIiJ9l/etKImf+7X9LuXd/2o55ubmMDc3z3d9VlYWzp8/j+DgYOUxURTRvn17nDx5ssBnnDx5EkFBQSrHAgICsG3bNq1i15bRJTMvXrwAALi4uOg4EiIiMhQvXryAjY1NsZRtZmYGR0fHIvkuWVpa5itn+vTp+Oqrr/Jd++zZM8jlcuVaiHkcHBxw8+bNAsuPi4sr8Pq4uDjtAteS0SUzTk5OePDgAaysrCAIgsblJCcnw8XFBQ8ePDCK5ipje1/A+N6Z71u68X01I0kSXrx4AScnpyKMTpWFhQWio6ORlZWldVmSJOX7thVUK1PaGF0yI4oiqlSpUmTlWVtbG8VfDHmM7X0B43tnvm/pxvdVX3HVyLzMwsICFhYWxf6cl9nZ2UEmkyE+Pl7leHx8PBwdHQu8x9HRUa3rSwpHMxERERkhMzMzNGzYEGFhYcpjCoUCYWFhaNq0aYH3NG3aVOV6ADhw4MBrry8pRlczQ0RERLmCgoIQGBgIX19fNG7cGIsXL0ZqaioGDx4MABg4cCCcnZ0REhICABg3bhz8/f2xcOFCdOnSBZs3b8a5c+ewevVqXb4GkxlNmZubY/r06UbRFgkY3/sCxvfOfN/Sje9LBenduzeePn2KadOmIS4uDt7e3ti7d6+yk29MTAxE8b9GnGbNmiE0NBRTpkzBl19+iZo1a2Lbtm2oW7eurl4BgBHOM0NERESlC/vMEBERkUFjMkNEREQGjckMERERGTQmM0RERGTQmMyo6e+//0a3bt3g5OQEQRB0vh5FcQsJCUGjRo1gZWUFe3t79OzZE5GRkboOq9isWLECXl5eyom2mjZtij179ug6rBIzb948CIKA8ePH6zqUYvHVV19BEASVzcPDQ9dhFatHjx6hf//+qFixIsqUKYN69erh3Llzug6r2Li5ueX7dywIAkaNGqXr0KgYMZlRU2pqKurXr49ly5bpOpQSceTIEYwaNQqnTp3CgQMHkJ2djQ4dOiA1NVXXoRWLKlWqYN68eTh//jzOnTuHtm3bokePHrh27ZquQyt2Z8+exapVq+Dl5aXrUIpVnTp1EBsbq9yOHTum65CKzfPnz9G8eXOYmppiz549uH79OhYuXAhbW1tdh1Zszp49q/Lv98CBAwCAjz76SMeRUXHiPDNq6tSpEzp16qTrMErM3r17VfbXr18Pe3t7nD9/Hq1atdJRVMWnW7duKvtz5szBihUrcOrUKdSpU0dHURW/lJQUfPzxx1izZg1mz56t63CKlYmJic6nXi8p8+fPh4uLC3788UflsWrVqukwouJXqVIllf158+ahRo0a8Pf311FEVBJYM0NqSUpKAgBUqFBBx5EUP7lcjs2bNyM1NVXnU3UXt1GjRqFLly5o3769rkMpdrdv34aTkxOqV6+Ojz/+GDExMboOqdjs2LEDvr6++Oijj2Bvbw8fHx+sWbNG12GVmKysLGzcuBGffPKJVgsLk/5jzQwVmkKhwPjx49G8eXOdz/ZYnK5cuYKmTZsiIyMDlpaW2Lp1Kzw9PXUdVrHZvHkzLly4gLNnz+o6lGLn5+eH9evXo1atWoiNjcWMGTPQsmVLXL16FVZWVroOr8hFRUVhxYoVCAoKwpdffomzZ89i7NixMDMzQ2BgoK7DK3bbtm1DYmIiBg0apOtQqJgxmaFCGzVqFK5evVqq+xgAQK1atRAREYGkpCT8/vvvCAwMxJEjR0plQvPgwQOMGzcOBw4cKPEVe3Xh5SZiLy8v+Pn5wdXVFb/++iuGDBmiw8iKh0KhgK+vL+bOnQsA8PHxwdWrV7Fy5UqjSGbWrl2LTp06wcnJSdehUDFjMxMVyujRo7Fz506Eh4ejSpUqug6nWJmZmcHd3R0NGzZESEgI6tevj++++07XYRWL8+fP48mTJ2jQoAFMTExgYmKCI0eOYMmSJTAxMYFcLtd1iMWqfPnyeOedd3Dnzh1dh1IsKleunC8Jr127dqluWstz//59HDx4EEOHDtV1KFQCWDNDbyRJEsaMGYOtW7fi8OHDpb7zYEEUCgUyMzN1HUaxaNeuHa5cuaJybPDgwfDw8MCkSZMgk8l0FFnJSElJwd27dzFgwABdh1Ismjdvnm8qhVu3bsHV1VVHEZWcH3/8Efb29ujSpYuuQ6ESwGRGTSkpKSo/xUVHRyMiIgIVKlRA1apVdRhZ8Rg1ahRCQ0Oxfft2WFlZIS4uDgBgY2ODMmXK6Di6ohccHIxOnTqhatWqePHiBUJDQ3H48GHs27dP16EVCysrq3z9n8qVK4eKFSuWyn5REydORLdu3eDq6orHjx9j+vTpkMlk6Nu3r65DKxYTJkxAs2bNMHfuXPTq1QtnzpzB6tWrsXr1al2HVqwUCgV+/PFHBAYGwsSEnzmjIJFawsPDJQD5tsDAQF2HViwKelcA0o8//qjr0IrFJ598Irm6ukpmZmZSpUqVpHbt2kn79+/XdVglyt/fXxo3bpyuwygWvXv3lipXriyZmZlJzs7OUu/evaU7d+7oOqxi9ddff0l169aVzM3NJQ8PD2n16tW6DqnY7du3TwIgRUZG6joUKiGCJEmSbtIoIiIiIu2xAzAREREZNCYzREREZNCYzBAREZFBYzJDREREBo3JDBERERk0JjNERERk0JjMEBERkUFjMkOkI4IgYNu2bboOo1hkZWXB3d0dJ06cKPZnNWnSBH/88UexP4eI9BeTGaIiNGjQIAiCAEEQYGpqCgcHB7z77rtYt24dFAqFyrWxsbEqqzi/iaElPitXrkS1atXQrFmzYn/WlClT8MUXX+T78yUi48FkhqiIdezYEbGxsbh37x727NmDNm3aYNy4cejatStycnKU1zk6OsLc3FyHkRYPSZKwdOlSDBkypESe16lTJ7x48QJ79uwpkecRkf5hMkNUxMzNzeHo6AhnZ2c0aNAAX375JbZv3449e/Zg/fr1yuterm3JysrC6NGjUblyZVhYWMDV1RUhISEAADc3NwDAe++9B0EQlPt3795Fjx494ODgAEtLSzRq1AgHDx5UicXNzQ1z587FJ598AisrK1StWjXfIoMPHz5E3759UaFCBZQrVw6+vr44ffq08vz27dvRoEEDWFhYoHr16pgxY4ZKUvaq8+fP4+7duyqrFR8+fBiCICAxMVF5LCIiAoIg4N69ewCA9evXo3z58ti5cydq1aqFsmXL4sMPP0RaWhp++uknuLm5wdbWFmPHjoVcLleWI5PJ0LlzZ2zevPmN/16IqPRiMkNUAtq2bYv69evjzz//LPD8kiVLsGPHDvz666+IjIzEpk2blEnL2bNnAQA//vgjYmNjlfspKSno3LkzwsLCcPHiRXTs2BHdunVDTEyMStkLFy6Er68vLl68iJEjR2LEiBGIjIxUluHv749Hjx5hx44duHTpEj7//HNlk83Ro0cxcOBAjBs3DtevX8eqVauwfv16zJkz57XvevToUbzzzjuwsrJS+88pLS0NS5YswebNm7F3714cPnwY7733Hnbv3o3du3fj559/xqpVq/D777+r3Ne4cWMcPXpU7ecRUenAtdGJSoiHhwcuX75c4LmYmBjUrFkTLVq0gCAIcHV1VZ6rVKkSAKB8+fJwdHRUHq9fvz7q16+v3J81axa2bt2KHTt2YPTo0crjnTt3xsiRIwEAkyZNwrfffovw8HDUqlULoaGhePr0Kc6ePYsKFSoAANzd3ZX3zpgxA1988QUCAwMBANWrV8esWbPw+eefY/r06QW+y/379+Hk5KTWn02e7OxsrFixAjVq1AAAfPjhh/j5558RHx8PS0tLeHp6ok2bNggPD0fv3r2V9zk5OeHBgwdQKBQQRf6MRmRsmMwQlRBJkiAIQoHnBg0ahHfffRe1atVCx44d0bVrV3To0OGN5aWkpOCrr77Crl27EBsbi5ycHKSnp+ermfHy8lL+XhAEODo64smTJwBym3p8fHyUicyrLl26hOPHj6vUxMjlcmRkZCAtLQ1ly5bNd096ejosLCzeGPvrlC1bVpnIAICDgwPc3NxgaWmpciwv/jxlypSBQqFAZmYmypQpo9GzichwMZkhKiE3btxAtWrVCjzXoEEDREdHY8+ePTh48CB69eqF9u3b52tOednEiRNx4MABLFiwAO7u7ihTpgw+/PBDZGVlqVxnamqqsi8IgrIZ6W0f/pSUFMyYMQPvv/9+vnOvS1js7Oxw5cqVN5YLQKXfy5tifVP8eRISElCuXDkmMkRGiskMUQk4dOgQrly5ggkTJrz2Gmtra/Tu3Ru9e/fGhx9+iI4dOyIhIQEVKlSAqalpvo//8ePHMWjQILz33nsAchOPvM60heXl5YUffvhB+ZxXNWjQAJGRkSpNT2/j4+ODFStWFFgTFR8fj/LlywMAoqKi1Ir1Ta5evQofH58iK4+IDAsbl4mKWGZmJuLi4vDo0SNcuHABc+fORY8ePdC1a1cMHDiwwHsWLVqEX375BTdv3sStW7fw22+/wdHRUfnhd3NzQ1hYGOLi4vD8+XMAQM2aNfHnn38iIiICly5dQr9+/dSea6Vv375wdHREz549cfz4cURFReGPP/7AyZMnAQDTpk3Dhg0bMGPGDFy7dg03btzA5s2bMWXKlNeW2aZNG6SkpODatWv5zn3xxRe4efMmTp48iRkzZgAAzp07l682SV1Hjx59a7McEZVeTGaIitjevXtRuXJluLm5oWPHjggPD8eSJUuwfft2yGSyAu+xsrLC119/DV9fXzRq1Aj37t3D7t27lZ1ZFy5ciAMHDsDFxUVZA7Fo0SLY2tqiWbNm6NatGwICAtCgQQO1YjUzM8P+/fthb2+Pzp07o169epg3b54yzoCAAOzcuRP79+9Ho0aN0KRJE3z77bcqHZRfVbFiRbz33nvYtGlTvnNeXl5o2rQpevbsiYkTJ8LHxwfjxo1DRkaGWnG/7NGjRzhx4gQGDx6scRlEZNgESZIkXQdBRKXL5cuX8e677+Lu3buwtLTE4cOH0aZNGzx//lxZ21RUJk2ahOfPn+ebP4eIjAdrZoioyHl5eWH+/PmIjo4u9mfZ29tj1qxZxf4cItJfrJkhomJXnDUzRERMZoiIiMigsZmJiIiIDBqTGSIiIjJoTGaIiIjIoDGZISIiIoPGZIaIiIgMGpMZIiIiMmhMZoiIiMigMZkhIiIig8ZkhoiIiAza/wNZLWUVRhYbIgAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Construct and diagonalize the Hamiltonian for a pair of Rydberg atoms at different distances\n",
    "distances = np.linspace(1, 7, 100)  # in micrometers\n",
    "pair_energy = 2 * (ket.get_energy(unit=\"GHz\") + energy_shift)\n",
    "pair_basis = pi.BasisPair(\n",
    "    [system, system],\n",
    "    energy=(pair_energy - 10, pair_energy + 10),\n",
    "    energy_unit=\"GHz\",\n",
    "    m=(2 * ket.m, 2 * ket.m),\n",
    ")\n",
    "print(f\"Number of states in the pair basis: {pair_basis.number_of_states}\")\n",
    "\n",
    "pair_systems = [\n",
    "    pi.SystemPair(pair_basis).set_distance_vector([0, 0, d], unit=\"micrometer\")\n",
    "    for d in distances\n",
    "]\n",
    "pi.diagonalize(\n",
    "    pair_systems,\n",
    "    diagonalizer=\"lapacke_evr\",\n",
    "    float_type=\"float64\",\n",
    "    energy_range=(pair_energy, pair_energy + 0.5),\n",
    "    energy_unit=\"GHz\",\n",
    "    sort_by_energy=False,\n",
    ")\n",
    "\n",
    "# Obtain the eigenenergies and the overlaps with the ket we are interested in\n",
    "pair_eigenenergies = [s.get_eigenenergies(unit=\"GHz\") for s in pair_systems]\n",
    "pair_overlaps = [s.get_eigenbasis().get_overlaps([ket, ket]) for s in pair_systems]\n",
    "\n",
    "# Plot the pair potentials\n",
    "distances_repeated = np.hstack(\n",
    "    [d * np.ones_like(e) for d, e in zip(distances, pair_eigenenergies)]\n",
    ")\n",
    "pair_eigenenergies_flattened = np.hstack(pair_eigenenergies)\n",
    "pair_overlaps_flattened = np.hstack(pair_overlaps)\n",
    "sorter = np.argsort(pair_eigenenergies_flattened)\n",
    "\n",
    "scat = plt.scatter(\n",
    "    distances_repeated[sorter],\n",
    "    pair_eigenenergies_flattened[sorter] - pair_energy,\n",
    "    c=pair_overlaps_flattened[sorter],\n",
    "    s=10,\n",
    "    cmap=alphamagma,\n",
    "    vmin=0,\n",
    "    vmax=1,\n",
    ")\n",
    "plt.colorbar(scat, label=f\"Overlap with 2x{ket}\")\n",
    "plt.xlabel(\"Distance (μm)\")\n",
    "plt.ylabel(\"Eigenenergy (GHz)\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In addition to these basic functionalities, the pairinteraction library comes with methods for the perturbative calculation of effective Hamiltonians and the calculation of properties of Rydberg atoms, such as lifetimes. For examples, see the corresponding tutorials."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": ".venv",
   "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.13.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}