{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Flux Qubit (3JJ + L)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  Introduction "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Calculating the flux-tuning eigenvalues for a three-junction flux qubit incorporating nonzero linear loop inductance. The goal is to compare SQcircuit's solution to a case-specific solution following the two papers mentioned below.\n",
    "\n",
    "Reference for circuit model and parameters: https://journals.aps.org/prb/abstract/10.1103/PhysRevB.73.174526"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"pics/fluxQubit_pic1.png\" width=\"230\" align=\"left\" />"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note:\n",
    "\n",
    "2D Brillouin zone for periodic coordinates used an incorrect unit cell of 2x the correct size, leading to parasitic doublets in the low spectrum (with both curves of a doublet close to the correct single curve), and large deviations in the high spectrum. This error was mentioned and corrected by Yamamoto et. al. NJP 2014: https://iopscience.iop.org/article/10.1088/1367-2630/16/1/015017/pdf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Setup\n",
    "import SQcircuit as sq\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First we calculate the value of the circuit components from paper:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Fundamental constants\n",
    "h    = 6.626e-34\n",
    "GHz = 1e9\n",
    "e0   = 1.602e-19\n",
    "phi0 = h/(2*e0)\n",
    "phi0_red = phi0/2/np.pi\n",
    "\n",
    "# Circuit parameters in paper's convention\n",
    "EJ = 50.0 * GHz * h\n",
    "EC = 1.0 * GHz * h\n",
    "αQ = 0.63\n",
    "βQ = 0.15\n",
    "κQ = 0.00\n",
    "σQ = 0.00"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Pre-convert to explicit circuit elements for checking values\n",
    "L_val   = βQ*(phi0_red**2/EJ)*(1/(1+κQ) + 1/(1-κQ) + 1/αQ)\n",
    "C1_val  = (e0**2/(2*EC))*(1+κQ)\n",
    "C2_val  = (e0**2/(2*EC))*(1-κQ)\n",
    "C3_val  = (e0**2/(2*EC))*(αQ)\n",
    "JJ1_val = EJ/h*(1+κQ)\n",
    "JJ2_val = EJ/h*(1-κQ)\n",
    "JJ3_val = EJ/h*(αQ)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Circuit Description"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Initialize loop(s)\n",
    "loop1 = sq.Loop(0.0) # \"Value\" corresponds to phiExt / phi0 threading the loop (can change later)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create circuit\n",
    "L  = sq.Inductor(L_val, 'H', loops=[loop1]) \n",
    "C1 = sq.Capacitor(C1_val, 'F')\n",
    "C2 = sq.Capacitor(C2_val, 'F')\n",
    "C3 = sq.Capacitor(C3_val, 'F')\n",
    "JJ1 = sq.Junction(JJ1_val, 'Hz', loops=[loop1]) \n",
    "JJ2 = sq.Junction(JJ2_val, 'Hz', loops=[loop1])\n",
    "JJ3 = sq.Junction(JJ3_val, 'Hz', loops=[loop1])\n",
    "\n",
    "elements = {\n",
    "    (0, 1): [JJ1, C1], # Only include C_j if not included in JJ_j object\n",
    "    (1, 2): [JJ2, C2],\n",
    "    (2, 3): [JJ3, C3],\n",
    "    (0, 3): [L]\n",
    "}\n",
    "\n",
    "cr1 = sq.Circuit(elements)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\hat{H} =~\\omega_1\\hat a^\\dagger_1\\hat a_1~+~E_{C_{22}}(\\hat{n}_2-n_{g_{2}})^2~+~E_{C_{23}}(\\hat{n}_2-n_{g_{2}})(\\hat{n}_3-n_{g_{3}})~+~E_{C_{33}}(\\hat{n}_3-n_{g_{3}})^2~~-~E_{J_{1}}\\cos(-0.63\\hat{\\varphi}_1+\\hat{\\varphi}_2+0.33\\varphi_{\\text{ext}_{1}})~-~E_{J_{2}}\\cos(0.63\\hat{\\varphi}_1+\\hat{\\varphi}_3-0.33\\varphi_{\\text{ext}_{1}})~-~E_{J_{3}}\\cos(\\hat{\\varphi}_1+\\hat{\\varphi}_2-\\hat{\\varphi}_3-0.33\\varphi_{\\text{ext}_{1}})$"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$------------------------------------------------------------$"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$\\text{mode}~1:~~~~~~~~~~~\\text{harmonic}~~~~~~~~~~~\\hat{\\varphi}_1~=~\\varphi_{zp_{1}}(\\hat a_1+\\hat a^\\dagger_1)~~~~~~~~~~~\\omega_1/2\\pi~=~51.63971~~~~~~~~~~~\\varphi_{zp_{1}}~=~2.33e-01$"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$\\text{mode}~2:~~~~~~~~~~~\\text{charge}~~~~~~~~~~~~~~~~n_{g_{2}}~=~0$"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$\\text{mode}~3:~~~~~~~~~~~\\text{charge}~~~~~~~~~~~~~~~~n_{g_{3}}~=~0$"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$------------------------------------------------------------$"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$\\text{parameters}:~~~~~~~~~~~E_{C_{22}}~=~2.886~~~~~~~~~~~E_{C_{23}}~=~2.231~~~~~~~~~~~E_{C_{33}}~=~2.886~~~~~~~~~~~E_{J_{1}}~=~50.0~~~~~~~~~~~E_{J_{2}}~=~50.0~~~~~~~~~~~E_{J_{3}}~=~31.5~~~~~~~~~~~$"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$\\text{loops}:~~~~~~~~~~~~~~~~~~~~\\varphi_{\\text{ext}_{1}}/2\\pi~=~0.0~~~~~~~~~~~$"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Ask about the ordering and Hilbert spaces for the coordinates (modes)\n",
    "cr1.description()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Circuit Spectrum"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To determine the size of the Hilbert space, we specify the truncation number for each circuit mode via `set_trunc_nums()` method. Note that this is a necessary step before diagonalizing the circuit."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Set truncation numbers, desired eigenvalues, and flux sweep\n",
    "n1 = 1 # harmonic [see cr.description()]\n",
    "n2 = 6 # charge [see cr.description()]\n",
    "n3 = 6 # charge [see cr.description()]\n",
    "\n",
    "cr1.set_trunc_nums([n1, n2, n3])\n",
    "\n",
    "n_eig = 7\n",
    "n_ext = 300\n",
    "phi_ext = np.linspace(0.0, 1.0, n_ext)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To generate the spectrum of the circuit, firstly, we need to change and sweep the external flux of `loop1` by the `set_flux()` method. Then, we need to find the eigenfrequencies of the circuit that correspond to that external flux via `diag()` method. The following lines of code find the `spec` a 2D NumPy array so that each column of it contains the eigenfrequencies with respect to its external flux. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Calculate eigenvalue spectrum\n",
    "spec = np.zeros((n_eig, len(phi_ext)))\n",
    "\n",
    "for i, phi in enumerate(phi_ext):\n",
    "    loop1.set_flux(phi)\n",
    "    spec[:, i], _ = cr1.diag(n_eig=n_eig)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEKCAYAAAAfGVI8AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAABftElEQVR4nO3ddZhc1f348fcZt92Zdd9sZLNxTwiBENwhUKQ4hQJtoe5O2++vQgsUqQBFixcPECAQAiEkQDaejev67qzvuJ3fH3djZN3lvJ5nnpEr87nZzP3cI/ccIaVEURRFGbl0Ax2AoiiKMrBUIlAURRnhVCJQFEUZ4VQiUBRFGeFUIlAURRnhVCJQFEUZ4QwDHUB3JCcny7y8vIEOQ1EUZchYt25djZQypbVlQzIR5OXlUVhYONBhKIqiDBlCiINtLVNVQ4qiKCOcSgSKoigjnEoEiqIoI5xKBIqiKCOcSgSKoigjnEoEiqIoI9yQ7D7aXf5ttSBA6AS0PMRRz8Ko0x5mPcKk114LMdBhK4oyTMmYRIaixIJRZDgG0RgyIpHRGES1ZxmVENWmCxAGHZbxCb0ex4hKBLXPbINYFzYQaMnApEdnMaCzGdDZjEee7dqz3m5E7zSjjzehizNpyUVRlBFFRmJEG4NEPWFinjAxb5ioN0zMEzryOhBFBiLaib/l0RV6p5mMX8zr9dhHVCJIvWMmSImMSYhJLdPGJEi0zBuJIUMxZCh6JEu3vI8FIsR8EaKNQcIVXmK+sJbBv0yALs50ODEYnGb0iRYMSRYMSVYMCRaEUdXIKcpQEwtFidT4idT6iTYEDz8ijdpzzBNudTth0qNzGNHZjeisBnQuM8KsXVwKsx6dWX+kFkIvEHodGMSR1/qWmguhfdYXRlQiMGU5enV/Mhwl6osQ84S1K4GmENGmINFG7Tni9hPc3YAMHZX1BejjzUcSQ4oVY7od82inShCKMgjEglGC+xqIuP3aib/lEW0KHbOeMOnRu0zoXRZMGQ70ThN6l1m7ELQb0Tm02gJh1A/QkXTeiEoEvU0Y9RicenCaoY0kI6Uk5g0TqQ0QqQsQrfUffu3fVkvMq11FCJMey8REbNOSsRQkIgwqKShKf4mFogS21eLbUkNgZz1EtNK+zm7AkGTFPM6FIcWqXbwlt5TsLfph04aoEkEfE0Kgd5jQO0yYR8UftzzqDRMubca/tRZ/UQ3+TW50NgO2manYZqdhyuzdUoyiKBopJaHiZnzrqvBtciODUXRxJhzz0rFMTsKUYUdnMw50mP1CJYIBprcb0RckYilIxHXJWIJ7GvCuq8LzWQWeT8sxZjmwz07DNisVnUX9uRSlp6KeEL511XgLK4m4/QijDuu0FOyzUzHlOUdkZw91ZhlEhF6HpSUpRL1h/JvceAsraViyl8Z3D2Cfk4Z9QSbGZOtAh6ooQ06owovn0zJ8G6shIjGNiifhsmys05LRmUf2qXBkH/0gprcbcSzIxLEgk1BpM57V5Xg+r8CzphxLQSKOkzIxj3MNmzpKRekLMiYJ7KjDs6qM4L5GhFGHfU46jgWZGFNtAx3eoCGklAMdQ5fNmTNHjsT5CKLNITyfVeD9vIKYJ4wxy0H86TlYJiaNyOKsorRFRiW+zW6aPywm4vajd5pxLMjEPjdtxNT7f5kQYp2Uck5ry1SJYAjRx5lwnjWK+FNz8G2spumjEmqf3o4x3Ubc6blYpySrhKCMaDIaw7de+21EawMY0+0kXj1B+230UR/84UAlgiFIGHXY56Zjm5V2+Kqn7rkdGFKtxJ85CuvUZFVlpIwoMirxra+i6cNiovVBjFkOXNdPVKXlTlKJYAgTeoF9Ziq26Sn4t9TQ1JIQjNkOnOeNxjLWNdAhKkqfklIS2F5H47sHiFT7MObE4Vo8DktBgroY6gKVCIYBoRPYpqdgnZqsFYvfP0DNf7ZgmZCI89w8jOn2gQ5RUXpdsLiJxqX7CR1owpBsJem6iVgmJ6kE0A0qEQwjQiewz0nDNj0Zz+pymlaUUHX/euxz04k/Jw+9fWQ2kinDS6QhSOPSffg316BzGHFdMg773DRtXB6lW1QiGIaEUU/cohzsc9Np+rAEz+oy/FtriD87D/u8dFVnqgxJMhKjeVUZzcuLkRLizsgl7pRsdObBP5bPYKcSwTCmsxlxXTgG+9w0Gt7YS8Pre/CurcS1eCzm3OOHu1CUwSqwu56GJXuJuP1YJiXhunAMhkTLQIc1bKhEMAIY0+wk3zoV/2Y3DW/vx/2vTdjnpeM8f7QatkIZ1KLNIRre3It/cw2GJAvJN03GUpA40GENO+osMEIIIbBNT8UyIZGm94vxfFpGYEcdrq/kY52gfljK4CKlxL/JTcOSvcSCUeLPGkXcomw1Km8fUYlghNGZDbguHINtegp1L++i9skibDNTcV44RjUmK4NCtDFI/Wt7COyow5QbR8Ll49VwEH1MJYIRypQTR9p3ZtK0ooTmFSUEdteTcMk4rFOSBzo0ZYSSUuJbV0XDm/sgJnFeOAbHgkzVuaEf9GsiEEIcAJqBKBCRUs4RQiQCLwJ5wAHgSillfX/GNVIJgw7nWaOwTk6i/uVd1D6zXWs7uHAMOpPqiaH0n5gvTP2ru/FvrcU02kni5fkYktQou/1lICrcTpNSzjhq8KOfA8ullPnA8pb3Sj8yZTpIvX0GcYuy8a6tpPrBDYTKPAMdljJCBPY2UHXfevzb63CeN5qUW6eqJNDPBkPLy2LgqZbXTwGXDFwoI5cw6HCeN5rkW6Yig1Gq/7WR5pWlyNjQG51WGRpkNEbju/upeXQLwqQn9VvTtQZhVRXU7/p1GGohxH6gHpDAw1LKR4QQDVJK11Hr1EspE1rZ9jbgNoDc3NzZBw8e7KeoR56oVyumB4pqMee7SLxqgmpIVnpVpCFA7bM7CJc0Y5+bjvMiVR3Z19obhrq/E0GmlLJcCJEKvA98B1jSmURwtJE6H0F/klLi/aKShjf3oneYSLpuIqbsuIEOSxkGArvrqXt+BzIqSbg8H9vUlIEOaURoLxH0a9WQlLK85bkaeA2YB1QJITIAWp6r+zMmpXVCCBwnZJD6zekAVP97E94vKgc4KmUok1LS9FEJNY9vRRdnIvXbM1QSGCT6LREIIexCiLhDr4Gzga3AEuDGltVuBN7or5iUjpmy40j9zkzMY5zUv7qbupd3IcOxgQ5LGWJigQi1T2+n6d0DWKcmk3r7DIwp6t6AwaI/u4+mAa+1DBFrAJ6TUr4rhFgL/E8I8XWgGLiiH2NSOkFvN5J80xSa3j9I84oSwpVekm+YhD7ePNChKUNA2O2j9qltROr82r0BJ2WqoaIHGTVnsdIl/m211L2wE51FT9KNkzFlOQY6JGUQC+xtoPaZ7QgdJF07EfMY10CHNGINmjYCZeizTkoi5VvTQQjcD2/Cv612oENSBinv2kpqHtuKPs5I6u0zVBIYxFQiULrMlGEn9Y4ZGFJt1D69TbvfYAiWLJW+IWOSxnf2U//KbsxjnaR+a4a6QWyQU4lA6RZ9vImU26ZhnZJM49L9NL69X918piAjMer/t5Pmj0uxn5BO8tcmo7OqIc0GO/UXUrpNZ9KTePUEGuP34VlVRswbJuHyfDVl4AgVC0WpfWY7wV31xJ+TR9yp2apReIhQiUDpEaET2iB1DhNN7x0g6g2TdN1EdZfoCBP1hql9sohQaTMJX8nHPi99oENSukBduik9JoQg/rQcEr6ST3B3PTX/2ULMFx7osJR+Em0M4n5oE6EKD0nXTVRJYAhSiUDpNfZ56SRdN5FQuQf3o1uIelUyGO4iDQGqH9lMtClEys1TsU5W81kMRSoRKL3KOjmZpBsmEa72UfPoFqKe0ECHpPSRSF0A98ObiXnCJH99CuYxzoEOSekmlQiUXmctSCT5hsmE3X7c/1HJYDiK1PpxP7KZmD9Kyq1TMefGD3RISg+oRKD0Ccv4BJK/NploXQD3I5tVMhhGDiUBGdKSgBqVduhTiUDpM5ZxLpJvmky0PkjNY1uJ+SMDHZLSQ9HGIO7HtiLDMZJvmaqGGBkmVCJQ+pR5jIuk61vaDJ4sIhaKDnRISjdFvWHcj20h5g2TfNMUTJkqCQwXKhEofc4yPoHEqyYQKm6i9ultyIgaxnqoiQUi1Dy+lUhdkOQbJ2HKUdVBw4lKBEq/sE1NJuGy8QR3N1DbMjuVMjTEQlFqniwiXOEl6To1guhwpBKB0m/sc9JwXjSGQFEtDW/uVQPVDQEyJql/cSehg00kfnU81gmJAx2S0gfUEBNKv4o7KYtoYwjPylIMSRbiFmYPdEhKOxrf2Y+/qBbnhWOwTU8d6HCUPqISgdLvnOfmEa0P0Lh0P3qXBdtUdTfqYOT5rBzPJ2XYT8zAcVLmQIej9CFVNaT0O6ETJF45HlNOHHUv7iRY3DTQISlf4t9RR8Mbe7FMSMR10Vg1iugwpxKBMiCEUU/SDZPQx5ta5rMNDHRISotQuYe653ZgzHSQePUEhE4lgeFOJQJlwOgdJpJvmoyMxqh9epu6x2AQiHrD1D69DZ1FT/KNk9CZ1XDiI4FKBMqAMqbYSLxqAuFKLw2v7lY9iQaQjErqnt9BtClE0vWT0MebBzokpZ+oRKAMOOuEROLPHIVvoxvPp+UDHc6I1bjsAME9DSRcMk7dMDbCqESgDApxp+VgmZRE49J9BPc1DnQ4I45vixvPx6XY52dgn6smlhlpVCJQBoVDPYkMiVbqXtihJrXpR5FaP/Uv78aUE4frwjEDHY4yAFQiUAYNncVA4tUTiHrD1L+0S7UX9AMZiVH7/A4QaD2EDOqUMBKpv7oyqJiyHLjOH01gRx2eVaq9oK81LjtAuNRDwmXjMSRaBjocZYD0ayIQQuiFEBuEEG+1vE8UQrwvhNjd8pzQn/Eog5N9QabWXvDufkKlzQMdzrAV2FmHZ2UZ9vkZ6u7uEa6/SwTfA7Yf9f7nwHIpZT6wvOW9MsIJIUi8PB+9w0jd/3Yiw+r+gt4W9Yape3kXhjQbrgtGD3Q4ygDrt0QghMgGLgAePerjxcBTLa+fAi7pr3iUwU1nM5Jw2Xgi1X4alx0c6HCGnYYle4l5IyReWYAwqpvGRrr+LBHcB/wUOHpWkjQpZQVAy3ObwxsKIW4TQhQKIQrdbnefBqoMDpbxCdhPSMezqozgftWltLf4Nrvxb3ITf0aummpSAfopEQghLgSqpZTrursPKeUjUso5Uso5KSkpvRidMpg5zx+DPsFC3Uu7iAVVFVFPRZtDNLy+B2O2g7hT1RDgiqa/SgQnARcLIQ4ALwCnCyGeAaqEEBkALc/V/RSPMkTozHoSLx9PtD5A03sHBjqcIa9hyV5ioSiJV4xH6FWnQUXTL/8TpJS/kFJmSynzgKuAD6WU1wFLgBtbVrsReKM/4lGGFvMYJ/YTMvCsKVe9iHrAv6MO/5Ya4k/PxZhmH+hwlEFkoCem+QvwPyHE14Fi4Iq+/LIDm9YjhA6h06HT6xA6PTq9Dp1Oj9Dp0BuNGM0WjGYzRrMFvdGoxmEfJJzn5OEvqqH+tT2k3jFDDY3cRbFQlIY39mBItRJ3iqoSGixkLEYkFCISDh3zHA2FiEYiyFgMKWPEYjGklOgNBnImTe31OPo9EUgpPwI+anldC5zRX9/9yp/uBLpwt6oQ6A0m9EYzJosVs92BNS4OW3wc1vh4LHYHFocDa1w8dlcC9oRE7K4ErI44hE4Vu3uTzmrAdeFY6p7fgWdNOXEnZQ10SENK8/JiovVBUm6bpu4e7gMyFsPf3IS3sQFfQwO+xnrtdVMjAU8zAa8Xf1MzAY+XoNdLKOAjHPARi3ZtKBWDKZ7vPf1cr8c/0CWCfpUanYVEIIUgJkAKgRQ6ZMvrqE5HTAiiAm05YZBhopEwvuYQvqYAdeUVIA8gZQBkkNYSi9DpsTqc2FwJxCUl4UpLw5maRnxKKvEpaThT0jDb7aq00UXWacmY1yXQtOwgtinJ6J1qmOTOCFd6af6kDNvsNMxjnAMdzpAUjYRpclfTWFVJQ3UVDZUV1FVU0FhVia+xnoC3CVodEkWHEBYQZhBmhDCDiEOIZITBjMFoBPToAJ0EXUxiiMUQMYlOSgQgZAwhJUJKzKG+Ob4RlQguuGgMMhJFxqIQjSGjEe05FoVoFBmNQTRCLBQi5gsQCkQIByKEAlGCYQiFTQSiNoJRA8GogZDRTtBoIWQ0EtQLInqISR9SegkGfATKPdSW7kXGNgDHZn690YItPhlHYjKu9CyScrJIy8slKTsbR2KSShKtEEKQsHgslX9fT8PS/SRdPWGgQxr0pJTUv7EXnUWP83x141hbpJT4m5uoLS2j+kAptaXlNFRV0lxbha+xhpCvgWMv+vQInbPlkYPebEcIG0JYMcR0mKMxTOEI5nAAc9iLSYQw68NYDDEsZonZHMRsCWM06jCY9OhMRoTJ1PKwIoxGhEEPOr1WDSp0oNOhc/RN286ISgSJN97Y8UqdJGMxYs3NRBsatEdjI+HqGryl1XjKQ3hqdPgaTPi9ToL6TPwmGwGTnqBZT0gXJSJ9eJua8NSXULF7C3Cka6QQRky2JGyuNJwp6aTk5ZIzcTRpo3OxxjtHdJIwJFmJOyWL5g9LCJ6UiTk3fqBDGtQCRbWE9jfiumQsertxoMMZUJFQiEZ3FdX7S6nYW0JtaRlN1VX4mtyEA/XI2Jcut4UdoXOi16dhMedhjugwB2NYw1EsoSDmiA+rQ2BL1mFPt2BPd2JLc2FMScaQnIIhOQm904mwWAb9b1YMxREe58yZIwsLCwc6jE6R4TChklJC+/YS3Lvv8LN/335CEQNBUzwhRxKBnDz8Lic+ox5f1I/fW0/QV4OMNnL0PXg6vQWbMxVXegapebmkjcklMSMLV0YmFvvIuDkoFoxSeXchhgQzKd+aPuh/ZANFRmJU/n0dwqAj7buzEPrh/+8U8vuoryinrqKMhopy3MVl1JWV01xXTcjXyLFX9QZ0BhdmeyIWqwurzoQ1FMFWV4ul5AA2Xz2mUBNGg8SSn4+5YDyW8eMxjR6NKTcXY2Ymwjh0kqsQYp2Uck5ry7pcIhBC2IGAlFLd3dMJwmjEPGY05jGjiTvzyOcyGiW0fz+BoiIC27bhLyoi+OkyYj4fAPrERCwzZ8Lk06lzZeL2hHCXltNYXYGvyY23YTul2z4/5rtMtjgSM7NIysoiISMLV3omCRmZJKRnYrQMn5EldWY9znNGUf/ybvyb3dimt3lD+ojmWV1OtDZA8s1ThlUSiEYiNLmrqCsvo76ijPqW59qyUnyN9ceuLBwIvROdPov4tOkkpGeQOiqTRIMkvnIPctMXBDasItasdUvWORxYp03FctkiLBMnYC4owJSbi9AP72E4OiwRCCF0aH3/rwXmAkHADLiBpcAjUsrdfRznMYZSiaArZCxGaN8+fOvX41+/Ad/69YSLiwEQFgu22bOxL1iAed58PPZMqg82ULbzAFX7Smh0VyAj9cRi9QgaiEU8x+zbkZCoJYeWxJCQkUViVg4J6RlDsoeTjEmq/7GBmC9C+o9mq/FyviTqCVF5dyHmUfEk3zRloMPpFhmL0VBdSU3JQWqKD1BTUkxN8QEaKsuJRY9ch+qNNoQ+ARlzIXQJ6E2JJGRmkTFuFGl5CSRnWrHV7Se4bi3ezz/Hv2EjMhAAITCPH491xgys06ZhnT4N05gxQ/L30BntlQg6kwg+Bj5Au9lrq5Qy1vJ5InAacA3wmpTymV6Nuh3DNRG0JuJ241u/AV9hId41qwnt2QuAPikJ+/z52BcswHHKQmR8ItUHm6ja30TlvkYq9tbga6hGxhoQukYsVg/QSNDrJug9clOWwWwmJSePlFGjj3rkYbLaBuiIOy+4rwH3I1uIPyeP+NNyBjqcQaX+jT14P68g7fuzMaYO/r+lr7EB98ED1JQcwF18gNqSg9SUFhMJBg+vY09IxWxPIxZz4ffYkVI78TsSXWSMdZI+xknGOCdJ2Q5iNTV4P1mJZ+UneFevJubRLozMBQXYTpiHfd48bHPmoHe5BuiI+19PE4FRStluZ9fOrNObRlIi+LJwVRXe1Wvwrl6Nd80aojU1AFimTyPu9DOIO/00TOPGAdBcG6ByfyMVuxsp21VPfaVW7WQwhUlID2OP8wA1eGpLcRfvJ+j1Hv4eV1rG4cSQOnoMaWPycSQk9vvxdqTmqSKC+xvJ+OlcdLahU1/blyJ1ASrvKcQ+J42ES/MHOpzjeOpqqdq/h6p9e6nav4fq/Xvx1NUeXm5zukjKHoXNlQEk4m2Mo67STCxqQAhIynaQMbbl5D/WSVyiRevHv2kTnhUf4Vm5kuCOHQAY0tJwnHIK9pNPxjZvLoaEkTvlSY8SwVE7WQ7cI6VcetRnj0gpb+udMDtvJCeCo0kpCe7ciWfFCpo/XEFgyxYAjDk5xJ1+GnFnn4115szDRV1vY5Dy3Q2U72o4JjFYHEayJySQmgsmSwPN7lLcB/fjLt5PfWXF4f7RjsQk0sfmkzYmX3sem4/VETcwB98iVOGl+oH1xC3KwXlu3oDGMljUvbQL3yY3GT+ZM6D3Wkgpaa6tOXyyr9qnPXsbWurxhSAxI4u0MeNIHT0We2IWvqY4KveFKd9VTyigVf8kZdnJKkgguyCBzHwX5paEL6NRfOvW0bzsfZrff59IVRXo9dhmzcJ+ykIcpyzCPD5fdSZo0VuJYB9QgjZO0O9bPlsvpZzVa5F2kkoErQtXVeP56CM8H36Id80aZCiEIT2d+PPOI/7887FMmXzMj8LbGKR0Rz3F22op2VaHv1kr1CVlO8idlMioyUkkZZuoLTlI5d7dVO3bTeXe3dRXlB3ehystg7Qx40gfm0/6uPGkjc3HaOrfk0/tCzsIFNWS/tO56ONM/frdg0242kfV39fhOCmr3yeiDwcDVO3bQ8XunS2PHXjq6wAQQkdSdg6po8cePvEnZeVRUxrkYFEtxVtrD1+YxCdbyJmYSFZBAlnjE7DFH/mbymgU39q1NC19h+bly4nW1iLMZuwLTyb+7LNxnHoq+njVpbg1vZUI1gPzgAeAHOA6YIVKBINT1OPFs2IFTUuX4lm1CsJhjLm5xJ9/Hs6LLsI8duwx68uYpKbUczgpVOxpJBaTmG0GRk1NYvS0FHInJ2KyGAh4PVTt20PVvj1U7t1F5d7dNNdoc0To9HpSR48lc/xE7VEwgbjEvp0GMVLjp/LeQhzzM3FdPLbjDYax2me3E9hZT/pP56B39F1SlFLSUFVx+IRfsXsn7oP7DzfiOtPSycyfQPq4AtLH5pMyKg+j2ULAE2b/5hr2b3JTsqOeSDCKziDIyneROzmJUVOScKXZjruKD+zaRdOSJTS++RaRqiqEzUbcqYuIO/tsHAsXorOrQfQ60luJYIOUcmbL668BPwISpJT9PoKVSgRdE21spPmDD2h6eynezz6DWAzrjBm4Lr+MuHPPQ9/K3YqhQISS7XXs31TDgS01BL0RdAZB9vgERk9PJm9aCo6EI1f+vsYGKvbspHzXDsp3badyz24iIa2hLz4ltSUxTCBz/ERSRo1G18vd8epf3Y13XRXpP56DIWH4dJXtilCZh+oHNxB3eg7Os/N6dd/RSJiqfXso3V5E2Y4iKnbvxN/cBIDRbCF93Hgyx08gI7+AjHEF2Jyuw9s21wXYt9HN/o1uyvc0ImMSR4KZvKnJ5E5JIrsgAaP5+P8PEbebxrfepnHJEoLbt4Nej+Pkk3EuvhjH6aejG0ZdovtDbyWCb0gpHz7q/WzgDinlzb0TZuepRNB9EbebxjffouGVVwjt3YuwWok/91xcl30F6+zZrdanxqIxKvc1sn9TDfs31dDo9gOQMdbJ2NmpjJuVit11bHVQNBLBfWAf5bu2U7ZrB+U7tx1uEDSaLWTkjydrwhRyJk8lY1wBBlPPrl4jjUEq/7YW2/RUEq8Y36N9DVU1TxYRPNBExs/morP2bNCAcCBA+e4dx5z4DyX2hIwsMsdPJCO/gMzxE0jKyUWnO/ZE3lDtY09hNfs2unEXN7dsZ2fMjGTGzEghJTeu1f9rMhbD99ln1L/wIs0ffgiRCJYpU3BefDHxF5yPISmpR8c1kvVKIhhMVCLoOSklgU2baHjlVZqWLiXm9WIaM4aEa6/BufiSVksJh7arr/Cxb6ObPeuqqC3zgoDMcS7GzU5lzMwU7G00UDbVVFO+czvlu3ZQtmMb1Qf3gZTojUYyxhWQPWkqOZOmkJFfgNHc9au9hrf34VlVRtoPhkaXyd4UPNCI+6HNxJ+bR/ypXe9KG/B4KNtZpJ34txdRtX8PsWgUIXSkjBpN9sTJZE2cTFbBJOyu1nveeOoD7C6sZk9hFdUHtZN/2uh4xsxIYcyMFFxpbf9NIvX1NL76GvX/e5HwwWL0TifOr3wF1+WXHVeNqXRPT7uPPkg7YzdLKb/bs/C6TiWC3hXz+Wh6bxn1zz9PYPNmdA4HzksvJfHaazDl5bW7bV2Fl73rq9ldWE19hZYUssa7KDghnbEzUzG1c2Ua8Hoo27GNkm1bKN22ler9e5Eyhk5vIH3ceHImTSV70hSyxk/s1J3RUW+YyrvWYilIIOnaiV39ZxiypJS4H95MpDZA+k/moDN1XO0W8vso3VFE8ZZNFBdtxn1wv5aUDQbSxo4ne+JksidMJrNgImZb2/Xv/uYQe9dXs2ttFRV7tHmlU3LjyJ+bxrjZqcQltv938xcVUffUUzS/8y4yHMY6ezYJV32VuLPPRmdWo8v2pp4mgqNHavs9cOfRy6WUT/U4wi5SiaDv+Ddvpu6ZZ2h6510Ih7GfspDE62/AfvJJHXbDqy33sGddNbu/qKLR7cdg1DF6RgoTTkwne0Iiug4mkwn6vJTt3Ebptq2UbNtC1b49yFgMnV5PRn4Bo6bNZNTUmaSPzW+zjaHx/YM0Ly8m9TszR8zE7P4dddQ+WYTrkrE45me2uk4kFKJ81w5KijZxcOsmKvfsQsZi6A0GMsdPJGfyNLInTSF93PgOe31FIzEObKlhx5pKDm6tRcYkCek28uemkT8nrd0rf9CqfzwrV1L3xJP4Pv8cnd2O85JLSLjqq5jzB999D8NFr1UNHd1gPJBUIuh7Ebeb+v/9j/oXXiDqrsE8YQJJt95C/DnnIAzt1z9LKana38SOzyrZU1hF0BfB7jQxfl46BSemk5TZuRN0yO+jfOd2irdtoXjLRqr27wUpMdvs5EyeRt50LTG40jMObxMLRKj861oMKTZtEpZhNMZOa2KhKNUPbkBGpTbURss8xFJKaksOsm9DIQc3b6B853Yi4RBC6Egfm0/OlGnkTplOZsHETnf3rSltZvvqCnZ9UUXAE8bmNFFwQjoFJ6STmNnx/BqxYJDGJUuoe+JJQvv2YUhPJ/H663FdeQX6uIG9H2Uk6M1EMCD3DXyZSgT9R4ZCNL71NrWPPkpo3z6Mubkk3Xwzzksv6VTRPRKOcnBLLTs+q6R4ay2xmCRtdDyTF2Yxbk4qxk5UYxzia2qkeOsmDm7eyMHNG2iu1bqsOlPTGDNrHvknLCBrwiQCm2qpe3EnjgWZOC8cM2yntZSRGHUv7cK/2U3yzVMwjY2neMsmdn+xmv0b1h3+90nOGUXu1BnkTplG9sQp7Vb1fFnAG2bn55XsWFNBTYkHnUEweloKExdkkDMxAZ2+43F5Yj4f9c+/QO0TTxCtqcE8aSJJN91M/LnnDKnRO4c6lQiUHpOxGM3Ll1P7yH8IbNmCISWFpG98A9eVV6DrZI8ff3OInZ9XUvRJOQ1VPkxWAwUnpDN5YSZJXazG0Rqtyzi4eQMHNq3n4JaNRMNh4lNSmXvRZYwKFOBdU4Ex20HcSVlYpyYPmykaY4EI/qJamleWEqny4Tgzh13eQja89xae2hqMFiujps5g9Mw5jJ45u1v3cVQdaGLrx6XsLqwmGo6RkhvHxAUZ5M9Jw+Lo3Mk75vVS//zz1D7+BNG6OuwLTiTpttuwnXCCutt3APS0jaCZI43FNsB3aBEgpZT9fhufSgQDR0qJ77PPqPnnv/AVFmLIyCD5m9/E9ZVLO311J6WkYk8DW1eWs3dDNbGIJH2Mk8mnZDJudiqGbowkGvL72Lt+LZuWvU3Zjm2MmTWXs864Dc9HZURq/OhsBqxTk7HNSMU0Kn7IlRJiwSiB7bX4NtcQ2FUHEYkhxYr9jCzefOUeKnbtYNS0mUw74xzGzD6hZQrErgmHouxeW8XWj8twFzdjMOspmJfG5FOySMnpfNVNzOej7plnqXv8caINDdhPPpnk22/HNmvAa5VHNNV9VOl1Ukp8a9bgvv8B/Js2YczOJvmOO3BefFGXxm73e0LsWFPJtlVaKcEaZ2TyKVlMOSWrzW6oHcW14d03WfHkI8w6fzGnXn8Lwd31eNdXE9hWiwzH0DvNWKenYJuRgjFj8M4dLWOS4N4GvF9U4t9eB5EYungTtqnJWKelYMqJ4+0H/squzz7l/O/8iAknLerW9zRU+9jyUSk7P6sk6IuQmGlnyilZFJyQ3m6vr+PijURoePkV3P/8B1F3DfZTFpJy++1YZ8zoVlxK71KJQOkzUkq8K1fivv8BAtu2YZ4wgbSf/RT7iSd2eT+lO+vZ/GEpB7bUoNMJ8uemMf2MnC5djR7ywaP/ZPMH73Hj3f8gKTsXOOqqeqObwK56iEmM2Q4SFo/D1I3v6EuBPfU0vLGXiLulNDM9Bdv0FEy5R0ozJUWb+d8ffsmCK67lxMuv7tL+tVJZIxs/KGb/Zu3fe+zMFKYsyiJjnKtLyVFKSfP77+O+9++EDhzAOmsWqT/+sSoBDDI9rRpaDGRLKf/Z8v5zIKVl8c+klC/1ZrCdoRLB4COlpPmdd6i+517CZWU4Tj2V1J/+BPOYrg981lDlY/OKUravqSASjJKZ72L6GTnkTUvusAvqIb6mRh777q3kTpnG4h//+rjlUW8Y/2Y3TStKiPnCJF03CeuEwTHMtnddFfUv78KQaCH+zFFYpyQjjMe2b0gpee7XP8JbX89N9z3U6Z4/0WiMveur2fRBCdUHm7HYjUxZlMWURd0rgfk3baLqT3/Gv2kTprFjSf3RD3GcdtqgLWWNZD1NBJ8CV0kpS1rebwTOAOzAE1LKM3o33I6pRDB4xYJB6p9+mpqHHibm95Pw1StJ+e53uzUBSNAXZtunFWxeUYKnLogrzcasc0Yx/oQ09J3orfLpi0/z2Wv/46Z7HyIxM6v1eH1h3I9tJVLtI+37szAkWbscZ28KlTZT/a+NmMe4SLphUps3hxVv3cRL//crzrzlDqafdV6H+w2Homz7pJyNHxTjqdf+LaefkUPB/PQu9dw6JFJbS/W999L4yqvoU5JJ+e53cV16aYddi5WB09NEsFZKOfeo9/+QUn675fVnUsr5vRptJ6hEMPhF6uqo+cc/qH/hRfQuF6k//QnOxYu7daUYi8bYu8HN+vcOUlPiwZFoZtbZo5i4IANDOycxb0M9/7njJqacdjZn3nJ727E2BKn6+zpMOXEkf33KgF3Nyqik+sH1RH0R0n8wu93xgl7502+pPrCPW//xeLvjNIX8EbZ8XMqm5SX4m8Nk5ruYcVYueVOSutVgLiMR6p97HveDDxLz+0m88QaSv3V7m0OSKINHTyevP2ZgkUNJoEUKitIKQ2Ii6b/9La4rrqDyd7+n4ue/oPHV10i/87ddHjtGp9eRP0cbsuDg1lrWvXOQlS/sYu3b+5lxZi5TFmVhshz/X9nuSmDiwtMo+ng5J331OqxxrXdwM7jMOM8eRcOb+wjuacCSPzCzWPnWVxGu9JF03cR2k0BdeRkHNq1nwZXXtpkEAt4wmz8sYfOKUoK+CLmTE5l9Xh6Z41zdj2/DBirv/B3BXbuwL1hA2q9/1a2qP2Xw6Uwi+FwIcauU8j9HfyiE+AbwRWe/SAhhAVaiTXxvAF6WUt7ZMvfxi0AecAC4UkpZ39n9KoObZeJERj3/HA0vvUz1vfey75JLSbrpJpJv/1aXhxEWQpA3NZlRU5Io393AuncPsua1vWxYVsysc0Yx9dSs40oIsy+4hK0r3mfTsqXMv+yqNvdtPyGD5pWlNH1QjLmLjaW9QUZiNH1QjCknDsvk9kfY3LTsbXR6A9POOPe4ZUFfmA3Litm8opRwMMro6cnMOT+P1FHd7+Ud8/movu8+6p9+BkNGOlkPPkDcmWeqdoBhpDOJ4AfA60KIa4D1LZ/NRjuhX9KF7woCp0spPUIII7BKCPEO8BVguZTyL0KInwM/B37Whf0qg5zQ6Uj46pXEnXkG1X+7m9pHHqH5gw/I/POfsE6f3vX9CUHWeG32qqoDTXyxZB+rX93DxuXFzDkvj0knZ6JvuXksOWcUo6bNZNP7S5m7+HL0bdRhC4OOuNNyaHh9L8HdDVjG92+pwLehmmhjkITL2p9aMRwIUPTxcsbPP+mYUUDDoSibPyxhw7Jigr4I42anMvu8PJKzezbeknfNGip+81vCpaUkXHMNKT/8oaoGGoY6c6ulW0q5APg/tCv2A8AfpJQnSimrAEQnLg2kxtPy1tjykMBi4NDAdU/RteSiDCGGpCQy//Jnch57lJjfz4Grr6H6nnuIBYPd3mdaXjwXfXcGl/5oJs4UKytf2MWzv/2MHZ9VIGNa+9eMcy7EU1/HvnXtF2Dtc9LRO800fXCQ/uxWLWOS5o9LMWY5MOe72l13+6qPCPq8zDj7AkAbAG7LR6U88+s1fPb6PtLHOLnyV3M559YpPUoC0eZmKn7zG4pvuhmh1zPqmadJ/+1vVBIYpjpTIlghhHgFeENK+eGhD4UQJuBk4EZgBfBkRzsSQuiBdcA44J9Sys+FEGlSygoAKWWFECK1jW1vA24DyM3N7UTYymDlOOkkxry5hOq77qL2P4/SvGIFWX/7G5aJ3R86OjM/gUt/NIuS7XV8/sY+lj+5nS0rSjnp8nzGzJpDXFIKG5e9Tf4JC9rcx5FSwR6Cexux9KA+vSv8W2uI1PhJvHZiu6UBKSUb33uLlFGjyRg/gf2ba/j0pd00uv1kjHNyzm1TetQGcIhv/QbKf/ITwpWVJN16C8l33KFmAxvmOtNryALcDFwLjAYaAAugB5ahndA3dulLhXABrwHfAVZJKV1HLauXUrZbLle9hoYPzyefUPGrXxOtryf1Jz8h4frrelz3LGOSXWur+Oz1vXjqg4ydmYLFtol1b73ATX9/uM2upAAyHKPir19gTLOTcsvUHsXRqVilpPqBDchIjLQfzG63J0/pjiJevPNnLLjyNtxluZRsqyMh3caCy8YxakpSz//dIhFq/v0QNf/+N8bMTLLu/pu6K3gY6VGvISllAPgX8K+Wuv1kwC+lbOhuQFLKBiHER8C5QJUQIqOlNJABVHd3v8rQ41i4kNFvvE7FL35J1Z/+hHfNGjL+9EcMCd2voxc6QcEJ6YyZmcLG94tZ/95BohEXQuhYt3QJZ93yrba3NeqIW5hN49L9hEqa+/yO4+CuesIVXhIuz++wO2fhm2+gN1pY/4EFs7WJk6/MZ8qirE7dU9GRcFkZZT/+Cf4NG3Auvpi03/wGvWNkzOegdK6N4DApZVhKWdGdJCCESGkpCSCEsAJnAjuAJWjVS7Q8v9HVfStDmyEhgex//4u0X/4C76pV7F98Cb716zvesANGk565F4zmuj+cSMG8MeiME9n8/lKKPl7e7nb2E9IRVgNNK0p6HENHmlaUoHeasM1otUYU0EoNn7/+MnsLVyP005iycBTX/mE+00/P6ZUk4PlkFfu/chnB3bvJvPtuMu+6SyWBEaY/x+XNQGtv2AysBd6XUr4F/AU4SwixGzir5b0ywgghSLzhBvJefAFhtXDwhhupf+GFXmm0tbvMnPG1SZx64zfQGbJ4919/Z+mDd+P3NLe6vs5swLEgk8C2WsJV3h5/f1uCBxoJHWjCsTC7zSGym9zVvPzH37Dq+SfRGfM5+5tfZ9E1BVgdnRv6uz0yFqPm3/+m5LbbMKSlMfqVl3FeeEGP96sMPf12P7iUcjNw3ChUUspatCErFAXLpEmMfuklyn7yEyp/93sCRUWk/eY3nZ7zoD1TT81j4wfXEA1+wc41n3Bwy0bOuvXbjJt7/M3xjgWZeD4ppfmjUhK/WtDj725N04cl6OwG7PPSj1smpWTrivf56L//QUpIHXsx4fBExs/LaGVPXRdtaqL8pz/D89FHxF90ERm//x06W/tTTCrD1/CYqUMZVvTx8eT8618kffMbNLz0Mgevv55wVc+bjvQGHdPPyCPgn8UF3/8jdlcCb9z9/3jnn/cSDgSOXdduxD4vA9+maiJ1gTb22H2hMg/BXfU4Ts46bjwhf3MTr//1Dyx7+AHSRo/jK7+4m6b6cUw+pXfaA0KlpRy46mo8q1aR9utfk/nXu1QSGOFUIlAGJaHXk/r975P1wP0Ed+/hwNVXEdy7t8f7nXBiOnqDjsp9Jq79073Mv+xqtn/yEc/+6ofUlh3bJhC3MAuEoHllaY+/98uaVxQjLHocJx472XzFnp08/fPvcXDzBk678Vau+M0fKd4WRgjB5JPb7u3UWf5Nmzjw1auI1NaS+/hjJF53rbpDWFGJQBnc4s8+m1FP/xcZCnPgmmvxrVvXo/1ZHSbGzkph52cVxKI6TrryWi771R/wNTXy7C9/yM41qw6vq3easc9Ow1tYSbQp1NNDOSxc7cNfVIvjxEx0R42RtHHZUl747c8QQnDVH/7GrPMXE41Ktn9awZjpyTgSuj5M9NGali3j4A03orPZyHv+eezz5vX0UJRhotOJQAihuhEoA8I6eTJ5LzyPITGR4ptupum9ZT3a3+RTsggFouwurAJg1NQZXH/X/aTk5vHWfX+h8K3XDq8bd0o2RCXNq8p69J1Ha15RgjDocLRc4Usp+fiZx1n+2L8YNXU61/3lftLH5gOwZ101AW+YKadm9+g7655+hrLvfR/LhAnkvfgC5jGje3wcyvDRlRLBBiHEJCHE4UsYIcT4PohJUY5jys5m1HPPYpk0ibLvf5/6F17o9r4yxjpJzLRTtPLIyT0uMZkrfvsnxs8/mY+ffoyVzz2JlBJDshXrtBS8n1UQ84V7fByRugC+TdXY56WjtxuJxaK89+/7KXzzVaaffQGX/Oy3WB1H7l3Y8lEZCek2ssa7uvV9UkpqHnqIqj/+kbgzzyD3qScxJA6OCXiUwaMriSAR+CdQIoTYIoT4H7C0b8JSlOMZEhLIffIJHIsWUfm731P39DPd2o8QgskLs6g+2Ez1waYj+zcaueB7P2H6Weex9o2X+eKNlwGIOzUHGYriWVPR42No/rgEhCDulGyklHz01KMUffwBJ15+NWfc/E10uiMNx9UHm6g+0MSURdndqseXUuK+917c992Pc/HFZP3972qoCKVVXUkExVLK06SUGcB5wH+Ar/VJVIrSBp3FQvYD9xN31plU/fGP1D7xZLf2UzA/HYNJR9En5cfuX6fnjJu/xYSTFrHq+afYtvJDTBl2LBMT8XxaRiwY6Xbs0cYg3sIq7LPT0DvNFL75KhvefZPZF1zCgiuOb7Td8nEZBrOegvnHdy/tiJSS6r/8hdr/PIrrqq+S8ec/q9nDlDZ1JRE4hRALhBBOKWWplPJ9KeWqjjdTlN4lTCay7r2XuHPO0Qau60YyMFsN5M9NY9faKkL+Y0/uQqfjnG99n9wp03jvoQeo2r+X+NNzifkieFaVt7HHjjV9WAxoJYwDmzew8tknKDhxIYuuu/m4dQPeMLvXVlEwLw1zO5PUtMV93/3UPfVfEq6/nvQ770ToVL8QpW1d+d/hAH4MfCGEOCCEeEcI8bc+iktR2iWMRrLuuZu4c8+l+q67aHj1tY43+pIpp2QRCUbZ9UXlccsMRiMXfv/n2OLjWfrg3ejSLVgmJdG8srRbbQWRWj/etVXY56YTMUV4919/JzErh3Nu/36rJ+ntqyuIhmNMWdT1RuKahx6m9uGHcV15JWm//IXqHqp0qCuJ4EQp5VeklAXABOBXwOa+CUtROiYMBjL/ehf2BQuo+M1vaF7e/hhCX5Y6Kp6U3Di2rixrdSgLa1w8Z3/ze9SVlbDu7ddxnj0KGYrS/HHX7yto+qAYdIL403NY/dKz+BoaOP/bP8JoOr5LqIxJtq4sI2Ocs8tzCtQ98yzu++4j/qKLSL/ztyoJKJ3SqUQghJgA5B3qQiqlDEgp1wPuvgxOUTqiM5nIfvABLJMnU/aDH+Jbu7ZL209emEltmZeq/U2tLh89Yzbj5p7IZ6++QNgaxjY9Bc/qcqLNnb+vIFzlxbexGseCTBqaq9n43ttMO+s80saMa3X9ku11NLn9TO1iaaD5wxVU/fGPOE4/ncw//wmh13e8kaLQiUQghPgu2oig3wG2CiEWH7X4T30VmKJ0ls5uJ+fhhzBmZ1P67e8QKi7u9Lb5c9MwWvRsXdn2fQILr/ka0VCYwrdeI/7MUcho7HB9f2c0vX8QYdITtyibz1/7H3qTkQWXX93m+ls+LsMaZ2TMzJROf0dg2zbKfvxjLJMnk3XP3aphWOmSzpQIbgVmSykvAU4FfiOE+F7LMlXuVAYFQ0ICOQ/9GwmU3H47UY+nw20ATBYDBSeks6dQu3GrNYmZWRQsWMjG994mYolin5OO94tKIjX+DvcfLG7Cv7WWuIVZNHtq2bHqY6afdT42p6vV9Ztq/BzYUnPMvMsdCVdWUvLNb6F3Osn+1z/RWa2d2k5RDunM/zT9obmGpZQH0JLBeUKIe1GJQBlETLm5ZN9/H6H9Byj/yU+R0Wintpu8MItoJMaOdu4TmHPRVwgHAxR9/AHxZ45C6HU0vLWv3f3KmKRhyV50cSYcC7PY9P5SEDD7/MVtbrPlo1KEEEw5pXPVQrFQiNI7vk3M6yXnoYcwprY9r4GitKUziaBSCDHj0JuWpHAh2kxlfT+Xn6J0gX3+fNJ++Qs8K1bgfvDBTm2TnO0gfYyTok/KD094/2Vpo8eSOX4iG5e9jc5hIP7MXAI76vDvqGtzv77CKsKlHlwXjCYqomxd8T75c08kLim51fVDgQjbPq1g7KyUTo8rVP2XuwgUFZH5t79iKVA3+ivd05lE8Dug6ugPpJQRKeUNwCl9EZSi9ETCNdfgvPwyah9+BM+nn3Zqm6mnZtFQ5aN4W9sn9ulnn09DZQVlO7fhWJCJIcVKw5t7keHjSx4xX5jG9/ZjyovHOj2FfevWEvA0M+3M89rc/64vtHsapp2W06mYm955h/rnniPxppuIO/30Tm2jKK3pTCK4BlgqhHhBCPE1IcTh2xyllJ37lSlKPxJCkP6rX2EaO4byn/6MiLvjzm1jZ6Vid5rY9GHb01Pmzz0Ro9nC9k8+Qhh0uBaPJVoboPHdA8esJ6Wk/vU9xPxRXBePRQjBztUrsbsSyJnSeiFaSsnmD0tIyY0jfUx8h/GGDhyg4te/wTpjBqk//EGH6ytKezpMBFLKb0opZ6KVDBKAJ4UQa4QQfxJCnCKEUH3UlEFHZ7WSde+9xDweyn/xyw6nvNQbdEw9LZuSbXXUlrfe0Gy0WBg3dz67PltFJBzGMi4B+4kZeD4tx7f5SLLxri7Hv7mG+LNyMWU6CPq87NuwlvEnnnzMWEJHK91RT32lj2mndzyukAyHKfvhjxAGA1n33oMwGjv411CU9nX6hjIp5Q4p5d+llOcCpwOrgCuAz/sqOEXpCcv48aT+9Cd4V62i4eWXO1x/8slZ6I06Nn/Y9g1jExeeRsDrYf/GQgBc54/BlBdP3Ys7afrgIA1L99Pw1j4sk5OIW6RV8exZ+xnRcJgJC9quSd28ohRrnJH82WkdxlnzyCMEtm0j/f/9H8bMzA7XV5SOdGsAEimlX0q5VEr5HSnlnN4OSlF6S8LVV2M74QSq/3IX4fL2xwmyOIwUzE9n5+eV+D2t3zA2auoMbE4XOz75CABh1JF8wyQsBYk0fVCMZ2UptllpJF1VgNBpV/Y7V68kPiWVjPwJre6z0e3jwJYabSpKY/s/ycCOHdT8+yHiL7yQ+LPOav/gFaWTupUIWtoLnhZC/FcI8dfeDkpReovQ6cj44/9DSknFb+/ssIpo+mk5RMMxila2njR0ej0FCxayd/0XBH1e7TObkaTrJ5L+07mk/2IeiVeMRxi1KiBfUyMHt2ykYMEpbVb5bFlRhk4IppzS/lSUMhKh/Je/RO9ykfarX3Z06IrSad0dknCNlPL6lp5Dd/VmQIrS20zZ2aR+//t4V62iedn77a6bmGknZ1IiWz4uJRqJtbrOhAWLiIbD7C08UisqhMCQaMHgPLbb5+7PVxOLRtusFgr5I2xfXc7Y2anYne13Ga1/4UWC27aT/utfYUhIaHddRemK7iaCxUKI7wghxkspa3s1IkXpAwnXXI15wgSq/vxnYj5fu+tOPyMHX2OIPeuqW12ekV9AXHIKO9d80uH37ly9koTMbFJGtT41ZNEn5YQCUaaf0X6X0UhNDe7778e+YAFx55zT4fcqSld0NxFcD+wFLhNCPNqL8ShKnxAGA+m//Q2Rykpq/v3vdtfNnZhIQrqNTctLWq1KEkJQcOJCDmzaQKCdoSw8dbWUbN/KhDaqhaKRGJuWF5NV4CItr/0uo9V/u5tYIEDab36tRhRVel1XJq+/XwixTAjxHvA9oFxK+Wcp5S19F56i9B7brFk4Fy+m7smnCJe1Pcic0Ammn5GDu7iZsp31ra5TMP9kYtEIewo/a3M/O9esAikpWLCw1eW7vqjE2xhi1tmj2o3bX1RE4xtvkPS1r2EerSadV3pfV0oE24G/AfcD1cAzQohv90lUitJHUr73XRAC9wPtDz9RMD8dW7yJde8ebHV52th84lPS2NVO9dDO1StJyRtDUtbx1T4yJtmwrJikbAc5k9qfTN59zz3oXS6Sbru13fUUpbu6ch/BQy3TUy6VUt4NzAG+0dnthRA5QogVQojtQoiiQyOYCiEShRDvCyF2tzyrVjClzxgzM0m4/joalywhsGNHm+sZjHqmn5FD6Y76Yya4P0SrHjqZg1s24vc0H7e8sbqSij0722wkPrClhvpKH7POzm23qsfz6ad4V68h+VvfRB8X14kjVJSu63IbgRDim0KIu9FKBq3P5tG6CPAjKeVEYD5whxBiEvBzYLmUMh9Y3vJeUfpM8q23oouLo/rvf293vSmnZGG2GVjfRqmg4MSFxKJR9qxdc9yyHZ+uPLxOazYsKyYu0cK42W2PFiqlpPqeezBmZeG6uu35CxSlp7rTWLwUrZooG/hzZzeSUla0zGqGlLK5ZR9ZwGLgqZbVngIu6UZMitJpepeLpK9/He/HK/Fv2drmeiargSmLsti70U19pfe45amjx+JMS2fn6uOrh3auXknG+Ak4U4+/U7hiTwMVexuZcVYOOn3bP0HPio8IbttO8re/jc5k6uTRKUrXdaWx+H9CiIlSymIp5WPARcAfu/OlQog8YCba8BRpUsoK0JIF0OolkhDiNiFEoRCi0N2JQcQUpT0J116DLj6emocfane96afnYDDoWL/s+BnJDvUeKt66CV9T4+HPa0uLcRcfaLNaaP2yYix2IxMXtD08hJSSmoe0WdecF17QyaNSlO7pynx2zwAvCq1Ccx3gAFq/46YdLfMevwJ8X0rZ1NmucFLKR4BHAObMmdP+7aFt2HNgBTqhR+gM6HQGdDojQm9AJ/To9Np7nc6A0BnR6U3odHrtRiFhwKg3YhAG1XVvmNA7HCRefz01//wngZ272hzL3xpnYuJJmRStLGPehaOJS7Qcs3ziSYv44vWX2P7JCmZfcAkAWz/6QLsDuZVqodoyDwc21zD3gjyM5rbHa/R+uprA5s2k//73alC5YSYaixKOhYnEIoRjYaIyioyEiUX8yFiIWDSKlDFiMoKMRYjJKDKmPXQ6A3k5C3o9pk4nAinlEmCJEGIaMAOtNLG0K18mhDCiJYFnpZSvtnxcJYTIkFJWCCEy0Hok9YnLP/oO0R6cyIUEI2BEYBI6jAiMwoBRp8ekM2LUGTHpTRj1JkwGC1aDDaspDqvZicXixGq0YTVYsegtWA1WrEYrVr1V+8xgwWaw4TA5cBgd2Iw2dKK7t3konZF4/XXUPfkktQ8/RNa997a53oyzcihaWcamD0o4+cr8Y5Yl5+aRMa6ALR8uY9b5i4lFIxR9vJwxs+Zhdx3f72Ht2wcwWvRMO739G8hqHvo3hvR0nJde0q1jUzpPSokv4sMT8uANe/GEPXjCHgKRgPaIBvBH/IdfB8J+/MEG/MEmAqFmAiEvgYifQDRAOBYmHI0QlhHCsQhhGSNCjLCMEUYSBmI9uJbMiEqW3dx2dWZ3dXmGaynlZmBzV7drKUk8BmyXUh79q1sC3Aj8peX5ja7uu7PuGn0ZMhYlJqPEYlGkjBKTscPvYzLW8lnL8miIaDRMNBoiFAsRjoYIxcLaHziqPYdivpbnCGEhWh7gFYIKIfALHX6dICAEfl3nT+wCgd1oP5wYHEYHDpODOGMcdpOdOGPc4WXx5nhcZhcuswunyYnT4iTOGKdKLx3Qu1wkXH0VtY8/QcoPyzBltz7WT3ySlfx5aRStKmPWuaOwxR9bXz/1jHNY9vADHNi4Dk9DHf6mRqadcfzdv7XlHvZuqGb2OaOw2Nu+yvdv2YK/cB1pv/i5ahvohJiM0RxqpjHYSGOwkYZgA40h7XVTqOm4E/zR770h7VnS+UoGayyGRUrtEdOerTKGXUrMMYkRtFoEnQGjzohBb8SoN2HUGTHqDBh0Rox6I0adCYPOiEFvQuhNLbUUBoTQIYQOXcuzEHp0QodO6LGa+6bnWJcTQQ+chHZH8hYhxMaWz36JlgD+J4T4OlCMNrR1nzhn0e/7atcQi0GoGYLNEGiCYJP22t8Avlrw1RLz1hD0ufH7a/H76wj46/EHG/ELiV8IfDodXiHw6HR4zDY8ZgPNYYHXEKZZ30S9TlAiY3hkGE/ETzAabDMcvdDjNDuJNx1JEsckDLOTBEsCSZYkEi2JJFmTcBgdIy55JFx3HbVPPkX900+T9ou2O6zNPncUuz6vZMOyg5x0+bGlgkmnnMYXb7zEiqf+Q8jvIyO/gLwZs4/bR+HSAxhNemacmdtuTHX/fRqd3Y7zssu6d1BDmJSS5nAztf5a6gJ11AXqqPXXUh+spynYpJ3kg42HT/SHTvYx2XYttdVgxaG34NCZcCCwxyQp0TCOcBBHyI896MMRDuCQMRwxiT0Wwx6TWNFhsSZitSZgsSRgsSZhtiUhbElgSwRrYstzAlicYHKA2QFGO3Thom8w6LdEIKVcRduT3Z/RX3H0GZ1O+89gcYKzjVUAa8vjMCnBXw/NldBc0fJc3vJcCU3lUFMJniqQR0+JKAjHZ9LsyqIxPoNGexKNNicNljgajFYaCdN41A+nwlvBjvodNAYb8Uf8rcZn0plItCZqieGoBHHMsyWJJGsSCeYE9G1MsjKUGNPTiT/nHBpefpnk73wbvcPR6noJ6XbGn5DOlo/LmHFW7jEDxOkNRk7/2jdYcs+f0BkMnP61bxyXUOvKvexZV82sc0ZhcbRdGghXV9P07rskXHVVm7EMNZFY5PAJvTbQcoL3t5zkAy2f+esOL4vEIq3ux2F0HHNxk2nPxGl24jQ5cQkDzpAfV6CZeE8NruZqnM1VxDWUYfB+qaFf6CE+E5zZkDABHOkQl/al53TtBD9CLoz6s0SgtEYI7arClghpk9peLxrREkRDMdQfhIZijA3FJDYcJLFkPTSVwdFXRRYnJOVDcj4kTYDscdrrxLEEdYKGQAP1wfpjfoCHfqi1gVpq/DXsrN/Z5g9TJ3QkWZJIsaWQak0lxZZCijVFe29LPfw60ZI46Ns6Er92I01vv03jK6+QeOONba4394I8dn1Rxbp3D3LKV49tXB49cw53PP4CCIGhlcbdwqX7MZj0zDiz/baBhhdehEiExOuu7d7B9KOYjFEXqMPtc+P2u6n2VeP2uan2tzz7qnH73dQF6lq9YjfpTIcvMFJsKRQkFhy5CDnqgiTJmoTT7MQYCUPNTqjZDbV7wL0Xar+Aun1aCfwQnRFcudqjYCo4c8GVA84c7TkuE/Tq1Hc09a8xVOgNR/5z5518/PJoGBpLoHYf1O5u+bHshn0fw6bnj1pRYE7IIy1tMmlpU7TkkzYFEvKglSt8KSVNoabjruhq/DXU+Guo9lVT4a1gc81m6gLHT/xuEAaSrEnHJIejX6dYtfcus2vAqqWsU6dinTWLuqefIeG66xD61ks6zhQbE09Mp+iTMmaelXtcDyJDG/X5teUedq+rZtbZo7A62q7zj4VC1L/4Io5FizCNan/8ob4kpaQh2HD4RH70Sf3Qyd7td1PjryF6TClVk2hJPPw3npQ06fDfOcmadExJ02awtf43j4S0E331NqheDu4d2uu6/XC4Ll9ov4WksZAzD5LGQeJY7b0zR53ou0j9aw0XeiMkjtEe+Wceuyzo0X5YtXu0BOHeDlVFsHPpkVKEwQqpE7XEkD4NMmdB+hSE0aoVv81ORjvbH/AsHA1ryeGoK8JDycLtd1PcXMy66nU0BhuP29aoM7aaLFJtqceUOvqqHSPx+uso+8EP8XzyCXGnntrmerPPz2PHZ5UUvnOA0649dsaxcDAKAoymYxPJmlf3YrIYmHlW+20DTW8vJVpbS+IN13f7ODrij/ip9lUf96jyVR1zwg/Hwsdt6zK7Dv8txrrGHve3SbWlkmRJwqjvQnfXkBcqt0D5RqjYpD1qdsKhUqjQayf59Gkw7SpInQDJBZA4Ggztz9+gdJ5KBCOB2QGZM7TH0UK+I1dbVUUtyeFd2PCMtlzoIXUSZM2EzJlackidBIbWr2qNeiMZjgwyHBnthhOMBnH7tCvKKl/V4WRx6GpzT8Me1pSvwRM+fohnq8F6JFkcdQL6cgKxGW1d+ieKO/NMDCkp1D/7XLuJID7JyqSTM9n2STkzz8zFlWajodrH52/sY8+6aoSAiQsymHvhGBwJZkp21HFway0nfmVsu20DUkrqnv4vpnFjsZ14YpdiB61vel2g7vBJvbWTfbW/mubQ8eMi2Qw2Um2ppNpSmZU2ixRbCmm2tGMScbI1GbO+hyfesB/KN7Sc9DdqzzW7OHyVb0+BjBkw/mxInaxdmCTnqxN+P1CJYCQz2SBrlvY4REqt0bpsfcuPdj1sfxPW/1dbrjdD+lTInQ+5J2rP9uQufa1ZbyY7LpvsuOx21/OFfcdVRxxdD11UW0R1STWBaOC4bR1Gx3HJItGSeLjH1NG9qOJN8RiMRhKuvQb3fffjKyzENqftqbhnn5bGzjXlLHtoPdbEOEq316EzCGacmUMkHGPbqnJ2fl5F7uRE6ksacMQbmDy3/RFGm997j+C27WT8v/9DCEEwGqQ51ExzqBlPyEN9sP5IG85RVXSHetfUB+qP6wKpF3qSrEmk2dLIc+YxL2Pe4RP+4Yc1FYepjxqlm8qh5HMo+UJ7rth05Eo/LgMypsPkS7ULlIzp2mcjpHF2sBEdzeE6GM2ZM0cWFhYOdBgjh5TQcPBIcihdq70+1H01adxRieFErXqqn37QUko8YQ9un1ur3mgncbTVGwXAbrTjiln5wz/qiZdmTFlZCL1Bay+Qkpjfrz2amog2NlKccwZ7xn6FOIdg3IIcpp+Rc7gnUaPbz6blJez9rARfQDB52+OkVa9DZ7djSElBWCwInQ4MBoROh4xGCezfR3V8jN/dGkdTxEMoFmozVpvBdlwvrkRLIsnWZFJtqaTZ0g4nvn7r2SUlVG+HA6ug5DPt5N9Yoi0zWCBrtlaXnz1Pex13/BhMSt8SQqyTUrZ6haMSgdI9kaBWtC9eA8WfaT9+f8skLvZUGLMIxpyqPZztX/n3h5iM4Ql7aAxofdAbgg1H+qQHG7Ur77CH/Z+8w4WfRzkxdR4iJiGqNYYKqxWd1YouzoExPQNDRgalD/0Xc7iRsUveQGc9plMwUY+HveddQCQtj6wbLydSVUW4qppIjRsZCkMkgozFIBohooNPawt559Q4xs89i3hzPHHGOOJMRx4us+vwCd9qsLZ2iP1LSq23zv6V2uPAJ+BtGQMsLhNyT4CcE7STf9rUNqsTlf7TXiJQVUNK9xjM2o899wTtfSym1fcWr9GuCvd9BFte0pYl5R9JCnkng9XV7+HqhI54UzzxpnhyaLsL58pRZ3FHxh38bdGlnJt3brv7NKanUXzDjdT885+k/vjHxyyrvusuorU1jP7ng1inTWt3P89uf5a/frGWFy58hMlJkzt/UP3NXw97lmuP/SuhqVT7PC4Dxp4Oo0+BvIWQMHA9npTuUYlA6R06ndajI3UCzLmppapgG+xdoSWFjc/C2v9oDdC582H8uVBwntYYOIicnHUyWY4sXtr5UoeJwD5vHq4rLqf2iSexn3QS9pZG3sa33qbhpZdJuvWWDpOAlJKXd73M5KTJgy8JSKl1INi9THuUfK71MrMmaif90T+E0Yu0Lpuqbn9IU1VDSv+IhLS2hb3LYdd7UNUycFbi2JakcC7kLhgU/b8f3fIo96+/nzcueYMxzjHtrhv1eDh49TWEKypI+/WvkH4/VX/6M5bp08h9/PEOxwraWL2R69+5njtPvJPLx1/em4fRPdGwdu/Jjrdg9/tHrvozpkP+2ZB/jta5YBjcVT7SqDYCZfBpKNYSwq53tWqGaAhsSTDxYq0nSd7JA3ayqfHXcNZLZ3H1xKv56dyfdrh+uKKC0m9/h0BREQC2+fPJvu/v6F2uDrf99apf8/7B9/nwyg+xG+09Db17IiGt1LbtddjxNgQatHFzxpwK48+BcWdBfPtdgpXBTyUCZXALerSSwrY3tPsYwl6tT/mhpDDqpH4fxOtHH/2Izys/54PLP8BisHS4voxE8BUWIgwGrLNmab2COtAcaub0/53OhWMv5M4T7+yNsDsvGtFO/ltfhh1LIdgIZidMOB8mLYYxp4Gx4+NWhg7VWKwMbmaHdvKZtFi7yW33Mih6DTY+B4WPQXw2zLgaZlyjdU3tB1cWXMmyg8t4/+D7XDT2og7XFwYD9vnzu/QdS/YuIRANcMX4Phtw93hVRdqQI5v/pw1kaHHBxItaTv6nqt49I5QqESiDV8gLO9/REsLeDwGplQ5mXKOVFEx9V5UipeTi1y/GZXbx9PlP98n+L3njEuxGO89d8Fyv7/8YgSbY/CKsf0obzkFn0Or6Z1ytPauT/4igSgTK0GSyw9TLtUdjGWx+ATY8C2/cAe/+EmZeB/Nu6ZNSghCCy8dfzt2Fd7OrfhfjE1qfyrK7CqsK2de4j/876f96db/HqCqCtY/Cphe16raM6XDeX2HK5WBP6rvvVYacwT0+sKIc4syChT+C76yDm97VBtb74mF4YBY8eyXs+UC7l6EXXTz2Ykw6Ey/tfKlX9wvw4s4XiTfFd9hFtctiUSh6HR4/F/69QCtNTb4Ebv0QvrESTviGSgLKcVSJQBlahIBRJ2qPpgpY9wQUPgHPXAYpE7VkMfnSXumGmmBJ4Oy8s3lr31v8YPYPujyQXVtq/DUsP7icqyde3amG6E6JBLW6/0/v1+74TciDs/5PKzXZ2h/nSFFUiUAZuuIz4LRfwg+2wqWPaJ+9egv8Y442SF6k7fF6OuuK8VfgCXt498C7Pd7XIa/ufpWIjHDl+Ct7vrOwH1b/A+6bBm9+D8zxcMVT8J31cNJ3VRJQOkUlAmXoM5hh+lfhW6vhq89os7Mt+Y6WELa+ot0h200zU2cyzjWOF3e+2CuhRmNRXt71MidknECeM68HO4rAuqfgwdmw7FeQMh5ueANu+0irClI3fCldoBKBMnzodFpXyNs+gmtfBnMcvHwzPHoGHFzdrV0KIbiy4Eq21W5jQ/WGHoe4vHg5Fd4Kriq4qns7kFIbFvxf8+HN72pz7974Ftz4ptb9Uw31oHSDSgTK8CME5J+lNY4u/pfWlvDEefDi9doY+V20eOxiXGYXj299vEdhSSl5YusTjIofxWk5p3V9B7V7tbaQF6/Trviveg6+/j6MXtijuBRFJQJl+NLpYea1Wk+j036t3aj2j3nwxX+61MPIZrRx9YSr+ajkI/Y17Ot2OIVVhWyt3cqNk2/s2jwB0TB89BetFFDyBZx7F3zzU5hwgSoBKL1CJQJl+DPZYNFP4PY1kD0blv4YnjwfGko6vYurJlyFRW/hyaInux3GY1sfI9GSyMVjL+78Ru5d8OiZ8NGftSE3vlMI8785KAbnU4YPlQiUkSNxDFz/Olzyb+0O24dOhu1vdW5TSyKXjLuEN/e9Sbmn69VL22q38WnZp1w78drOzf0rpVZyefgUbYC+K5+Gyx+DuPQuf7eidEQlAmVkEUIbouIbK7W+9i9eC+/8TOuF04GvT/06OnT8c+M/u/y196+/H5fZxdUTru545bAfXvuGVnLJOxlu/wwmdaEUoShdpBKBMjIljdUaWk/4Fnz+EDx3pTYmTzvS7elcM/Ea3tz7Jrvrd3f6qz6v+JzV5au5ZeotxJni2l/ZU63dFbz5f1q7xjX/U/P7Kn2u3xKBEOJxIUS1EGLrUZ8lCiHeF0LsbnlO6K94FAWDCc77C1x0P+z/WOtZ5K1pd5OvT/k6dqOdBzc82KmvkFJy37r7SLenc9WEDrqM1h+Ax88B906tR9Cin/T78NvKyNSf/8ueBL48sMrPgeVSynxgect7Relfs7+mXXnX7oGnLgKPu81VXRYXN0+5mRUlK1hTvqbDXb+570221m7l9um3t982UH9AKwn46uDGJdq8AIrST/otEUgpVwJ1X/p4MfBUy+ungEv6Kx5FOca4M+CaF6FuP/x3cbvVRDdMvoFR8aP4w5o/4I/421yvPlDP39b+jekp01k8bnHb391cBf+9BCIBuGkp5MzrwYEoStcNdLkzTUpZAdDynNrWikKI24QQhUKIQre77Ss2Rem2MafC1c9DzU7tjuQ2GpDNejO/nf9bSj2lPLzp4TZ3d3fh3XhCHu488U50oo2fWsir3STmqdbuhk4bZBPYKyPCQCeCTpNSPiKlnCOlnJOSkjLQ4SjD1djT4Py7Yc/7sOzXba42L2Mel4y7hKeKnmKze/Nxyz8u+Zgle5dw05SbyE/Ib30nUsJbP4SqrXDlfyG71TlDFKXPDXQiqBJCZAC0PFcPcDyKAnNuaulN9G/Y9V6bq/14zo9Js6fx449/TI3/SCNzcVMxv1z1SyYkTuAb07/R9vesf0qbbOfUX2jzKyjKABnoRLAEuLHl9Y3AGwMYi6IccdbvIW2KNoqp78tNWxqn2ck9i+6hIdjArctuZXf9bra4t3DLslvQCR33LLqn7Qbi2r2w9KfaJPGn/LgPD0RROtaf3UefB9YABUKIUiHE14G/AGcJIXYDZ7W8V5SBZzDDpQ9pSeDtH7W52uTkyTxw+gNU+ar4ypKvcM3SawjHwjx81sPkxue2vpGU2s1ih75DDRmtDLB+G7BEStnWLZVn9FcMitIl6VNh0U9hxR9hzs1tjvI5P2M+b17yJu/sfwe9Ts+FYy5s/8axba/D3g+1+YPVkBHKICBkDybtGChz5syRhYWFAx2GMhKE/fDgHHCkwC0f9vwGr5BXm0zGngK3rlCDxyn9RgixTkrZao+EgW4jUJTBzWiF038N5Rtg22s939/ax6C5As7/m0oCyqChEoGidGTalVrD8Qe/1+YG6K6QD1Y/oDUQ587vvfgUpYdUIlCUjuj0Wqmg4SAU9aBUsO4J8LrhVDWSijK4qESgKJ2Rfw4kF2hX9N1pVwv74dP7YfQpqjSgDDoqEShKZ+h0sODb2oQ2+z/u+vZbXwVPFSxsuyuqogwUlQgUpbOmXgn2VFjduSGoj7H2UUiZAKMX9X5citJDKhEoSmcZLTDvNtjzgTZnQGeVrYPy9TD3FjXZvDIoqUSgKF0x+2ugM0Lh453fZu1jYLTDtK/2WViK0hMqEShKVzhSYNJi2PicdnNYR3x1sPUVmP5VsMT3fXyK0g0qEShKV829BYJNsOXljtfd+Kw24czcW/o+LkXpJpUIFKWrcudD6iStAbi9rqSxmFYtlLtATTijDGoqEShKVwkBc78OlZu1huC27F0O9fu1dRVlEFOJQFG6Y9pXweTQrvjb8sV/tO6mEy/uv7gUpRtUIlCU7jDHaclg6yutT1xTfwB2L9N6GRlM/R2donSJSgSK0l1zvw7RoNYg/GVrHwOh06a9VJRBTiUCRemutMmQe6J20o/Fjnwe9sOGp2HCBRCfOXDxKUonqUSgKD0x9xatQXj7kiOfrf8v+OvhhHYmrleUQUQlAkXpicmXaqOSfvj/tLkKgs2w8m+QtxBGnTTQ0SlKp6gpkhSlJ3R6OPN38MLV8MYdWsOxr1b7TI0rpAwRKhEoSk9NOB9O+Sms/Kv2/oJ7ILvVqWEVZVBSiUBResPpv4LpV2klhIS8gY5GUbpEJQJF6S1JYwc6AkXpFtVYrCiKMsKpRKAoijLCqUSgKIoywg2KRCCEOFcIsVMIsUcI8fOBjkdRFGUkGfBEIITQA/8EzgMmAVcLISYNbFSKoigjx2DoNTQP2COl3AcghHgBWAxs6+0v+v2bRWwrb+rt3SqKovSLSZnx3HlR709yNOAlAiALKDnqfWnLZ8cQQtwmhCgUQhS63e5+C05RFGW4Gwwlgtbuwz9u/j8p5SPAIwBz5sxpZ37AtvVFJlUURRnqBkOJoBTIOep9NlA+QLEoiqKMOIMhEawF8oUQo4UQJuAqYEkH2yiKoii9ZMCrhqSUESHEt4H3AD3wuJSyaIDDUhRFGTEGPBEASCmXAksHOg5FUZSRaDBUDSmKoigDSCUCRVGUEU4lAkVRlBFOJQJFUZQRTkjZrXuzBpQQwg0c7ObmyUBNL4YzFKhjHv5G2vGCOuauGiWlTGltwZBMBD0hhCiUUo6oCWXVMQ9/I+14QR1zb1JVQ4qiKCOcSgSKoigj3EhMBI8MdAADQB3z8DfSjhfUMfeaEddGoCiKohxrJJYIFEVRlKOoRKAoijLCDctEIIQ4VwixUwixRwjx81aWCyHEAy3LNwshZg1EnL2pE8d8bcuxbhZCrBZCTB+IOHtTR8d81HpzhRBRIcTl/RlfX+jMMQshThVCbBRCFAkhPu7vGHtbJ/5vO4UQbwohNrUc800DEWdvEUI8LoSoFkJsbWN575+/pJTD6oE2lPVeYAxgAjYBk760zvnAO2izo80HPh/ouPvhmBcACS2vzxsJx3zUeh+ijW57+UDH3Q9/ZxfafN+5Le9TBzrufjjmXwJ3tbxOAeoA00DH3oNjPgWYBWxtY3mvn7+GY4lgHrBHSrlPShkCXgAWf2mdxcB/peYzwCWEyOjvQHtRh8cspVwtpaxvefsZ2kxwQ1ln/s4A3wFeAar7M7g+0pljvgZ4VUpZDCClHOrH3ZljlkCcEEIADrREEOnfMHuPlHIl2jG0pdfPX8MxEWQBJUe9L235rKvrDCVdPZ6vo11RDGUdHrMQIgu4FHioH+PqS535O48HEoQQHwkh1gkhbui36PpGZ475H8BEtClutwDfk1LG+ie8AdHr569BMTFNLxOtfPblPrKdWWco6fTxCCFOQ0sEJ/dpRH2vM8d8H/AzKWVUu1gc8jpzzAZgNnAGYAXWCCE+k1Lu6uvg+khnjvkcYCNwOjAWeF8I8YmUsqmPYxsovX7+Go6JoBTIOep9NtqVQlfXGUo6dTxCiGnAo8B5Usrafoqtr3TmmOcAL7QkgWTgfCFEREr5er9E2Ps6+3+7RkrpBbxCiJXAdGCoJoLOHPNNwF+kVoG+RwixH5gAfNE/Ifa7Xj9/DceqobVAvhBitBDCBFwFLPnSOkuAG1pa3+cDjVLKiv4OtBd1eMxCiFzgVeD6IXx1eLQOj1lKOVpKmSelzANeBm4fwkkAOvd/+w1goRDCIISwAScA2/s5zt7UmWMuRisBIYRIAwqAff0aZf/q9fPXsCsRSCkjQohvA++h9Th4XEpZJIT4Zsvyh9B6kJwP7AF8aFcUQ1Ynj/m3QBLwr5Yr5IgcwiM3dvKYh5XOHLOUcrsQ4l1gMxADHpVSttoNcSjo5N/5/4AnhRBb0KpNfialHLLDUwshngdOBZKFEKXAnYAR+u78pYaYUBRFGeGGY9WQoiiK0gUqESiKooxwKhEoiqKMcCoRKIqijHAqESiKooxwKhEoiqKMcCoRKEorhBCPCSFu74X9PCyEOKm396sovUndR6AoRxFCLAb+DsQBYaAM+JqUsqib+9uINvbPhb25X0XpTcPuzmJF6S4hxFjg32h3dX4XKASagJeFEFOklNEu7m8i2hg/eb25X0XpbapqSFGOOBt4/eixmKSUr6IN1ZDfmR0IIRKOense8G5v7FdR+pJKBIpyrNbqSgWtD/3bmr8f9foctETQ7n6FEHYhxFNCiP8IIa7tfKiK0jtUIlCUIz4ALm2pIgIOtxkYgH1CiD8JIe4TQvyrZVmKEOIJIUR2yzyzFwEThBA/bhn50yWlLO9gv7uArwAvSylvBS7ur4NVlENUG4GitJBS7m4Z6XI52qQui9GmDLwMuKXlswa0+XORUrqFEMXAPWiT/UwAnpFS/kMIcQGwoqP9tkyak402sxaAai9Q+p1KBIpyFCnlq0KI14DHgfVSygcBhBDfA+6QUgYPrSuEcKAlhYiU0tMy8c+mlsXnoc2B0O5+W5SiTS6yEVVKVwaA6j6qKJ3QUu1zDdpcsR+iVfc8AvweuBJtApV4tNLDX4BngROklOFO7NuONu9uAFglpXy2L45BUdqiEoGiKMoIp4qhiqIoI5xKBIqiKCOcSgSKoigjnEoEiqIoI5xKBIqiKCOcSgSKoigjnEoEiqIoI5xKBIqiKCOcSgSKoigj3P8Hq8P4Ux4adzgAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Display eigenvalue spectrum\n",
    "plt.figure()\n",
    "for i in range(n_eig):\n",
    "    plt.plot(phi_ext, (spec[i, :] - spec[0, :]))\n",
    "\n",
    "plt.xlabel(r\"$\\Phi_{ext}/\\Phi_0$\")\n",
    "plt.ylabel(r\" $\\omega_n / 2\\pi$  (GHz)\")\n",
    "plt.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}