{ "cells": [ { "cell_type": "markdown", "id": "c96b5dc6", "metadata": {}, "source": [ "# Dissipative spin bath\n", "In this tutorial we will go over the steps needed to simulate the decoherence of a central spin coupled to a dissipative, interacting spin baths governed by Lindblad Master equation using CCE method (ME-CCE) within the **PyCCE** module. Two methods of interest include:\n", "\n", "* Master Equation CCE (ME-CCE).\n", "* Master Equation gCCE (ME-gCCE).\n", "\n", "$\\newcommand{\\bra}[1]{\\left\\langle {#1} \\right|}$\n", "$\\newcommand{\\ket}[1]{\\left| {#1} \\right\\rangle}$\n", "$\\newcommand{\\Tr}{\\mathrm{Tr}}$" ] }, { "cell_type": "markdown", "id": "bf087276", "metadata": {}, "source": [ "The Lindblad master equation of the overall system can be written as:\n", "\n", "\\begin{equation}\n", "\\frac{d}{dt} \\hat \\rho=-\\frac{i}{\\hbar}[\\hat H, \\hat \\rho]+\\sum_i{\\gamma_i (\\hat L_i \\hat \\rho \\hat L_i^\\dagger - \\frac{1}{2}\\{\\hat L_i^\\dagger \\hat L_i, \\hat \\rho \\})},\n", "\\end{equation}\n", "\n", "where $\\hat \\rho$ is the density matrix of the system, $\\hat H$ is the Hamiltonian, and $\\hat L_i$ are jump operators with corresponding dissipation rates $\\gamma_i$.\n", "\n", "Within the conventional CCE framework, the coherence of the central spin is recovered from the trace of the partial inner product $\\hat \\rho_{01} (t) =\\bra{0}\\hat \\rho (t)\\ket{1}$ as $\\mathcal{L}(t)=\\Tr [\\hat \\rho_{01} (t)]/\\Tr [\\hat \\rho_{01} (0)]$. The evolution of $\\hat \\rho_{01}$ by solving the following:\n", "\n", "\\begin{equation}\n", "\\frac{d}{dt} \\hat \\rho_{01} (t) = \\mathfrak{I} \\cdot \\hat \\rho_{01} (t)= -\\frac{i}{\\hbar} \\hat H^{(0)} \\hat \\rho_{01}(t) + \\frac{i}{\\hbar} \\hat \\rho_{01}(t) \\hat H^{(1)} + \\sum_i{\\gamma_i (\\hat L_i \\hat \\rho \\hat L_i^\\dagger - \\frac{1}{2}\\{\\hat L_i^\\dagger \\hat L_i, \\hat \\rho \\})},\n", "\\end{equation}\n", "\n", "Within the generalized CCE framework, one needs to solve full Lindbladian for each cluster including central spin." ] }, { "cell_type": "code", "execution_count": null, "id": "00c24e3e", "metadata": { "tags": [] }, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "import pandas as pd\n", "import pycce as pc" ] }, { "cell_type": "markdown", "id": "02b520b1", "metadata": {}, "source": [ "## Setup the simulator properties\n", "\n", "As an example we consider the dissipative electron spin bath. First, we generate electron spin bath with concentration $\\rho=10^{16}$ cm$^{-3}$ " ] }, { "cell_type": "code", "execution_count": null, "id": "821cf053", "metadata": { "tags": [] }, "outputs": [], "source": [ "electrons = pc.random_bath('e', [5e3, 5e3, 5e3], density=1e16,\n", " density_units='cm-3', seed=2) # Density in 1/cm^3" ] }, { "cell_type": "markdown", "id": "b8ead82a", "metadata": {}, "source": [ "The properties of the central spin-1/2 are stored in the `CenterArray`." ] }, { "cell_type": "code", "execution_count": 3, "id": "18ec8abe", "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "CenterArray\n", "(s: [0.5],\n", "xyz:\n", "[[0. 0. 0.]],\n", "zfs:\n", "[[[0. 0. 0.]\n", " [0. 0. 0.]\n", " [0. 0. 0.]]],\n", "gyro:\n", "[[[-17608.59705 -0. -0. ]\n", " [ -0. -17608.59705 -0. ]\n", " [ -0. -0. -17608.59705]]])\n" ] } ], "source": [ "center = pc.CenterArray(spin=1/2, alpha=[1,0], beta=[0,1])\n", "print(center)" ] }, { "cell_type": "markdown", "id": "014e8d43", "metadata": {}, "source": [ " To run the simulations we need to setup the `Simulator` object." ] }, { "cell_type": "code", "execution_count": 4, "id": "ab333eb1", "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Simulator for center array of size 1.\n", "magnetic field:\n", "array([ 0., 0., 300.])\n", "\n", "Parameters of cluster expansion:\n", "r_bath: 1400.0\n", "r_dipole: 700.0\n", "order: 2\n", "\n", "Bath consists of 104 spins.\n", "\n", "Clusters include:\n", "104 clusters of order 1.\n", "539 clusters of order 2.\n", "\n" ] } ], "source": [ "calc = pc.Simulator(center, bath=electrons, order=2, r_bath=1.4e3, r_dipole=0.7e3, pulses=1, magnetic_field=300, n_clusters=None)\n", "print(calc)" ] }, { "cell_type": "markdown", "id": "91273f59", "metadata": {}, "source": [ "As a reference we compute the coherence with conventional CCE assuming no dissipation:" ] }, { "cell_type": "code", "execution_count": 5, "id": "64528c63", "metadata": { "tags": [] }, "outputs": [ { "data": { "text/plain": [ "Text(0.5, 0, 'Time (ms)')" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiMAAAGwCAYAAAB7MGXBAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABAQElEQVR4nO3deXhU9cH+//vMTCYLZCEkBEISCPu+hTWIFsUgKsrjhtWyiVasS5HWVupTt+f3LV2spS5gKYu1RcSq1A2XtAphVyCALLIGEkhCCEsSkpBtzu+PQGpMgExIcmYm79d1TWFOPjPn/nQ4V27PnMUwTdMUAACARWxWBwAAAM0bZQQAAFiKMgIAACxFGQEAAJaijAAAAEtRRgAAgKUoIwAAwFIOqwPUhcvlUmZmpoKDg2UYhtVxAABAHZimqYKCAkVHR8tmu/j+D68oI5mZmYqNjbU6BgAAqIeMjAzFxMRc9OdeUUaCg4MlVU4mJCTE4jQAAKAu8vPzFRsbW/V7/GK8ooxc+GomJCSEMgIAgJe53CEWHMAKAAAsRRkBAACWoowAAABLUUYAAIClKCMAAMBSlBEAAGApyggAALAUZQQAAFiKMgIAACxFGQEAAJZyu4ykpKRo/Pjxio6OlmEY+te//nXZ16xevVoJCQkKCAhQp06d9Nprr9UnKwAA8EFul5HCwkL1799fr7zySp3Gp6Wl6cYbb9SoUaOUmpqqX/3qV3rsscf07rvvuh0WAAD4HrdvlDdu3DiNGzeuzuNfe+01xcXFae7cuZKknj17avPmzXrhhRd0++23u7v6BnXybImKSitkGJU38TGkyr/LOP+nJEOyG4aCnA4F+Nkue7MfAADgnka/a++GDRuUlJRUbdnYsWO1aNEilZWVyc/Pr8ZrSkpKVFJSUvU8Pz+/UbI9++Fufbg9s87jDUNq4XQoyGlXC//zfzodCvK3q6W/QzGtgtShdZA6hAepQ0QLtQ0JkN1GeQEA4FIavYxkZ2crKiqq2rKoqCiVl5crNzdX7dq1q/GaOXPm6LnnnmvsaPKzG/J32GRKkimZMmWakinJNM3zf/53vGlKZ0vKdbakXCooqf1Nv8NptykmPLCynLRuoY6tg9Q3JlS9o0MV4GdvpFkBAOBdGr2MSKrx1YZ5/jf8xb7ymD17tmbNmlX1PD8/X7GxsQ2e68W7BujFuwZcdlyFy1RxWYWKSspVVFqhwtLzf5b898+84jJlnCrSkVNFOnKySEdPF6m0wqVDJwp16EShpBNV7+ewGerRLlj9Y8I0ILby0TmypWzsRQEANEONXkbatm2r7OzsastycnLkcDjUunXrWl/j7+8vf3//xo5WZ3aboZb+DrX0r/v/XRUuU5lnipV+qkiHTxYq/WSRDp44q20Zeco9W6Kdx/K181i+lm5KlyS19HeoX0yoBsSGqd/5ktI2NKCxpgQAgMdo9DIyYsQIffjhh9WWff755xo8eHCtx4v4CrvNUGx4kGLDgzSyS0TVctM0lZl3TtvSz2j70TPaln5G3xzL09mScq0/eFLrD56sGtsm2F/9YsLUPyZU/WPD1C8mVGFBTiumAwBAo3G7jJw9e1YHDhyoep6WlqZt27YpPDxccXFxmj17to4dO6Y33nhDkjRjxgy98sormjVrlh544AFt2LBBixYt0rJlyxpuFl7EMAy1DwtU+7BA3dSv8niZ8gqX9uec1baMM1UlZX/OWeUUlOjfe47r33uOV72+Q+sg9YsJU4+2werSpqW6RQUrLjyIA2UBAF7LMM3vHqJ5eatWrdLo0aNrLJ8yZYpef/11TZ06VYcPH9aqVauqfrZ69Wo9/vjj2rVrl6Kjo/XLX/5SM2bMqPM68/PzFRoaqry8PIWEhLgT12sVl1ZoV2aeth/N0/aMM9px9IwOnyyqdazTYVOniBbqGhWsrm1aVj6ighUf0YKSAgCwTF1/f7tdRqzQHMtIbc4UlWrH0TztzMzT/uNntT+nQAdyzupcmavW8YF+dvVsF6ze0aHq0z5EvaND1S0qWE4HdwEAADQ+ykgz4XKZOnamWPtzCrTv+FntP35WB87/vbisosZ4P7uhblHB6hMdqj4xoRrTs43ahQZakBwA4OsoI81chctUWm6hdmXmaVdmvnYey9POY3nKP1debZxhSCM7R+i2Qe01tndbtXDjjCEAAC6FMoIaTNPU0dPF2pWZp53H8rUp7aS+Pny66udBTrtu6N1Wtw2K0YjOrTneBABwRSgjqJOMU0VakXpM7209Wu0A2bYhAZowsL1uH9ReXaOCLUwIAPBWlBG4xTRNbU0/o/e2HtWH2zOrfZ1zU992euaWXmoTzEXYAAB1RxlBvZWUV+iLPTl6d+sxfbk3RxUuUyEBDs2+sacmDo7lsvUAgDqhjKBB7DyWp9nvfaNvjuVJkobGh2vObX3VObKlxckAAJ6urr+/ueAELqlP+1Ct+Emi/vemngr0s+urtFMaN3eNXvrPfpWW1359EwAA3EEZwWU57DbdP6qTPn/8al3TLVKlFS69mLxPN7+8RluOnL78GwAAcAmUEdRZbHiQXp82RH++e4Bat3Bq3/GzuuO19Xr6/Z06W1J++TcAAKAWlBG4xTAM3Tqgvf496xrdkRAj05Te2HBEE/+yQScKSqyOBwDwQpQR1EurFk69cGd/Lb1/mFq3cGpXZr7ueG290i9yMz8AAC6GMoIrMrJLhN55KFExrQJ15GSRbpu/Xrsy86yOBQDwIpQRXLH4iBZ676FE9WgbrNyzJbr7Lxu14eBJq2MBALwEZQQNok1IgJY/OEJD48NVUFKuKUu+0qc7s6yOBQDwApQRNJjQQD+9cd9QJfWKUmm5Sz9ZulVvbkq3OhYAwMNRRtCgAvzsmnfvIN09JFYuU/rVim/00n/2ywsu9AsAsAhlBA3OYbdpzm199cjoLpKkF5P36ZkPdsnlopAAAGqijKBRGIahn4/trmfH95JUeS2Sp/61kz0kAIAaKCNoVFNHxuvPdw+QYUjLvkrXgpRDVkcCAHgYygga3a0D2ut/b6rcQ/LbT7/lLBsAQDWUETSJ+0Z21KThHWSa0szl27Q944zVkQAAHoIygiZhGIaeGd9LP+geqXNlLk3/22YdPc2l4wEAlBE0IYfdppd/OLDqSq3TX9+s/HNlVscCAFiMMoImFRzgp8VTh6hNsL/2Hi/QI2+mqrzCZXUsAICFKCNoctFhgVo0ZYgC/exK2XdCz3ywi1N+AaAZo4zAEn1jQqtO+V26KV2L1qZZHQkAYBHKCCyT1LutnrqxpyTp/63co892ZVucCABgBcoILDX9qnjdOyyu8pTft7bpm6N5VkcCADQxyggsZRiGnrult67uFqnisgo9umyrikrLrY4FAGhClBFY7sIpv+1CA3T4ZJHmrPzW6kgAgCZEGYFHCA300x/u6C9J+vvGI1q974TFiQAATYUyAo9xVdcITU3sKEn6xTvbdaao1NpAAIAmQRmBR/nlDT3UKbKFjueX6Nfv77I6DgCgCVBG4FECnXb96a4BstsMfbg9Ux9sz7Q6EgCgkVFG4HH6x4bpkdFdJEm//tdOHc8/Z3EiAEBjoozAIz1ybRf1bR+qvOIyPfHODi4XDwA+jDICj+Rnt+lPE/vL32FTyr4TWrop3epIAIBGQhmBx+rSJli/vKGHJOn/fbxHabmFFicCADQGygg82tTEjhrRqbWKyyr0s7e3qbzCZXUkAEADo4zAo9lshl64q7+C/R3amn5Gf0k5ZHUkAEADo4zA47UPC9Szt/SWJM399z7tPMbN9ADAl1BG4BVuG9ReY3tHqazC1FP/2imXi7NrAMBXUEbgFQzD0P9N6KMWTru2Z5zRhzu4GBoA+ArKCLxGm+AA/eT8xdB+98m3Ki6tsDgRAKAhUEbgVaZfFa/2YYHKzDunhWs4mBUAfAFlBF4lwM+uX9zQXZI0f/VBLhUPAD6AMgKvc0v/aA2MC1NRaYVe+Gyv1XEAAFeIMgKvYxiGfn1zL0nSO1uPcqovAHg5ygi80qC4Vrqlf7RMU/r/Pt7NjfQAwItRRuC1fjmuh/wdNm08dEqf7z5udRwAQD1RRuC12ocF6v5R8ZKkOSv3qLSc+9YAgDeijMCrPfSDLooM9tfhk0V6Y8Nhq+MAAOqBMgKv1tLfoZ8ndZMk/fk/+3WqsNTiRAAAd1FG4PXuSIhVz3YhKjhXrj//e5/VcQAAbqKMwOvZbYZ+fVNPSdI/NqXrQE6BxYkAAO6gjMAnJHaJ0JieUapwmfp/H++xOg4AwA31KiPz5s1TfHy8AgIClJCQoDVr1lxy/NKlS9W/f38FBQWpXbt2mjZtmk6ePFmvwMDF/OrGHnLYDH2594RS9p2wOg4AoI7cLiPLly/XzJkz9dRTTyk1NVWjRo3SuHHjlJ6eXuv4tWvXavLkyZo+fbp27dqlf/7zn/r66691//33X3F44Ls6RbbU5BEdJUl/+GwvF0IDAC/hdhl58cUXNX36dN1///3q2bOn5s6dq9jYWM2fP7/W8Rs3blTHjh312GOPKT4+XldddZUefPBBbd68+aLrKCkpUX5+frUHUBcPj+6sQD+7vjmWp1XsHQEAr+BWGSktLdWWLVuUlJRUbXlSUpLWr19f62sSExN19OhRrVy5UqZp6vjx43rnnXd00003XXQ9c+bMUWhoaNUjNjbWnZhoxlq39NePhsdJkl7+z372jgCAF3CrjOTm5qqiokJRUVHVlkdFRSk7O7vW1yQmJmrp0qWaOHGinE6n2rZtq7CwML388ssXXc/s2bOVl5dX9cjIyHAnJpq5B67uJH+HTVvTz2j9QY5NAgBPV68DWA3DqPbcNM0ayy7YvXu3HnvsMT399NPasmWLPv30U6WlpWnGjBkXfX9/f3+FhIRUewB11SY4QD8cWrl35KX/7Lc4DQDgctwqIxEREbLb7TX2guTk5NTYW3LBnDlzNHLkSD3xxBPq16+fxo4dq3nz5mnx4sXKysqqf3LgEh68ppP87IY2pZ3SV2mnrI4DALgEt8qI0+lUQkKCkpOTqy1PTk5WYmJira8pKiqSzVZ9NXa7XZL4Ph+Npl1ooO5IqDzW6OUv2DsCAJ7M7a9pZs2apYULF2rx4sXas2ePHn/8caWnp1d97TJ79mxNnjy5avz48eP13nvvaf78+Tp06JDWrVunxx57TEOHDlV0dHTDzQT4np/8oLPsNkNr9ucqNf201XEAABfhcPcFEydO1MmTJ/X8888rKytLffr00cqVK9WhQwdJUlZWVrVrjkydOlUFBQV65ZVX9LOf/UxhYWG69tpr9bvf/a7hZgHUIjY8SP8zsL3e2XJUL39xQIunDrE6EgCgFobpBd+V5OfnKzQ0VHl5eRzMCrek5Rbquj+uksuUPnr0KvVpH2p1JABoNur6+5t708CnxUe00Pj+lV8HvvLFAYvTAABqQxmBz3tkdBcZhvTprmztzeaOvgDgaSgj8Hldo4I1rk9bSdIrX7J3BAA8DWUEzcIjo7tKkj7akamDJ85anAYA8F2UETQLvaJDNKZnlExTepW9IwDgUSgjaDYeu66LJOn9bZlKP1lkcRoAwAWUETQb/WLCdE23SFW4TM1bxd4RAPAUlBE0K49eW7l35N2tR3XsTLHFaQAAEmUEzczgjuEa0am1yipM/WX1QavjAABEGUEzdGHvyNubM3S6sNTiNAAAygianRGdW6tXuxCdK3Np6aYjVscBgGaPMoJmxzAM/fjqTpKkv204opLyCosTAUDzRhlBs3RTv3ZqFxqgEwUlen9bptVxAKBZo4ygWfKz2zQ1saMkaeGaQ/KCm1cDgM+ijKDZuntonFo47dp3/KxW7zthdRwAaLYoI2i2QgP9NHFInCRp4Zo0i9MAQPNFGUGzNm1kR9kMae2BXO3OzLc6DgA0S5QRNGux4UEa17edJGnh2kMWpwGA5okygmbvgVGVp/l+uD1Tx/PPWZwGAJofygiavQGxYRraMVxlFaZeX3/Y6jgA0OxQRgBJ94+KlyQt3XhEhSXlFqcBgOaFMgJIGtMzSvERLZR/rlxvb86wOg4ANCuUEUCSzWbovqsq944sXpemChcXQQOApkIZAc67Y1CMWgX5KeNUsT7blW11HABoNigjwHmBTrt+NLyDJOmvazjNFwCaCmUE+I5JIzrIabcpNf2Mthw5ZXUcAGgWKCPAd7QJDtCEgdGSpL+mcIl4AGgKlBHge+4/fxG0z3Zn68jJQovTAIDvo4wA39MtKlg/6B4p05QWrWXvCAA0NsoIUIsLl4h/Z8tR5RWXWZwGAHwbZQSoRWLn1uoeFayi0gr9k4ugAUCjoowAtTAMQ1NHdpQkvb7+MBdBA4BGRBkBLmLCgPYKC/LT0dPFSt593Oo4AOCzKCPARQQ67frh0DhJ0pJ1HMgKAI2FMgJcwuQRHWS3GdqUdkq7MvOsjgMAPokyAlxCu9BAjevTVpL0+rrD1oYBAB9FGQEuY9rIyrv5vr8tU7lnSyxOAwC+hzICXMaguDD1jwlVaYVLb25KtzoOAPgcyghwGYZh6L6rKveO/H3jEZWWuyxOBAC+hTIC1MG4Pu3UJthfJwpKtPKbLKvjAIBPoYwAdeB02DRpeAdJ0uJ1aTJNLoIGAA2FMgLU0T3D4uR02LTjaJ62pp+2Og4A+AzKCFBHrVv669b+0ZKkxZzmCwANhjICuOHCab6f7sxW5plii9MAgG+gjABu6BUdouGdwlXhMvX3jUesjgMAPoEyArjpwt6RZV+lq7i0wuI0AOD9KCOAm8b0jFJseKDOFJVpReoxq+MAgNejjABustsMTRnRUZL0+npO8wWAK0UZAerhriGxauG0a9/xs1p34KTVcQDAq1FGgHoICfDTHQkxkqQl69IsTgMA3o0yAtTTlMSOkqQv9uYo/WSRtWEAwItRRoB66hTZUld3i5RpSn/feNjqOADgtSgjwBWYmlh5v5rlX2eoqLTc4jQA4J0oI8AV+EG3NooLD1L+uXK9vy3T6jgA4JUoI8AVsNkMTR5RuXfkb+sPc5ovANQDZQS4QncmxCrQz65vswu0Ke2U1XEAwOvUq4zMmzdP8fHxCggIUEJCgtasWXPJ8SUlJXrqqafUoUMH+fv7q3Pnzlq8eHG9AgOeJjTITxMGtpckvbHhsLVhAMALuV1Gli9frpkzZ+qpp55SamqqRo0apXHjxik9Pf2ir7nrrrv0n//8R4sWLdLevXu1bNky9ejR44qCA55kyvkDWT/bdZy7+QKAmwzTzS+5hw0bpkGDBmn+/PlVy3r27KkJEyZozpw5NcZ/+umnuvvuu3Xo0CGFh4fXaR0lJSUqKSmpep6fn6/Y2Fjl5eUpJCTEnbhAk7l7wQZtPHRKD4/urCfGUrYBID8/X6GhoZf9/e3WnpHS0lJt2bJFSUlJ1ZYnJSVp/fr1tb7mgw8+0ODBg/X73/9e7du3V7du3fTzn/9cxcUX/6/HOXPmKDQ0tOoRGxvrTkzAEhfuV7PsqwydK+NuvgBQV26VkdzcXFVUVCgqKqra8qioKGVnZ9f6mkOHDmnt2rXauXOnVqxYoblz5+qdd97Rww8/fNH1zJ49W3l5eVWPjIwMd2IClri+V5SiQwN0qrBUH+/IsjoOAHiNeh3AahhGteemadZYdoHL5ZJhGFq6dKmGDh2qG2+8US+++KJef/31i+4d8ff3V0hISLUH4OkcdpvuHX7+NN8NnOYLAHXlVhmJiIiQ3W6vsRckJyenxt6SC9q1a6f27dsrNDS0alnPnj1lmqaOHj1aj8iA57p7SKycDpt2HM3TtowzVscBAK/gVhlxOp1KSEhQcnJyteXJyclKTEys9TUjR45UZmamzp49W7Vs3759stlsiomJqUdkwHO1bumv8f2iJVVeBA0AcHluf00za9YsLVy4UIsXL9aePXv0+OOPKz09XTNmzJBUebzH5MmTq8bfc889at26taZNm6bdu3crJSVFTzzxhO677z4FBgY23EwADzH1/N18P/4mSzkF56wNAwBewO0yMnHiRM2dO1fPP/+8BgwYoJSUFK1cuVIdOlR+V56VlVXtmiMtW7ZUcnKyzpw5o8GDB+vee+/V+PHj9dJLLzXcLAAP0jcmVAPjwlRWYeqtrzj4GgAux+3rjFihrucpA57i/W3H9NO3tikqxF9rf3mt/OzceQFA89Mo1xkBUDfj+rRTREt/Hc8v0ac7az/tHQBQiTICNAKnw6Z7hsVJ4n41AHA5lBGgkdw7LE4Om6GvD5/Wrsw8q+MAgMeijACNJCokQOP6tpPEab4AcCmUEaARTRlReZbZB9szlVdUZnEaAPBMlBGgESV0aKUebYN1rsyld7ZyxWEAqA1lBGhEhmFo0vm9I//YeEQul8efSQ8ATY4yAjSyCQPaq6W/Q2m5hVp/8KTVcQDA41BGgEbWwt+h2wa1lyT9feNha8MAgAeijABN4EfDK7+q+feeHGXlFVucBgA8C2UEaALdooI1LD5cFS5Ty7hfDQBUQxkBmsiFA1mXfZWusgqXxWkAwHNQRoAmktSrrSKD/XWioESf7zpudRwA8BiUEaCJOB023T0kVlLlab4AgEqUEaAJ/XBonGyGtOHQSR3IKbA6DgB4BMoI0ISiwwI1pmeUJOkfG9MtTgMAnoEyAjSxCweyvrvlqApLyi1OAwDWo4wATWxk5wjFR7RQQUm5PtieaXUcALAcZQRoYjaboXuHxUmS/r7hiEyT+9UAaN4oI4AF7kiIkb/Dpt1Z+dqafsbqOABgKcoIYIGwIKdu6R8tidN8AYAyAljkwoGsH+/I0smzJRanAQDrUEYAi/SLCVP/mFCVVrj0zy1HrY4DAJahjAAWuvf83XyXbjqiChcHsgJonigjgIXG94tWaKCfMk4VK2XfCavjAIAlKCOAhQKddt2ZECNJ+jsHsgJopigjgMUufFXz5d4cHT1dZHEaAGh6lBHAYvERLTSyS2uZpvTWVxlWxwGAJkcZATzAj4ZV7h156+sMlVW4LE4DAE2LMgJ4gDG9ohQZ7K/csyX6fNdxq+MAQJOijAAewM9u091DYiVVnuYLAM0JZQTwEHcPjZPNkNYfPKmDJ85aHQcAmgxlBPAQ7cMCNbp7G0nSsk3pFqcBgKZDGQE8yI/On+b7zy1Hda6swuI0ANA0KCOAB7m6W6TahwUqr7hMH+/IsjoOADQJygjgQew2Q/cMi5PEgawAmg/KCOBh7hwcI4fN0Nb0M9qdmW91HABodJQRwMO0CQ7Q2D5tJbF3BEDzQBkBPNC957+q+VfqMZ0tKbc4DQA0LsoI4IFGdGqtThEtVFhaofe3HbM6DgA0KsoI4IEM478Hsv5jY7pM07Q4EQA0HsoI4KHuSIiRv8OmPVn5Ss04Y3UcAGg0lBHAQ4UFOXVzv2hJ0tKNXJEVgO+ijAAe7N7hlV/VfLQjU2eKSi1OAwCNgzICeLCBsWHq2S5EJeUuvbuVA1kB+CbKCODBDMOoOs136aYjHMgKwCdRRgAPN2Fge7Vw2nXoRKE2HDppdRwAaHCUEcDDtfR3aMLA9pI4kBWAb6KMAF7g3mEdJEmf7cpWTv45i9MAQMOijABeoFd0iBI6tFK5y9RbX2dYHQcAGhRlBPASk4ZX7h15c1O6yitcFqcBgIZDGQG8xLi+bdW6hVPZ+ef07z05VscBgAZDGQG8hL/DrolDYiVJ/9h4xOI0ANBwKCOAF7lnWJwMQ1p7IFcHT5y1Og4ANAjKCOBFYloF6boebSSxdwSA76CMAF7mR+cPZH1ny1EVlZZbnAYArly9ysi8efMUHx+vgIAAJSQkaM2aNXV63bp16+RwODRgwID6rBaApKu7RiouPEgF58r1wbZMq+MAwBVzu4wsX75cM2fO1FNPPaXU1FSNGjVK48aNU3r6pa8MmZeXp8mTJ+u6666rd1gAks1m6Efn7+b7xgbuVwPA+7ldRl588UVNnz5d999/v3r27Km5c+cqNjZW8+fPv+TrHnzwQd1zzz0aMWJEvcMCqHRnQqycDpt2Z+UrNeOM1XEA4Iq4VUZKS0u1ZcsWJSUlVVuelJSk9evXX/R1S5Ys0cGDB/XMM8/UaT0lJSXKz8+v9gDwX61aODW+X7Qk6R8bOJAVgHdzq4zk5uaqoqJCUVFR1ZZHRUUpOzu71tfs379fTz75pJYuXSqHw1Gn9cyZM0ehoaFVj9jYWHdiAs3C5BGVB7J+tCNLpwpLLU4DAPVXrwNYDcOo9tw0zRrLJKmiokL33HOPnnvuOXXr1q3O7z979mzl5eVVPTIyuBcH8H39Y8PULyZUpRUuLed+NQC8mFtlJCIiQna7vcZekJycnBp7SySpoKBAmzdv1iOPPCKHwyGHw6Hnn39e27dvl8Ph0BdffFHrevz9/RUSElLtAaCmC6f5Lt10RBUuDmQF4J3cKiNOp1MJCQlKTk6utjw5OVmJiYk1xoeEhOibb77Rtm3bqh4zZsxQ9+7dtW3bNg0bNuzK0gPN3Ph+0QoN9NPR08VavY/71QDwTnU7iOM7Zs2apUmTJmnw4MEaMWKEFixYoPT0dM2YMUNS5Vcsx44d0xtvvCGbzaY+ffpUe32bNm0UEBBQYzkA9wU67bozIUYL16bp7xuO6NoeNfdQAoCnc7uMTJw4USdPntTzzz+vrKws9enTRytXrlSHDpW7i7Oysi57zREADefe4R20cG2aVu07ofSTRYprHWR1JABwi2F6wRWT8vPzFRoaqry8PI4fAWoxadEmrdmfqwev6aTZ43paHQcAJNX99zf3pgF8wOQRHSVJb3+doXNlFdaGAQA3UUYAH3BtjzZqHxao00VlWvlNltVxAMAtlBHAB9hthu4ZVnm/mr+tP8z9agB4FcoI4CMmDomV027T9qN52pp+xuo4AFBnlBHAR0S09NetAyrvV7N4XZrFaQCg7igjgA+ZNjJekvTpzmwdO1NscRoAqBvKCOBDekWHaESn1qpwmXpjw2Gr4wBAnVBGAB9z31WVe0eWbUpXUWm5xWkA4PIoI4CPua5HG3VoHaT8c+V6d+sxq+MAwGVRRgAfY7MZmpbYUZK0ZF2aXNzNF4CHo4wAPuiOwbEK9nfo0IlCrd5/wuo4AHBJlBHAB7X0d2jikFhJ0uK1nOYLwLNRRgAfNSWxo2yGtGZ/rvYdL7A6DgBcFGUE8FGx4UFK6tVWUuWxIwDgqSgjgA+7cJrve1uP6VRhqcVpAKB2lBHAhw3p2Ep92oeopNylZV+lWx0HAGpFGQF8mGEYmn5+78gbGw6rtNxlcSIAqIkyAvi4m/pGKzLYX8fzS/TJziyr4wBADZQRwMc5HTZNHt5BkrRobZpMk4ugAfAslBGgGbhnWJycDpt2HM3T1vTTVscBgGooI0Az0Lqlv/5nQHtJ0uK1h60NAwDfQxkBmolpV3WUJH2yM0tHTxdZGwYAvoMyAjQTPdqGaGSX1nKZ0hsbjlgdBwCqUEaAZuS+kZWn+S7blK684jKL0wBAJcoI0IyM7t5G3aJaqqCkXP/YyN4RAJ6BMgI0IzaboZ/8oIukytN8i0srLE4EAJQRoNm5uV87xYYH6lRhqZZ/zSXiAViPMgI0Mw67TQ9e3VmStCDlEJeIB2A5ygjQDN2REKPIYH9l5p3Tv7YdszoOgGaOMgI0QwF+dj0wqvLMmtdWH1SFi0vEA7AOZQRopu4Z1kGhgX46dKJQn+3KtjoOgGaMMgI0Uy39HZqS2FGS9OqXB7iBHgDLUEaAZmxaYkcFOe3alZmvlP25VscB0ExRRoBmrFULp344NE5S5d4RALACZQRo5h4Y1Ul+dkNfpZ3S5sOnrI4DoBmijADNXNvQAN0+KEaSNG/VQYvTAGiOKCMA9OA1nWUzpC++zdGuzDyr4wBoZigjABQf0UI39YuWJM1n7wiAJkYZASBJeuiaykvEr/wmS2m5hRanAdCcUEYASJJ6RYfo2h5t5DKlv6xm7wiApkMZAVDl4dGVe0fe3XpUWXnFFqcB0FxQRgBUSegQrqHx4SqrMPXXlDSr4wBoJigjAKp5eHQXSdLSTUeUnXfO4jQAmgPKCIBqru4aocEdWqmk3KWXvthvdRwAzQBlBEA1hmHoFzf0kCQt/zqDM2sANDrKCIAahsaH6wfdI1XhMvVi8j6r4wDwcZQRALV6Ymx3SdKH2zO5KiuARkUZAVCr3tGhGt+/8qqsL3y21+I0AHwZZQTARc26vpvsNkNf7j2hr7mjL4BGQhkBcFHxES101+BYSdLvP/1WpmlanAiAL6KMALikn17XVf4Om74+fFqr9p6wOg4AH0QZAXBJbUMDNCWxoyTp95/tlcvF3hEADYsyAuCyHrqms4L9HdqTla8Pd2RaHQeAj6GMALisVi2ceuDqTpKkF5P3qazCZXEiAL6EMgKgTqZfFa/WLZw6crJIb2/OsDoOAB9CGQFQJy38HXrk2sqb6L30n/06V1ZhcSIAvqJeZWTevHmKj49XQECAEhIStGbNmouOfe+993T99dcrMjJSISEhGjFihD777LN6BwZgnXuGxal9WKCO55fob+sPWx0HgI9wu4wsX75cM2fO1FNPPaXU1FSNGjVK48aNU3p6eq3jU1JSdP3112vlypXasmWLRo8erfHjxys1NfWKwwNoWv4Ou2aO6SpJmrfqoPKKyyxOBMAXGKabVzEaNmyYBg0apPnz51ct69mzpyZMmKA5c+bU6T169+6tiRMn6umnn67T+Pz8fIWGhiovL08hISHuxAXQwCpcpsbOTdGBnLN6ZHQX/fz8PWwA4Pvq+vvbrT0jpaWl2rJli5KSkqotT0pK0vr16+v0Hi6XSwUFBQoPD7/omJKSEuXn51d7APAMdpuhnyd1kyT9dc0hZZwqsjgRAG/nVhnJzc1VRUWFoqKiqi2PiopSdnZ2nd7jj3/8owoLC3XXXXdddMycOXMUGhpa9YiNjXUnJoBGNrZ3Ww3vFK6Scpee+3C31XEAeLl6HcBqGEa156Zp1lhWm2XLlunZZ5/V8uXL1aZNm4uOmz17tvLy8qoeGRmcRgh4EsMw9H+39pHDZujfe47rP3uOWx0JgBdzq4xERETIbrfX2AuSk5NTY2/J9y1fvlzTp0/X22+/rTFjxlxyrL+/v0JCQqo9AHiWrlHBmn5VvCTp2Q93caovgHpzq4w4nU4lJCQoOTm52vLk5GQlJiZe9HXLli3T1KlT9eabb+qmm26qX1IAHuex67qqbUiAMk4Va/6qg1bHAeCl3P6aZtasWVq4cKEWL16sPXv26PHHH1d6erpmzJghqfIrlsmTJ1eNX7ZsmSZPnqw//vGPGj58uLKzs5Wdna28vLyGmwUAS7Twd+jXN/eSJM1ffVBHThZanAiAN3K7jEycOFFz587V888/rwEDBiglJUUrV65Uhw4dJElZWVnVrjnyl7/8ReXl5Xr44YfVrl27qsdPf/rThpsFAMvc2LetruoSodJyl579YJfcvFoAALh/nRErcJ0RwLMdPHFWN8xNUVmFqQWTEpTUu63VkQB4gEa5zggA1KZzZEs9MKryrr7PfbhbxaUczAqg7igjABrEI9d2UfuwQB07U6xXvzxgdRwAXoQyAqBBBDn/ezDrgpRDOnTirMWJAHgLygiABjO2d5Su6Rap0gqXnuFgVgB1RBkB0GAMw9Bzt/SW027Tmv25+nRn3W4TAaB5o4wAaFAdI1poxjWVB7M+/9FuFZWWW5wIgKejjABocD8Z3UUxrQKVlXdOf/7PfqvjAPBwlBEADS7Az65nx/eWJP015ZA2Hz5lcSIAnowyAqBRjOkVpdsGtpfLlH761jblnyuzOhIAD0UZAdBonru1t+LCg3TsTLH+d8VOzq4BUCvKCIBGExzgp7l3D5DdZuiD7ZlakXrM6kgAPBBlBECjGhTXSjOv6ypJevr9XUo/WWRxIgCehjICoNH9ZHQXDe0YrrMl5frp8lSVVbisjgTAg1BGADQ6u83Qn+4eoOAAh1LTz+hlTvcF8B2UEQBNon1YoH7zP30lSa98eUBfpXG6L4BKlBEATWZ8/2jdPihGLlN6fPk25RVzui8AygiAJvbcrb3VoXXl6b5PrfiG030BUEYANK2W/g79+e6BctgMfbQjS+9u5XRfoLmjjABocgNiw/T49d0kSc+8v1OHcwstTgTASpQRAJaYcU1nDYsPV2FphR5ZtpW7+wLNGGUEgCXsNkN/mjhA4S2c2nksX4+8mapyrj8CNEuUEQCWiQ4L1MIpg+XvsOmLb3P06/d3cUAr0AxRRgBYalBcK730w4EyDGnZV+mat+qg1ZEANDHKCADLje3dVs+O7y1J+sNne7Ui9ajFiQA0JcoIAI8wJbGjfnx1J0nSL97ZoXUHci1OBKCpUEYAeIwnb+ihm/u1U1mFqRl/36Jvs/OtjgSgCVBGAHgMm83QH+/qr6Hx4SooKdfUxV8rK6/Y6lgAGhllBIBH8XfY9ddJg9WlTUtl55/TtCVfK/8c97ABfBllBIDHCQ3y0+vThigy2F/fZhfooX9sUWk51yABfBVlBIBHimkVpCVTh6iF0651B07qZ//crjIuigb4JMoIAI/Vp32o5v0oQQ6boQ+3Z+r+v21WYQmXjQd8DWUEgEe7plukFkxOUICfTav3ndA9f92ok2dLrI4FoAFRRgB4vGt7RGnZA8PVKshP24/m6Y7XNijjVJHVsQA0EMoIAK8wMK6V3nkoUe3DApWWW6jb5q/XzmN5VscC0AAoIwC8RufIlnrvJ4nq0TZYJwpKdPeCjVrPlVoBr0cZAeBVokIC9PaMERreKVxnS8o1ZclX+nB7ptWxAFwByggArxMS4Ke/3TdUN/WtvHT8o8tStXhtmtWxANQTZQSAV/J32PXSDwdqyogOkqTnP9qt//toNxdHA7wQZQSA17LbDD17S289Mba7JGnR2jT9z7x12n+8wOJkANxBGQHg1QzD0MOju+i1HyWoVZCfdmXm6+aX12rJujS5XKbV8QDUAWUEgE+4oU9bfTbzal3TLVIl5S499+FuTVnylbLzzlkdDcBlUEYA+Iw2IQF6fdoQPX9rbwX42bRmf67Gzk3RxzuyrI4G4BIoIwB8imEYmjyioz56dJT6tg9VXnGZHn5zq2Yt36b8c2VWxwNQC8oIAJ/UpU3lBdIevbaLbIb0XuoxjZu7Ruu4SBrgcSgjAHyWn92mnyV11z9njFBceJCOnSnWvQs3aeqSr7Qrk0vJA56CMgLA5yV0CNfKn47SlBEd5LAZWrX3hG56aa0eXZaqw7mFVscDmj3DNE2PP/ctPz9foaGhysvLU0hIiNVxAHixIycL9WLyPr2/rfIS8g6boYlDYvXYdV0VFRJgcTrAt9T19zdlBECztCszTy98tldf7j0hSQrws2nayHjNuLqzQoP8LE4H+AbKCADUwVdpp/T7T7/V5iOnJUkhAQ5NGxmvu4fGql1ooMXpAO9GGQGAOjJNU198m6Pff7pXe89fSt5mSNd0i9TEIbG6rmeU/OwcYge4izICAG6qcJn6+JssvbnpiDYeOlW1PKKlU7cPitFdQ2LVObKlhQkB70IZAYArkJZbqLc3Z+idLUd1oqCkavmQjq00cUicbuzbVkFOh4UJAc9HGQGABlBW4dKX3+bo7c0Z+uLbHF24957TbtOQ+Fa6umukru4WqR5tg2UYhrVhAQ9DGQGABpadd07vbj2qtzdn6MjJomo/iwz216iuEbq6a6Su6hqhiJb+FqUEPAdlBAAaiWmaOniiUCn7TmjN/hPaeOiUissqqo3pHR2i4Z1aq1e7EPVuH6LOkS05CBbNDmUEAJpISXmFthw+rdX7T2jNvlztzsqvMcZpt6lb25bq1S6k8hEdqh7tghUSwDVN4LsatYzMmzdPf/jDH5SVlaXevXtr7ty5GjVq1EXHr169WrNmzdKuXbsUHR2tX/ziF5oxY0ad10cZAeBNcgrOad2BXG3PyNPurHztycxXQUl5rWPbhQYoplWg2ocFKqZVUOXfW1X+PTosQP4OexOnBxpOo5WR5cuXa9KkSZo3b55Gjhypv/zlL1q4cKF2796tuLi4GuPT0tLUp08fPfDAA3rwwQe1bt06/eQnP9GyZct0++23N+hkAMATmaapo6eLtSszX7uz8rU7M197svJ17EzxZV/bJthfbUMDFN7CqfAgp1q1cFb+vYVTrYL++/eQQIda+jsU6GfnQFp4jEYrI8OGDdOgQYM0f/78qmU9e/bUhAkTNGfOnBrjf/nLX+qDDz7Qnj17qpbNmDFD27dv14YNG+q0TsoIAF90pqhUh08W6ejpIh09Xayjp4t07HTx+b8X1zgOpS4MQwrys6uFv0Mt/B0KctrVwulQkL9dgX52OR02+TtscjpsctqrP/d32ORnt8luM+SwGXLYbXLYjKrndpshh92Qzaj8u80wZBiSzbiwTDLO/91mSIYqf34h14XxF5Yb55dX/k3/HVs1F6Pa8++OqVxu1Lr8SjTnHhfewtngp6vX9fe3W2stLS3Vli1b9OSTT1ZbnpSUpPXr19f6mg0bNigpKanasrFjx2rRokUqKyuTn1/N70tLSkpUUvLf8/rz82t+/woA3i4syKkBQU4NiA2r8TPTNHWqsFRHTxfrREGJThWV6nRhqU6df5wuuvBnmU6eLVFBSblMUzJNqbC0QoWlFdJ3ro8CXM5LPxyoW/pHW7Jut8pIbm6uKioqFBUVVW15VFSUsrOza31NdnZ2rePLy8uVm5urdu3a1XjNnDlz9Nxzz7kTDQB8imEYat3SX63reIqwaZoqLqtQYUmFikrLdbakXEWlFSo8/+fZknKdK6tQablLJecfpRceFRUqKXOptMKlsgqXKlymKlymyi/8WXHhuatqucusXGeFy5TLNGWaUoVZ+XeXS1XLTF34s7IoqdpzU2ZV/srnlSMu/E/VH9XmWWPuNcbU8v9PjVH14/mnfNSf3cLdQvXaH/P97yNN07zkd5S1ja9t+QWzZ8/WrFmzqp7n5+crNja2PlEBoFkwDENBTsf53exc4wTexa0yEhERIbvdXmMvSE5OTo29Hxe0bdu21vEOh0OtW7eu9TX+/v7y92djAgCgOXDrCjxOp1MJCQlKTk6utjw5OVmJiYm1vmbEiBE1xn/++ecaPHhwrceLAACA5sXtywHOmjVLCxcu1OLFi7Vnzx49/vjjSk9Pr7puyOzZszV58uSq8TNmzNCRI0c0a9Ys7dmzR4sXL9aiRYv085//vOFmAQAAvJbbx4xMnDhRJ0+e1PPPP6+srCz16dNHK1euVIcOHSRJWVlZSk9PrxofHx+vlStX6vHHH9err76q6OhovfTSS3W+xggAAPBtXA4eAAA0irr+/uauTQAAwFKUEQAAYCnKCAAAsBRlBAAAWIoyAgAALEUZAQAAlqKMAAAAS1FGAACApSgjAADAUm5fDt4KFy4Sm5+fb3ESAABQVxd+b1/uYu9eUUYKCgokSbGxsRYnAQAA7iooKFBoaOhFf+4V96ZxuVzKzMxUcHCwDMNosPfNz89XbGysMjIyfPaeN74+R+bn/Xx9jr4+P8n358j86s80TRUUFCg6Olo228WPDPGKPSM2m00xMTGN9v4hISE++Q/su3x9jszP+/n6HH19fpLvz5H51c+l9ohcwAGsAADAUpQRAABgqWZdRvz9/fXMM8/I39/f6iiNxtfnyPy8n6/P0dfnJ/n+HJlf4/OKA1gBAIDvatZ7RgAAgPUoIwAAwFKUEQAAYCnKCAAAsJTPlZF58+YpPj5eAQEBSkhI0Jo1ay45fvXq1UpISFBAQIA6deqk1157rcaYd999V7169ZK/v7969eqlFStWNFb8y3Jnfu+9956uv/56RUZGKiQkRCNGjNBnn31Wbczrr78uwzBqPM6dO9fYU6mVO/NbtWpVrdm//fbbauM86fOT3Jvj1KlTa51j7969q8Z40meYkpKi8ePHKzo6WoZh6F//+tdlX+NN26C78/PGbdDdOXrbduju/LxtG5wzZ46GDBmi4OBgtWnTRhMmTNDevXsv+zqrt0OfKiPLly/XzJkz9dRTTyk1NVWjRo3SuHHjlJ6eXuv4tLQ03XjjjRo1apRSU1P1q1/9So899pjefffdqjEbNmzQxIkTNWnSJG3fvl2TJk3SXXfdpU2bNjXVtKq4O7+UlBRdf/31WrlypbZs2aLRo0dr/PjxSk1NrTYuJCREWVlZ1R4BAQFNMaVq3J3fBXv37q2WvWvXrlU/86TPT3J/jn/+85+rzS0jI0Ph4eG68847q43zlM+wsLBQ/fv31yuvvFKn8d62Dbo7P2/bBiX353iBt2yH7s7P27bB1atX6+GHH9bGjRuVnJys8vJyJSUlqbCw8KKv8Yjt0PQhQ4cONWfMmFFtWY8ePcwnn3yy1vG/+MUvzB49elRb9uCDD5rDhw+ven7XXXeZN9xwQ7UxY8eONe++++4GSl137s6vNr169TKfe+65qudLliwxQ0NDGyriFXF3fl9++aUpyTx9+vRF39OTPj/TvPLPcMWKFaZhGObhw4erlnnSZ/hdkswVK1Zccoy3bYPfVZf51caTt8Hvq8scvXE7vKA+n6E3bYOmaZo5OTmmJHP16tUXHeMJ26HP7BkpLS3Vli1blJSUVG15UlKS1q9fX+trNmzYUGP82LFjtXnzZpWVlV1yzMXes7HUZ37f53K5VFBQoPDw8GrLz549qw4dOigmJkY333xzjf9qawpXMr+BAweqXbt2uu666/Tll19W+5mnfH5Sw3yGixYt0pgxY9ShQ4dqyz3hM6wPb9oGG4Inb4NXylu2wyvlbdtgXl6eJNX4N/ddnrAd+kwZyc3NVUVFhaKioqotj4qKUnZ2dq2vyc7OrnV8eXm5cnNzLznmYu/ZWOozv+/74x//qMLCQt11111Vy3r06KHXX39dH3zwgZYtW6aAgACNHDlS+/fvb9D8l1Of+bVr104LFizQu+++q/fee0/du3fXddddp5SUlKoxnvL5SVf+GWZlZemTTz7R/fffX225p3yG9eFN22BD8ORtsL68bTu8Et62DZqmqVmzZumqq65Snz59LjrOE7ZDr7hrrzsMw6j23DTNGssuN/77y919z8ZU3yzLli3Ts88+q/fff19t2rSpWj58+HANHz686vnIkSM1aNAgvfzyy3rppZcaLngduTO/7t27q3v37lXPR4wYoYyMDL3wwgu6+uqr6/WeTaG+eV5//XWFhYVpwoQJ1ZZ72mfoLm/bBuvLW7ZBd3nrdlgf3rYNPvLII9qxY4fWrl172bFWb4c+s2ckIiJCdru9RkvLycmp0eYuaNu2ba3jHQ6HWrdufckxF3vPxlKf+V2wfPlyTZ8+XW+//bbGjBlzybE2m01Dhgxp8kZ/JfP7ruHDh1fL7imfn3RlczRNU4sXL9akSZPkdDovOdaqz7A+vGkbvBLesA02JE/eDuvL27bBRx99VB988IG+/PJLxcTEXHKsJ2yHPlNGnE6nEhISlJycXG15cnKyEhMTa33NiBEjaoz//PPPNXjwYPn5+V1yzMXes7HUZ35S5X+NTZ06VW+++aZuuummy67HNE1t27ZN7dq1u+LM7qjv/L4vNTW1WnZP+fykK5vj6tWrdeDAAU2fPv2y67HqM6wPb9oG68tbtsGG5MnbYX15yzZomqYeeeQRvffee/riiy8UHx9/2dd4xHbYIIfBeoi33nrL9PPzMxctWmTu3r3bnDlzptmiRYuqo56ffPJJc9KkSVXjDx06ZAYFBZmPP/64uXv3bnPRokWmn5+f+c4771SNWbdunWm3283f/va35p49e8zf/va3psPhMDdu3Ojx83vzzTdNh8Nhvvrqq2ZWVlbV48yZM1Vjnn32WfPTTz81Dx48aKampprTpk0zHQ6HuWnTJo+f35/+9CdzxYoV5r59+8ydO3eaTz75pCnJfPfdd6vGeNLnZ5ruz/GCH/3oR+awYcNqfU9P+gwLCgrM1NRUMzU11ZRkvvjii2Zqaqp55MgR0zS9fxt0d37etg2apvtz9Lbt0N35XeAt2+BDDz1khoaGmqtWrar2b66oqKhqjCduhz5VRkzTNF999VWzQ4cOptPpNAcNGlTtdKYpU6aY11xzTbXxq1atMgcOHGg6nU6zY8eO5vz582u85z//+U+ze/fupp+fn9mjR49qG1lTc2d+11xzjSmpxmPKlClVY2bOnGnGxcWZTqfTjIyMNJOSksz169c34Yyqc2d+v/vd78zOnTubAQEBZqtWrcyrrrrK/Pjjj2u8pyd9fqbp/r/RM2fOmIGBgeaCBQtqfT9P+gwvnOZ5sX9z3r4Nujs/b9wG3Z2jt22H9fk36k3bYG1zk2QuWbKkaownbofG+fAAAACW8JljRgAAgHeijAAAAEtRRgAAgKUoIwAAwFKUEQAAYCnKCAAAsBRlBAAAWIoyAgAALEUZAXBRzz77rAYMGGDZ+n/961/rxz/+caO9f05OjiIjI3Xs2LFGWweAy+MKrEAzdblbf0+ZMkWvvPKKSkpKqu7c2ZSOHz+url27aseOHerYsWOjrWfWrFnKz8/XwoULG20dAC6NMgI0U9+9Hfjy5cv19NNPa+/evVXLAgMDFRoaakU0SdJvfvMbrV69Wp999lmjruebb77R0KFDlZmZqVatWjXqugDUjq9pgGaqbdu2VY/Q0FAZhlFj2fe/ppk6daomTJig3/zmN4qKilJYWJiee+45lZeX64knnlB4eLhiYmK0ePHiaus6duyYJk6cqFatWql169a69dZbdfjw4Uvme+utt3TLLbdUW/aDH/xAjz76qGbOnKlWrVopKipKCxYsUGFhoaZNm6bg4GB17txZn3zySdVrTp8+rXvvvVeRkZEKDAxU165dtWTJkqqf9+3bV23bttWKFSvq/38mgCtCGQHgli+++EKZmZlKSUnRiy++qGeffVY333yzWrVqpU2bNmnGjBmaMWOGMjIyJElFRUUaPXq0WrZsqZSUFK1du1YtW7bUDTfcoNLS0lrXcfr0ae3cuVODBw+u8bO//e1vioiI0FdffaVHH31UDz30kO68804lJiZq69atGjt2rCZNmqSioiJJlced7N69W5988on27Nmj+fPnKyIiotp7Dh06VGvWrGng/6cA1BVlBIBbwsPD9dJLL6l79+6677771L17dxUVFelXv/qVunbtqtmzZ8vpdGrdunWSKvdw2Gw2LVy4UH379lXPnj21ZMkSpaena9WqVbWu48iRIzJNU9HR0TV+1r9/f/3v//5v1boCAwMVERGhBx54QF27dtXTTz+tkydPaseOHZKk9PR0DRw4UIMHD1bHjh01ZswYjR8/vtp7tm/f/rJ7agA0HofVAQB4l969e8tm++9/x0RFRalPnz5Vz+12u1q3bq2cnBxJ0pYtW3TgwAEFBwdXe59z587p4MGDta6juLhYkhQQEFDjZ/369auxrr59+1bLI6lq/Q899JBuv/12bd26VUlJSZowYYISExOrvWdgYGDVnhQATY8yAsAtfn5+1Z4bhlHrMpfLJUlyuVxKSEjQ0qVLa7xXZGRkreu48DXK6dOna4y53PovnCV0Yf3jxo3TkSNH9PHHH+vf//63rrvuOj388MN64YUXql5z6tSpi2YB0Pj4mgZAoxo0aJD279+vNm3aqEuXLtUeFztbp3PnzgoJCdHu3bsbJENkZKSmTp2qf/zjH5o7d64WLFhQ7ec7d+7UwIEDG2RdANxHGQHQqO69915FRETo1ltv1Zo1a5SWlqbVq1frpz/9qY4ePVrra2w2m8aMGaO1a9de8fqffvppvf/++zpw4IB27dqljz76SD179qz6eVFRkbZs2aKkpKQrXheA+qGMAGhUQUFBSklJUVxcnG677Tb17NlT9913n4qLixUSEnLR1/34xz/WW2+9VfV1S305nU7Nnj1b/fr109VXXy273a633nqr6ufvv/++4uLiNGrUqCtaD4D646JnADySaZoaPny4Zs6cqR/+8IeNtp6hQ4dq5syZuueeexptHQAujT0jADySYRhasGCBysvLG20dOTk5uuOOOxq17AC4PPaMAAAAS7FnBAAAWIoyAgAALEUZAQAAlqKMAAAAS1FGAACApSgjAADAUpQRAABgKcoIAACwFGUEAABY6v8HRRJXl6EOdxIAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "ts = np.linspace(0, 2, 51)\n", "lcce = calc.compute(ts, method='cce')\n", "plt.plot(ts, lcce.real)\n", "plt.xlabel('Time (ms)')" ] }, { "cell_type": "markdown", "id": "a5fc0706-fda4-42d7-8a0e-c5e0acbb3011", "metadata": {}, "source": [ "## Dissipation in the bath" ] }, { "cell_type": "markdown", "id": "26645225", "metadata": {}, "source": [ "In this example we assume that each electron spin in the bath decays into completely random state with the characteristic decay time $T_1=0.5$ ms.\n", "To add dissipation we use the method `.add_single_jump` method of the `calc.bath` object. Note because we have a single bath type (same name of all spins) the `which` keyword of the method is unnecessary." ] }, { "cell_type": "code", "execution_count": 6, "id": "af8f360d", "metadata": { "tags": [] }, "outputs": [], "source": [ "et1 = 0.5 # in ms \n", "decay_rate = 1 / et1 / 2 # in rad / ms\n", "calc.bath.add_single_jump('p', rate=decay_rate) # p for creation operator S+ \n", "calc.bath.add_single_jump('m', rate=decay_rate) # m for annihilation operator S-" ] }, { "cell_type": "markdown", "id": "22ecf3a7-dc81-4121-9720-4026438d96d3", "metadata": {}, "source": [ "The dissipators are stored in the `.so` attribute of the given spin type in the units of $kHz^{1/2}$ \n", "(square root to match how one sets up the simulations in the Qutip, but note that it's not in radial frequencies):" ] }, { "cell_type": "code", "execution_count": 7, "id": "b3b03225-f87b-42b8-8680-27e36402a5cd", "metadata": { "tags": [] }, "outputs": [ { "data": { "text/plain": [ "{'p': 0.3989422804014327, 'm': 0.3989422804014327}" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "calc.bath['e'].so" ] }, { "cell_type": "markdown", "id": "bd68f2a2", "metadata": {}, "source": [ "We can compute the coherence with ME-CCE to account for these dissipative dynamics:" ] }, { "cell_type": "code", "execution_count": 8, "id": "453e48ed-3ab4-4570-850f-455e7bd897e0", "metadata": { "tags": [] }, "outputs": [], "source": [ "lmecce = calc.compute(ts, method='mecce')\n", "\n", "calc.order = 1\n", "lmecce_1 = calc.compute(ts, method='mecce') # Coherence at first order" ] }, { "cell_type": "markdown", "id": "b1ed9e4b-d5d6-449f-8819-4f5c24a7123b", "metadata": {}, "source": [ "By comparing the full calculation with ME-CCE to the product of ME-CCE calculation of first order (ME-CCE1) and CCE2 we can directly see the interplay between single spin incoherent dynamics and the coherent flips:" ] }, { "cell_type": "code", "execution_count": 9, "id": "b3941aac", "metadata": { "tags": [] }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAGwCAYAAABVdURTAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAB/KUlEQVR4nO3dd1xV5R/A8c+5l3vZGwQHAiriXjhx71WOyixNczT8WZnZtJ0N29nSMleluXJmapl7m4ob90AFRED2vvf8/rhKEUO4wGV936/XfXk55znP8z3gla/PeYaiqqqKEEIIIUQloSnrAIQQQgghSpIkN0IIIYSoVCS5EUIIIUSlIsmNEEIIISoVSW6EEEIIUalIciOEEEKISkWSGyGEEEJUKlZlHYClGY1GwsPDcXR0RFGUsg5HCCGEEIWgqiqJiYnUqFEDjabgvpkql9yEh4fj4+NT1mEIIYQQwgxXr16lVq1aBZapcsmNo6MjYPrmODk5lXE0QgghhCiMhIQEfHx8sn+PF6TKJTd3HkU5OTlJciOEEEJUMIUZUiIDioUQQghRqUhyI4QQQohKRZIbIYQQQlQqVW7MjRBCiNwMBgOZmZllHYao4vR6/V2neReGJDdCCFGFqapKZGQkcXFxZR2KEGg0Gvz9/dHr9cWqR5IbIYSowu4kNtWqVcPOzk4WNxVl5s4iuxEREdSuXbtYfxcluRFCiCrKYDBkJzbu7u5lHY4QeHp6Eh4eTlZWFjqdzux6ZECxEEJUUXfG2NjZ2ZVxJEKY3HkcZTAYilWPJDdCCFHFyaMoUV6U1N9FSW6EEEIIUamUaXKzY8cO7r33XmrUqIGiKKxevfqu12zfvp2goCBsbGyoU6cO3333XekHKoQQQogKo0yTm+TkZJo3b84333xTqPKXLl1iwIABdO7cmZCQEF599VUmTZrEihUrSjlSIYQQQlQUZTpbqn///vTv37/Q5b/77jtq167NjBkzAGjYsCEHDx7k008/5f777y+lKAvHkJVFxNVzxGTEomh1aDVWaDVaNBodWq0OjaJFb2WD3soWvY0ddtY6dFp5KiiEEOaKjIzk/fff5/fff+f69etUq1aNFi1aMHnyZHr27AlASEgIH3zwATt27CA+Pp7atWvTtWtXXnzxRerXr8/ly5fx9/fPs/69e/fSvn17S96SKCEVair43r176dOnT45jffv2Ze7cuWRmZuY5bSw9PZ309PTsrxMSEkolttib4bj+1IH+fj75lumblMynN2MwqAqxOPJEDWf0qhW2Rj02Rlv0OKJTnNFr3bC3DcTDpwc+Hq74uttR290OJxvzp8UJIURlcvnyZTp27IiLiwsff/wxzZo1IzMzkz/++IOnnnqK06dPs27dOu6//3769u3LokWLqFu3LlFRUSxfvpw33niDpUuXZtf3119/0bhx4xxtyPT4iqtCJTeRkZF4eXnlOObl5UVWVhbR0dFUr1491zXTp0/nnXfeKfXYFCBZtUZvVDEooALG/4z6vtNPo1VU7JREztu4AAYg9fYrNrtsUOpGXtz+FidVP/YbA/jMqhqRShus3f3x9bDH182Oel6OtKjlgo+brcx2EEIUm6qqpGYWbwquuWx12iL9OzZx4kQUReHAgQPY29tnH2/cuDHjxo0jJSWFsWPHMmDAAFatWpV93t/fn3bt2uVakdnd3R1vb+9i34coHypUcgO5p4mpqprn8TumTp3KlClTsr9OSEjAxyf/3hVzeVSvDW/f4NB/YlONBozGLIzGTBRVBUVDRkoCabeu81HEfmISr3MrJYq49FjiMhKIN6QQp6bRJh2slSxaKeepa3WRGb61cDBuplGaAZvrLlw758+2lDYcU+vjZm9N81rONPdxoYWPC81rueBqX7ylq4UQVU9qpoFGb/5RJm2fmtYXO33hfiXFxsayceNG3n///RyJzR0uLi6sWrWK6OhoXnrppTzrcHFxKU64opyrUMmNt7c3kZGROY5FRUVhZWWVb/ehtbU11tbWlggvF0VRULRWaLRWgE32cb21A+6uNRhQp03+F6sq3LoE1w5y+dImHG7tJ0mj4YCdBuySgRM4GY7ySKKKPq452872ZOuZGtmX+7nb0fx2otPcx4XGNZyw0WlL72aFEMJCzp8/j6qqNGjQIN8y586dAyiwzL8FBwfn2rAxPj4erVb+3ayIKlRy06FDB3777bccx/78809at25drGWayyVFAbc64FaH5s0eZJfRwNmbxzl84XcOReznUNJVYrWwxgWmZu3kTXUDV20CWWvsyIKEIC7HwOWYFNYcCQfASqMQ6O1Is1ou2b08AdUcsJJBzUKI22x1Wk5N61tmbRfW3Xrs/12msJYuXUrDhg1zHJPEpuIq0+QmKSmJ8+fPZ3996dIljhw5gpubG7Vr12bq1Klcv36dn376CYAJEybwzTffMGXKFB5//HH27t3L3LlzWbx4cVndgsVoNVoaerWgoVcLRgIGo4HdYVtYdWQ292gNkLQNn7Qz1LW/ysM1VtHHuj7hDv1Yk96K/eFZRCelczI8gZPhCSw+YKrTRqehSQ1nmtR0JsDLgYBqjgRUc5BHWkJUUYqiFPrRUFkKCAhAURRCQ0MZMmRInmXq168PwOnTp+nQocNd6/Tx8aFevXolGaYoQ4pa1PS2BG3bto3u3bvnOv7oo4+yYMECxowZw+XLl9m2bVv2ue3bt/Pcc89x8uRJatSowcsvv8yECRMK3WZCQgLOzs7Ex8fj5ORUErdRPiRHw8lVPHRyJic1WQB0SE3loeRMujUYRlTj8YQkunD0WjzHrsVx7Fo8SelZeVbl4aA3JTpeDgRUc6BeNUcCvR1xk6RHiEolLS2NS5cu4e/vj42Nzd0vKEf69+/P8ePHOXPmTK5xN3Fxceh0Ovz8/OjUqVOOAcX/LuPi4pI9FTwkJIQWLVpYKHqRn4L+Thbl93eZJjdlodImN5i6Ybdc3cLyEz+x5+Zh7vxgG6ZnMPFWPF1r90Dp8BT4BmNU4WJ0MseuxXE6MpFzNxI5eyOJ63Gp+dZfw9mGxjWdaVzDKbvHx8vJWmZqCVFBVeTk5tKlSwQHB+Pm5sa0adNo1qwZWVlZbNq0iVmzZhEaGsqaNWsYNmwY/fr1Y9KkSdSrV4/o6GiWLVtGWFgYS5YsyU5u8poK7uLiUuG+LxWdJDdmqszJzb9dS7zG8jPLWXJ6ESkG0zo/9ycm8XZ0LFRvDh2ehkZDwCpnb0xyehYXbiZx7kYS56KSTElPVCJXY/NOejwc9DSq4UyTGk609nOlc4CnLE4oRAVRkZMbgIiICN5//33WrVtHREQEnp6eBAUF8dxzz9GtWzcADh48yPTp09m5c2f2bNkePXrw4osvUq9evQIX8Vu8eDEPPfSQBe9ISHJjpqqS3NxxK+0WP578kV9CFzHDLpDgU5sgK41MwMqxOkrbxyFoLNi5FVhPYlomp26P2TkRHs/J6wmci0rE+J+/PW72egY1r8F9rWrStKaz9OoIUY5V9ORGVD6S3JipqiU3d8Snx+Okd0JJiYWD85h1ch57tAaeuhVHO6MOpd0ECH4abF0LXWdapoHTkYmcuB7P8WvxbD4dRXTSP6tB16vmwH2tajKkRU1quNiWxm0JIYpBkhtR3khyY6aqmtz8W6Yxkz7LexOdFgNAm9Q0Xo69RaBiBx2egvb/A5uif2+yDEZ2no9m5eHr/HkykvQsI2Ca1R5c1537WtaiXxNv7K3L/2wMIaoCSW5EeSPJjZkkuTGJSoli7vG5LD+7nExjJhoVhiUm8vSteFysnaHjs9D2CdDnXv2zMBLSMtlwPIIVh69z4NI/20o42+p4bUBDhrWuJY+shChjktyI8kaSGzNJcpNTRFIEnx36jD8um5Zcd1Jh+o0ouqSmgb0ndHoOWo8DnfmPla7GprAq5DorDl/jSkwKAB3quPPBfU3x9zAveRJCFJ8kN6K8KankRqa1VHHVHarzaddPmdd3HgGuAaRqdfh1exNc/SD5JvzxKnzVEv6eC0bzNtTzcbNjUs8ANk/pytT+DbDRadh7MYa+M3bw7dbzZBqMJXtTQgghqjTpuRHZsoxZnIg+QYtqLcCQCUd+4df9n9ApJhxvgwFqtoZBX4FX47vWVZCwmBReW32cneeiAWjg7cj0+5rSsnbhBzMLIYpPem5EeSM9N6LEWWmsTIkNgFbHCd8gptlrGOTryxx3TwzXD8L3XeCvdyAz/8X+7qa2ux0/jWvL5w82x9VOx+nIRO6btYe3157Md9VkIYQQorAkuRH5stHa0LJaS1LVLL50suXxOg2JUlTY9TnMCoZLO8yuW1EU7mtVi83Pd+O+ljVRVViw5zK9P9/OltM3SvAuhBBCVDWS3Ih81XOtx4J+C3iv43vYWdnxt5rMsDr12eNaHWIvwo/3wpqnICX27pXlw81ez+fDW/Dz+Lb4uNkSEZ/GuAUHmb3jQgneiRBCiKpEkhtRIEVRGFxvMEvvWUqgayCxWclMcNHzXaPbG56GLIRv28LxX6EYw7c6B3jy5+SujO7gC8AH60/zwfpQjP9dAlkIIYAxY8agKEqeGydPnDgRRVEYM2ZMjrL/ffXr1++u7WzdupUBAwbg7u6OnZ0djRo14vnnn+f69evZZVRVZfbs2bRr1w4HBwdcXFxo3bo1M2bMICXFNEP07bffzjOGBg0alMw3ROQgyY0oFD9nPxYNXMSD9R9ERcWr5aMw7g/wCDTNqloxHhY/DKlxZrdhq9cybXATpvY3fdhn77jIC78eldlUQog8+fj4sGTJElJT/xkDmJaWxuLFi6ldu3aOsv369SMiIiLHa/HixQXW//3339OrVy+8vb1ZsWIFp06d4rvvviM+Pp7PPvssu9yoUaOYPHkygwcPZuvWrRw5coQ33niDNWvW8Oeff2aXa9y4ca4Ydu3aVULfDfFvslSsKDRrrTVvdHiDQfUG0cyjmWnp4Qk7Sd7xCfa7v4SzG2Bub3h4CbjXNbudJ7vWxc1ezysrj7Py8HXiUjL5dkQrbPXaErwbIURF16pVKy5evMjKlSsZOXIkACtXrsTHx4c6derkKGttbY23t3eh67527RqTJk1i0qRJfPHFF9nH/fz86NKlC3FxcQAsW7aMRYsWsXr1agYPHpyj3KBBg0hISMg+ZmVlVaQYhPmk50YUWXPP5tmrC8dlpTI0dgdfdHmcTKeaEH0W5vSEy8X738iw1j7MHhWEjU7DltNRjJyzj7iUjJIIXwiRH1WFjOSyeZn5WHvs2LHMnz8/++t58+Yxbty4Yn8rli9fTkZGBi+99FKe511cXABYtGgRgYGBORKbOxRFwdnZudixiKKTnhtRLJvDNhORHMG85N84Ur8534R74hh+BH4aAvfOgJaPmF13z4ZeLHqsHWPn/83hsDiGfbeXn8a3pbqzbMIpRKnITIEPapRN26+Gm7Xdy6hRo5g6dSqXL19GURR2797NkiVL2LZtW45y69atw8HBIcexl19+mTfeeCPPes+dO4eTkxPVq1cvsP1z584RGBhYqFiPHz+eK4aHHnqIOXPmFOp6UXiS3Ihiub/+/TjqHXlrz1scjjnBU7Wa852zD3ahv5lmUkWfhZ5vg8a8TsIgXzeWTwhm9Lz9nItK4v6Ze/hpfFvqVXMs2RsRQlRIHh4eDBw4kB9//BFVVRk4cCAeHh65ynXv3p1Zs2blOObm5gbAhAkTWLhwYfbxpKQkVFUt1P53hS0HEBgYyNq1a3Mcc3SUf8tKgyQ3otj6+PXBx9GH8X+OJ+TmUSZ5t+Pbzs9jvfMz2P0lRJ+H+2aDtcPdK8tDoLcjK/4XzOh5B7h4M5kHvtvL/DFtZEVjIUqazs7Ug1JWbZtp3LhxPP300wB8++23eZaxt7enXr16eZ6bNm0aL7zwQo5j9evXJz4+noiIiAJ7b+rXr09oaGih4tTr9fnGIEqWjLkRJaKhe0Nm9ZqFnZUd+yP3M4UbZA75DrR6OPM7zO8H8dfvXlE+arna8euEYJrXciYuJZORc/Zz/Fp8Cd6BEAJFMT0aKotXIXs/8tKvXz8yMjLIyMigb9++Rb6+WrVq1KtXL/sF8MADD6DX6/n444/zvObOgOIRI0Zw9uxZ1qxZk6uMqqrEx8u/U2VBkhtRYpp7Nuebnt9go7UhLCGM+MDe8Og6sPOAyOPwQw+4fsjs+t3s9fzyeHuC67qTkmFg/I9/Ex5n/jYQQojKQavVEhoaSmhoKFpt3rMq09PTiYyMzPGKjo7Ot04fHx+++OILvvzyS8aPH8/27du5cuUKu3fv5sknn+Tdd98F4MEHH2T48OE8/PDDTJ8+nYMHD3LlyhXWrVtHr1692Lp1a3adWVlZuWK4cUNWZC8N8lhKlKg23m2Y2Wsm/s7+eNh6QG0PeHwL/DIcbobCgntg9BrwaWtW/fbWVnw3KogHZu3h7I0kxi34m1//F4yDtfxVFqIqu9tGihs3bsz1eCkwMJDTp0/ne83EiROpX78+n376KUOHDiU1NRU/Pz/uuecepkyZAphmRP3yyy/Mnj2befPm8d5772FlZUVAQACjR4/O0ZN08uTJXDFYW1uTlpZW1NsVdyG7gotSdyL6BI3tfVCWj4aL28DGBcZthGoNza7z2q0Uhny7h+ikdLoFejJndGustNIRKURRyK7goryRXcFFhbD6/GpGrh/Jx8e/Qx2+CGq1gbQ4+Pk+iAszu95arnbMebQ1NjoN287cZNq6U1SxPF0IIUQ+JLkRpUpVVYyqkYWhC/n65DwYsQw8G0BiOPw8FJLzf+Z9Ny18XPjiwRYoCvy09wrzd18uucCFEEJUWJLciFI1NGAor7V7DYAfjv/AnAsr4ZGV4OwDMedh0QOQnmh2/f2bVueVfqa9qN79/RR/nZLBeUIIUdVJciNK3UMNHmJKkGnw3ZeHv2R7wnkYtQrs3CE8BJaMhKx0s+t/oksdHm7rg6rCM4tDOHFdpl4KIURVJsmNsIixTcbySEPTVgxv7nmTaHtXGLkcdPZwaTusfAKMBrPqVhSFaYOb0DnAg9RM0xTxiHiZIi6EEFWVJDfCYiYHTSbANYDYtFjWXlgLNYPgoUWg0cGp1bD+BbM3z9NpNXw7shUB1Ry4kZDOuAUHSUrPKtkbEEIIUSFIciMsxlprzYedP2Ra8DTGNh5rOli3u2lrBhQ4OA+2fWh2/U42OuaNaYOHg57QiAQmLQ7BYJQZVEIIUdVIciMsqr5rfYYGDM250VyT+2Dgp6b32z+E/bPNrt/HzY4fRrfG2krDltNRzN11sZgRCyGEqGgkuRFlJj49ni8Pf0mmIRPaPAbdpppObHgJLmwxu96WtV15e1BjAD794yxnIs2fjSWEEKLikeRGlAlVVXnsz8eYc3wO3x65vYtv15eh5ShAhRWPQ0KE2fU/1MaHHg2qkWEwMnnpETKyjCUTuBBCiHJPkhtRJhRFYUKzCQDMOzGPvyP/Nu0KPOAT8GoCKdGw4jEwmDcoWFEUPry/Ka52OkIjEvhy89mSDF8IUYbGjBlj+jdkwoRc5yZOnIiiKIwZMyZX+f+++vXrd9e2tm7dyoABA3B3d8fOzo5GjRrx/PPPc/369ewyqqoye/Zs2rVrh4ODAy4uLrRu3ZoZM2aQkpICwNtvv51nDA0aNMiuZ+XKlfTt2xcPDw8UReHIkSOF/p5069at0GVXrFhBt27dcHZ2xsHBgWbNmjFt2jRiY2Ozy2RkZPDxxx/TvHlz7Ozs8PDwoGPHjsyfP5/MzEzg7t/X2NhYnnnmGQIDA7Gzs6N27dpMmjTJIjulS3IjykxP357cF3AfKiqv7nqVhIwE0NnCsB9B7wBXdsG26WbXX83Rhg+GNgVg1rYLHLpyq6RCF0KUMR8fH5YsWUJq6j/LPqSlpbF48WJq166dq3y/fv2IiIjI8Vq8eHGBbXz//ff06tULb29vVqxYwalTp/juu++Ij4/ns88+yy43atQoJk+ezODBg9m6dStHjhzhjTfeYM2aNfz555/Z5Ro3bpwrhl27dmWfT05OpmPHjnz4YeEmVpw9e5YlS5bkOHb48GF+++23fK957bXXGD58OG3atGHDhg2cOHGCzz77jKNHj/Lzzz8DpsSmb9++fPjhhzzxxBPs2bOHAwcO8NRTT/H1119z8uTJQn1fw8PDCQ8P59NPP+X48eMsWLCAjRs3Mn78+ELdX7GoVUx8fLwKqPHx8WUdilBVNTkjWe2/or/aZEET9cXtL/5z4thyVX3LSVXfclbVc5uK1cbkJSGq78vr1K4fb1GT0zOLF7AQlUhqaqp66tQpNTU1taxDKZJHH31UHTx4sNq0aVN14cKF2ccXLVqkNm3aVB08eLD66KOP5ipfFFevXlX1er06efLkPM/funVLVVVVXbp0qQqoq1evzlXGaDSqcXFxqqqq6ltvvaU2b968UG1funRJBdSQkJACy8XExKhPPPGEOmzYMLV58+bqm2++qfbr108NDQ3Ns/z+/ftVQJ0xY0aB9/TRRx+pGo1GPXz4cK4yGRkZalJSkqqq5n1fly1bpur1ejUzM+9/iwv6O1mU39/ScyPKlJ3Ojg87f4hW0bLh0gbWXVxnOtH0AWg9DlBNC/wlhJvdxtuDGlPd2YbLMSl8sD60ZAIXohJLyUzJ95VuSC902bSstEKVNdfYsWOZP39+9tfz5s1j3LhxZtf3b8uXLycjI4OXXnopz/MuLi4ALFq0iMDAQAYPHpyrjKIoODs7l0g8eXFzc8vuXTp69Cjnz59nw4YNOR51/duiRYtwcHBg4sSJeZ7/9z316tWLli1b5iqj0+mwt7c3O+Y7O3pbWVmZXUdhlG7tQhRCM89mPNn8SWYemcnc43Pp79cfrUYLfafDtb8h8jj8Oh4e/Q20Rf8r62yr45MHmvPI3P0s3BdG70bedK3vWQp3IkTl0O6Xdvme61yzMzN7zcz+utuybqRm5b0ieGuv1szv90/y0W9FP26l5348fPzR42bFOWrUKKZOncrly5dRFIXdu3ezZMkStm3blqvsunXrcHBwyHHs5Zdf5o033siz7nPnzuHk5ET16tULjOHcuXMEBgYWKt7jx4/niuGhhx5izpw5hbr+v27dusVrr71GdHQ0zZs3p27duvTv358ZM2bkGdO5c+eoU6cOOp2uwHrPnTtX6DE8Rfm+xsTE8O677/Lkk08Wqu7ikORGlAuPN30cVVUZ2XCkKbEB0NmYxt983xXC9sDW96HXW2bV3ynAgzHBfizYc5kXlx/lz+e64GKnL8E7EEJYmoeHBwMHDuTHH39EVVUGDhyIh4dHnmW7d+/OrFmzchxzc3MDYMKECSxcuDD7eFJSEqqq5lyPKx+FLQcQGBjI2rVrcxxzdHQs1LV5iYqKonPnzjz88MN069aNadOmcfjwYc6ePZtnclMa91TQ9/XfEhISGDhwII0aNeKtt8z7d7woJLkR5YKVxoqJLfLoKnWvC4O+gl/Hwq7PwbcjBPQyq42X+zVgx7mbXLyZzBtrTvL1w7m7XIUQsH/E/nzPZf/n47ZtD27Lt6xGyTnyYeP9G4sVV17GjRvH008/DcC3336bbzl7e3vq1auX57lp06bxwgsv5DhWv3594uPjiYiIKLD3pn79+oSGFu5xt16vzzcGcwQGBuZKYlq1akWrVq3yLF+/fn127dpFZmZmgb03Rbmngr6vdyQmJtKvXz8cHBxYtWrVXXuOSoKMuRHljqqqbA7b/M+z+Cb3mRb5A1j5OMRfz//iAtjqtXz+YAu0GoXfjoaz9qj543iEqMzsdHb5vqy11oUua2NlU6iyxdGvXz8yMjKyZ/iYo1q1atSrVy/7BfDAAw+g1+v5+OOP87wmLi4OgBEjRnD27FnWrFmTq4yqqhaZ9gzk+Sjuv0aMGEFSUhIzZ87M8/y/7+mvv/4iJCQkV5msrCySk5MLHVdCQgJ9+vRBr9ezdu1abGxs7n5RCZDkRpQ7r+9+nclbJzP3xNx/DvZ5H6o3h9RY+HWc2evftPBx4anupn+83lh9gsj4tLtcIYQoz7RaLaGhoYSGhqLVavMtl56eTmRkZI5XdHR0vuV9fHz44osv+PLLLxk/fjzbt2/nypUr7N69myeffJJ3330XgAcffJDhw4fz8MMPM336dA4ePMiVK1dYt24dvXr1YuvWrdl1ZmVl5Yrhxo0b2edjY2M5cuQIp06dAuDMmTMcOXKEyMjI4n6bAGjXrh0vvfQSzz//PC+99BJ79+7lypUrbN68mWHDhvHjjz8CMHnyZDp27EjPnj359ttvOXr0KBcvXmTZsmW0a9eOc+fOFer7mpiYSJ8+fUhOTmbu3LkkJCRklzEYDCVyT/kq0hyuSkCmgpd/my5vUpssaKIG/RykhieG/3Mi5oKqflDLNEX8zzfNrj8jy6De89VO1ffldeqouftVo9FYAlELUfFU9Kng+clrKjiQ6xUYGHjXtjZt2qT27dtXdXV1VW1sbNQGDRqoL7zwghoe/s+/TQaDQZ01a5bapk0b1c7OTnVyclKDgoLUL7/8Uk1JSVFV1TQVPK8YrK2ts+uZP39+nmXeeuutIn+PCrJ06VK1S5cuqqOjo2pvb682a9ZMnTZtWvZUcFVV1bS0NHX69Olq06ZNVRsbG9XNzU3t2LGjumDBguxp3Hf7vm7dujXP84B66dKlPGMrqangiqqqVWrb5ISEBJydnbOno4nyR1VVxv0xjoM3DjLAfwAfdfnon5MnV8PyR03vH10H/p3NauN8VCIDvtpFRpaRLx9qweAWNYsfuBAVTFpaGpcuXcLf399ijwuEKEhBfyeL8vtbHkuJckdRFF5s8yIKCusvrefYzWP/nGw8BILGmN6vfQYyzFsjo141R56+/Xjq441nSMss5S5SIYQQFiPJjSiXGrk3YlDdQQB8/PfH5Ohg7P0uONWEW5dM08PN9HjnOlR3tuF6XCpzd10qbshCCCHKCUluRLk1qdUkbK1sOXrzKBsv/2sKqY0T3POF6f2+mXDtkFn12+q1vNzPtJLnzK3niUqUwcVCCFEZSHIjyq1qdtUY12Qc9Vzq4WH7n4W56veFpg+CaoQ1T0FWet6V3MWg5jVoXsuZ5AwDn/8pO4cLIURlIMmNKNfGNRnH8nuX08a7Te6T/T4EOw+4GQo7Pzerfo1G4Y17GgGw9OBVToUnFCdcISqkKjavRJRjJfV3UZIbUa7ptXqsNPkspG3vDgM+Mb3f+SncOGlWG6393BjYrDqqCu/9fkr+oRdVxp2VYlNSzN+8UoiSlJGRAVDgmkWFIdsviAoh3ZDOwlMLiU6N5uW2L/9zovFQOLECTq8zPZ4a/5dZm2u+0q8Bm07dYM+FGP4KjaJ3I68SjF6I8kmr1eLi4kJUVBQAdnZ2hd5TSIiSZjQauXnzJnZ2dsXeNVySG1EhnIk9w4zDM1BQGFR3EA3dG5pOKAoM/Awu7YTwENj3LXR8tsj1+7jZMb6TP7O2XeCD9aF0re+J3ko6NkXl5+3tDZCd4AhRljQaDbVr1y52ki2L+IkK46XtL7Hh8gbaeLdhbp+5Of/yH/4Z1j4NVjYwYTd4FH1zusS0TLp/uo3opAzevKcR4zr5l2D0QpRvBoOBzMzMsg5DVHF6vR6NJu//WBbl97ckN6LCCE8KZ9DqQaQb0vmy+5f0qN3jn5OqCj8PhYtbTTuHP7oO8vmAFOSX/WG8uuo4zrY6tr/YDRc7fQnegRBCCHPJCsWiUqrhUIPRjUYD8NnBz8g0/Ot/mYoC934JOnu4shsOzTOrjeFtfGjg7Uh8aiYz/jp39wuEEEKUO5LciAplfNPxuNu4E5YYxm8Xf8t50tUXer1ler/pLYi7WuT6tRqF1weapoYv3HeFCzeTihuyEEIIC5PkRlQo9jp7xjYZC8Dc43MxqsacBdo8Dj7tISMJ1k02Pa4qok4BHvRoUI0so8r09aElELUQQghLkuRGVDjD6g/jgfoP8FWPr9Ao//krrNHAoK9Baw3n/4JTq81q49UBDbHSKPwVGsXu89HFD1oIIYTFlHlyM3PmzOytzYOCgti5c2eB5RctWkTz5s2xs7OjevXqjB07lpiYGAtFK8oDO50db3V4i7oudfMu4FkfOj1nev/nG2btHF6vmgOPtPcF4N11pzAYq9S4eyGEqNDKNLlZunQpkydP5rXXXiMkJITOnTvTv39/wsLC8iy/a9cuRo8ezfjx4zl58iTLly/n77//5rHHHrNw5KI8yTTmMX2147PgVAvir8Ker8yq99meATjZWHE6MpF1x8KLGaUQQghLKdPk5vPPP2f8+PE89thjNGzYkBkzZuDj48OsWbPyLL9v3z78/PyYNGkS/v7+dOrUiSeffJKDBw/m20Z6ejoJCQk5XqJyuJlyk9d2vcao9aNyb5mgt4M+75re75ph1uBiV3s9T3SpA8A3W85jlN4bIYSoEMosucnIyODQoUP06dMnx/E+ffqwZ8+ePK8JDg7m2rVrrF+/HlVVuXHjBr/++isDBw7Mt53p06fj7Oyc/fLx8SnR+xBlR6/V89eVvzgZc5Lt17bnLtB4KPh2gqxU2PSGWW2MDvbDycaKc1FJbDwZWcyIhRBCWEKZJTfR0dEYDAa8vHLu4ePl5UVkZN6/RIKDg1m0aBHDhw9Hr9fj7e2Ni4sLX3/9db7tTJ06lfj4+OzX1atF/x+8KJ+crZ15qMFDAMw+Njt3742iQL/poGjg5Cq4vKvIbTjZ6BjT0bRS8ddbzsummkIIUQGU+YDi/+4foapqvntKnDp1ikmTJvHmm29y6NAhNm7cyKVLl5gwYUK+9VtbW+Pk5JTjJSqP0Y1GY6O14Xj0cfZG7M1doHozaPWo6f2GV8BoKHIb4zr6Ya/XEhqRwF+hsv+OEEKUd2WW3Hh4eKDVanP10kRFReXqzblj+vTpdOzYkRdffJFmzZrRt29fZs6cybx584iIiLBE2KKccbd154H6DwDw/dHv8y7U4w2wcYYbx+Hwj0Vuw8VOz+hgPwC+3nJOem+EEKKcK7PkRq/XExQUxKZNm3Ic37RpE8HBwXlek5KSkmtDLa1WCyC/cKqwMY3HoNPoOBx1mIOReQwut3eHbq+a3m9+F1JvFbmN8Z38sdFpOHYtnu1nbxYzYiGEEKWpTB9LTZkyhTlz5jBv3jxCQ0N57rnnCAsLy37MNHXqVEaPHp1d/t5772XlypXMmjWLixcvsnv3biZNmkTbtm2pUaNGWd2GKGNe9l4MrTcUgJ9O/ZR3oTbjwbMBpMbCtg+L3IaHgzUj25nWvZGxN0IIUb5ZlWXjw4cPJyYmhmnTphEREUGTJk1Yv349vr6mXyIRERE51rwZM2YMiYmJfPPNNzz//PO4uLjQo0cPPvroo7K6BVFOjGs6Dk87T0Y0HJF3Aa0O+n0IPw+BAz9A0Bio1rBIbTzZpQ4/77vCoSu32HshhuB6HsWOWwghRMlT1Cr2X9CibJkuKqElI+H0OqjTDUatNs2oKoI315zgp71X6FDHncVPtC+VEIUQQuRWlN/fZT5bSoiSpqoqKZn5bLnQ5z3TvlMXt8Hp34tc94SuddFpFfZejOHg5djiBSqEEKJUSHIjKpWTMSd5ZMMjvL779bwLuPlD8NOm93+8CplpRaq/hostDwTVAuCrLeeLE6oQQohSIsmNqFT0Gj3Hbh5j05VNnL+VT/LRaQo4Voe4K7D3myK38b+u9dBqFHacvcmRq3HFC1gIIUSJk+RGVCoBrgH0qt0LgDkn5uRdyNoBek8zvd/5GSQUbY2k2u52DGlRE4BvtpwzO1YhhBClQ5IbUek83uxxADZc2sDVxHy222g6DGq1hcwU2F702XZPda+LRoG/QqM4GR5fnHCFEEKUMEluRKXTyL0RwTWCMapGlpxeknchRYHe75jeH/4JYi4UqY06ng7c08y0ttI3MvZGCCHKFUluRKU0suFIAFadW5X/zCnfYAjoC6oBtrxb5Dae7lEPgA0nIjl7I9HsWIUQQpQsSW5EpdSpZidqO9YmMTORjZc35l+w55uAYto1PDykSG3U93KkX2NvQHpvhBCiPJHkRlRKGkXDlKApfNHtCwbVHZR/Qe8mpvE3AJunFbmdO703646Fcyk62ZxQhRBClDBJbkSl1dO3J718e2GlucsuI91fBY0OLmyBi9uL1EaTms70aFANowpzd10sRrRCCCFKiiQ3okowqsb8T7r5Q+uxpvd/vQ1F3JHk8c51APj10DVuJWeYGaEQQoiSIsmNqNRUVeWHYz/Qd0VfribkMy0coMuLoLOH8MMQ+luR2mhfx40mNZ1IyzSycN+VYkYshBCiuCS5EZWaoigcjjpMZHIki88szr+gQzXo8JTp/eZpYMgqUht3em9+3HuFtExDcUIWQghRTJLciEqvUNPCwbTnlK0bxJyDo78UqY0BTatT3dmG6KR01h4JL064QgghikmSG1HpBdcIxtfJl6TMJNZeWJt/QRtn6Py86f3W6ZCZWug2dFoNYzv6AfDDzouoRRy3I4QQouRIciMqPY2i4eEGDwPwy+lfCh5c3OYxcKoFieFw4IcitfNQ29o4WFtxLiqJbWdvFidkIYQQxSDJjagSBtcdjL3Onkvxl9gXvi//gjob6D7V9H7nZ5AaV+g2nGx0PNTGB4A5O2VauBBClBVJbkSV4KB3YEi9IYCp96ZAzR4Cj0BIi4M9XxWpnbGd/NFqFHafj5ENNYUQooxIciOqjIcbPMyguoOY2GJiwQW1Vre3ZQD2zYLEyEK3UdPFlgFNqwMwd+clc0MVQghRDJLciCrD18mX9zu9TyP3Rncv3GAg1GwNmSmw45MitfN4Z38A1h4NJzI+zZxQhRBCFIMkN0LkRVGg19um94cWQGzhx9A0q+VCW383sowqC/ZcLo3ohBBCFECSG1HlnL91nrf2vMWvZ38tuKB/Z6jbE4xZsOPTIrVxZ1G/RfuvkJRe+AUBhRBCFJ8kN6LK+fvG36w8t5IfT/5Y8LRwMG2qCXB0CcRcKHQbPRtUo46HPYlpWSz7u4BtH4QQQpQ4SW5ElTOo7iAcdA5cTrjM3vC9BReu1Rrq9QbVUKSxNxqNwvjbY2/m7b5EluEuSZQQQogSI8mNqHLsdfbZ08IXhS66+wXdbq97c2wpRJ8vdDv3t6qFm72ea7dS+ePkDTMiFUIIYQ5JbkSV9HCDh1FQ2Hl9J2EJYQUXrhUEAX1BNcKOjwvdho1OyyPtfQHZkkEIISxJkhtRJdV2qk1wzWAAVp1fdfcLur1i+vP4cog+V+h2RnfwRW+l4cjVOA5duWVOqEIIIYpIkhtRZd0fcD8Aq8+vJst4lxlNNVtB/f6m3pvtHxW6DQ8Ha+5rWRMw9d4IIYQofZLciCqrW61uNHZvzLD6w8gwZNz9gjt7Th3/FW6eKXQ7j90eWPznqRtcjk42J1QhhBBFIMmNqLJ0Wh1L7lnCxBYTsdPZ3f2C6s2hwT2ACts+LHQ79ao50j3QE1VFFvUTQggLkORGiKK4M/bm5CqICi30ZeM6mXpvfj10jcS0zNKITAghxG2S3IgqL9OYyeawzWy6sunuhb2bQsN7KWrvTad6HgRUcyApPYtlB6+ZH6wQQoi7kuRGVHkbLm1g8tbJzDg0o3DTtbve7r05tRpunCxUG4qiMKajHwA/7rmMwSjTwoUQorRIciOqvF61e2FnZUdYYhgHbxy8+wXeTaDRYNP7IvTe3NeyFs62OsJiU9hyOsrMaIUQQtyNJDeiyrPT2dHfvz8AK8+tLNxFXV8BFAhdC5HHC3WJrV7LQ219AJi/+5I5oQohhCgESW6E4J81bzZd2URCRsLdL/BqBI2HmN4XofdmdAc/tBqFPRdiOB1ZiHaEEEIUmSQ3QgBNPJoQ4BpAuiGd9RfXF+6iO703p9dBxLFCXVLTxZZ+jb0BmL/rsnnBCiGEKJAkN0JgGvB7X737gCI8mqrWAJqYrilK783Y2wOLVx+5TmxyIRYPFEIIUSSS3Ahx2z117kGn0aFRNIV7NAXQ9WVAgTO/Q/iRQl0S5OtK05rOpGcZWXzgLpt2CiGEKDJJboS4zcXGhfX3rWfJPUtw0jsV7iLPQGj6gOn9jk8KdYmiKNm9Nz/vvUKmwWhGtEIIIfIjyY0Q/+Jt7130i7q8SPbYm0LOnBrYrDoeDtZEJqSx4URk0dsUQgiRL0luhMhDQkYCl+ILOV3bMxAaDzW93/5xoS6xttIyqr0vAPN2ybRwIYQoSZLcCPEfm8M202NZD6btnVb4i7q8aPozdC3cOFWoS0a0q41eq+HI1ThCwm6ZEakQQoi8SHIjxH80dm9MpjGTgzcOciXhSuEu8mr0z6rFOwrXe+PpaM29zWsAMH/3ZTMiFUIIkRdJboT4D297bzrW6AjAqnOrCn9hl5dMf55cDVGnC3XJnYHF649HEBmfVoQohRBC5EeSGyHycF+Aaf2aNRfWkGnMLNxF3k2gwT2AWuiZU01qOtPWz40so8rCfYXsJRJCCFEgSW6EyEPXWl1xs3EjOjWandd2FuHC2703J1bAzbOFumRcJz8AFu2/QlqmoYiRCiGE+C9JboTIg06rY3Bd0xiaIj2aqt4cAgcAKuz8tFCX9G7kTU0XW26lZLLmyHUzohVCCPFvktwIkY8hAUMA2B2+m/j0+MJfeKf35vhyiLlw1+JajcKjwaZp4fN3X0ZV1aKGKoQQ4l8kuREiH3Wc6/BBpw/YeP9GnK2dC39hjZYQ0BdUI+z8rFCXDG9dG1udltORiey9GGNmxEIIIUCSGyEKdG/de6lmV63oF3Z92fTn0SUQe/GuxZ3tdNzXqiYAP+2RgcVCCFEcktwIUUhFelxUKwjq9QLVADs/L9Qljwb7AfDnqUiux6WaEaEQQgiQ5EaIuzoYeZAn/nyCTw4Wbnp3tuzem8Vw6+69MfW9HAmu645RRaaFCyFEMUhyI8RdpGSlsDdiL79f/L3wa94A+LSFOt3BmAW7itZ7s+RAmEwLF0IIM0lyI8RddKjRATcbN2LTYtlzfU/RLr7TexOyCOKu3rV4zwbVsqeF/3Y03IxohRBCSHIjxF3oNDoG+A8AYO2FtUW72LcD+HcBYybs+uKuxa20Gh65vVv4j3tlWrgQQpijzJObmTNn4u/vj42NDUFBQezcWfBqsOnp6bz22mv4+vpibW1N3bp1mTdvnoWiFVXVoLqDANh2dVvR1ryBf/Xe/Azxd1+kb3gbH/RWGk5cT+BwWFzR2hJCCFG2yc3SpUuZPHkyr732GiEhIXTu3Jn+/fsTFhaW7zUPPvggmzdvZu7cuZw5c4bFixfToEEDC0YtqqIGbg0IcA0gw5jBn1f+LNrFfp3AtxMYMgrVe+Nmr2fw7d3Cf9xz2YxohRCiaivT5Obzzz9n/PjxPPbYYzRs2JAZM2bg4+PDrFmz8iy/ceNGtm/fzvr16+nVqxd+fn60bduW4ODgfNtIT08nISEhx0uIolIUhUF1TL03v134regVdLvde3P4R0i4+1iaOwOL1x+PICpBdgsXQoiiMDu5uXDhAq+//joPP/wwUVFRgCn5OHnyZKGuz8jI4NChQ/Tp0yfH8T59+rBnT96DNteuXUvr1q35+OOPqVmzJvXr1+eFF14gNTX/NUGmT5+Os7Nz9svHx6eQdyhETgPqDKCNdxuG1hta9LEwfp2hdvDt3psZdy3epKYzQb6uZBlVfjmQf0+mEEKI3MxKbrZv307Tpk3Zv38/K1euJCkpCYBjx47x1ltvFaqO6OhoDAYDXl5eOY57eXkRGRmZ5zUXL15k165dnDhxglWrVjFjxgx+/fVXnnrqqXzbmTp1KvHx8dmvq1fvPmNFiLxUs6vGvL7zGBowFEVRinaxovzTe3NoASRE3PWSO703i/aHkZFlLFp7QghRhZmV3Lzyyiu89957bNq0Cb1en328e/fu7N27t0h1/feXhKqq+f7iMBqNKIrCokWLaNu2LQMGDODzzz9nwYIF+fbeWFtb4+TklOMlRJnw7wo+7cGQDru/vGvxfo29qeZozc3EdDacuHsyJIQQwsSs5Ob48eMMHTo013FPT09iYgq36Z+HhwdarTZXL01UVFSu3pw7qlevTs2aNXF2/mcTw4YNG6KqKteuXSvCHQhhvpjUGBaeWsipmFNFuzBH7818SMy7h/IOvZWGEe1qA/DTXlmxWAghCsus5MbFxYWIiNz/kwwJCaFmzZqFqkOv1xMUFMSmTZtyHN+0aVO+A4Q7duxIeHh49mMwgLNnz6LRaKhVq1YR7kAI8315+Es++vsjlp1ZVvSL63QHn3aQlVao3psR7Wqj0yocunKLE9eLOAVdCCGqKLOSmxEjRvDyyy8TGRmJoigYjUZ2797NCy+8wOjRowtdz5QpU5gzZw7z5s0jNDSU5557jrCwMCZMmACYxsv8u74RI0bg7u7O2LFjOXXqFDt27ODFF19k3Lhx2NramnMrQhTZvXXvBeCPy3+QllXEmUyK8s+6NwfnQeKNAotXc7RhQNPqACyQaeFCCFEoZiU377//PrVr16ZmzZokJSXRqFEjunTpQnBwMK+//nqh6xk+fDgzZsxg2rRptGjRgh07drB+/Xp8fU0rtEZERORY88bBwYFNmzYRFxdH69atGTlyJPfeey9fffWVObchhFmCvIKoYV+DpMwktl3dVvQK6vaAWm1MvTd77v53d3QHPwDWHg0nNjmj6O0JIUQVo6jFWN/94sWLHD58GKPRSMuWLQkICCjJ2EpFQkICzs7OxMfHy+BiYbavQ75m9rHZdK7ZmZm9Zha9gnN/waL7wcoWJh8Dh2r5FlVVlUHf7Ob49Xhe6hfIxG71ihG5EEJUTEX5/V2sRfzq1KnDAw88wIMPPlghEhshSsq9dUyPpvaE7yE6NbroFdTrCTWDICv1rr03iqJkTwtfuPcKWQaZFi6EEAUxK7l54IEH+PDDD3Md/+STTxg2bFixgxKivPNz9qOZZzMMqoH1F9cXvQJFga6vmN4fmANJNwssfk+z6rjZ6wmPT+Ov0CgzIhZCiKrD7EX8Bg4cmOt4v3792LFjR7GDEqIiGFx3MFYaK2LSCrf8QS4BvaFGy0L13tjotDzUxrS6tuw3JYQQBTMruUlKSsqxeN8dOp1O9m4SVcY9de5h67CtPBf0nHkVKAp0m2p6//ccSC748dYj7X3RKLD3YgxnIhPNa1MIIaoAs5KbJk2asHTp0lzHlyxZQqNGjYodlBAVgZ3ODhcbl+JVEtDH1HuTmQJ7vi6waA0XW/o08gZg4T5Z1E8IIfJjZc5Fb7zxBvfffz8XLlygR48eAGzevJnFixezfPnyEg1QiIrgSsIVajrUxEpTxI/UnXVvFj8EB36A4Elg755v8VEdfNl4MpJVIdd5uX8DHKzN+ggLIUSlZlbPzaBBg1i9ejXnz59n4sSJPP/881y7do2//vqLIUOGlHCIQpRvz2x5hntW3cO+iH3mVVC/H1RvDpnJsPebAosG13Wnjqc9SelZrA65bl57QghRyZk9FXzgwIHs3r2b5ORkoqOj2bJlC127di3J2ISoEGrY1wDgtwu/mVfBv1ctPjAbUmILKKrwSDvTIpcL912hGMtUCSFEpVWsdW4yMjK4du0aYWFhOV5CVCV3tmPYEraFpIyku5TOR+AA8G4KGUl3HXtzf1AtbHQaTkcmcvDKLfPaE0KISsys5ObcuXN07twZW1tbfH198ff3x9/fHz8/P/z9/Us6RiHKtcbujfF39ifNkMamK5vufkFe/j1zav/3Bc6ccrbVMbi5aYNaGVgshBC5mZXcjBkzBo1Gw7p16zh06BCHDx/m8OHDhISEcPjw4ZKOUYhyTVGU7BWL111cZ35FgQP+GXtzl3VvRnUwPZpafzyC6KR089sUQohKyKzk5siRI3z//ff079+fFi1a0Lx58xwvIaqae+rcA8CByAOEJ4WbV4miQLdXTe8P/FDgqsVNajrTwseFTIPK0r+vmteeEEJUUmYlN40aNSI62oz9dISopKo7VKetd1sA1l8yYzuGO+r3hRqtTOve7J5RYNFR7U29N7/sD8NglIHFQghxh1nJzUcffcRLL73Etm3biImJISEhIcdLiKrosaaP8UnXT3ik4SPmV6Io0P12783fcyHxRr5FBzarjoudjutxqWw7I/tNCSHEHYpqxlxSjcaUEymKkuO4qqooioLBYCiZ6EpBUbZMF6JMqCrM6QXXD0L7idBver5FP1gfyuwdF+kW6MmCsW0tGKQQQlhWUX5/m7W86datW80KTAhRCHd6bxbeBwfnQcdnwdE7z6Ij29Vm9o6LbD97kysxyfi621s4WCGEKH/MSm5ksT4h8paalcrPp35m57WdzOs7D51WZ15FdXuATzu4uh92fQH9P8qzmK+7PV3re7L97E1+2R/G1AENixG9EEJUDmYv4rdz504eeeQRgoODuX7dtAz8zz//zK5du0osOCEqGp1Gx5LTSzhy8wg7r+80v6J/j705OB8S8p+BdWdg8bKDV0nLLL+PhIUQwlLMSm5WrFhB3759sbW15fDhw6Snm9bZSExM5IMPPijRAIWoSKw0VgzwHwAUYzuGO/y7Qu1gMKTDzs/zLda9QTVquthyKyWT9ccjitemEEJUAmYlN++99x7fffcdP/zwAzrdP93uwcHBsoifqPLubMew/dp24tPjza/o3703h3+E+Gt5FtNqFEa0qw3Az7JisRBCmJfcnDlzhi5duuQ67uTkRFxcXHFjEqJCC3QLJNA1kExjJn9c/qN4lfl3Br/OYMiAnZ/lW+zB1j7otAohYXGcuF6MhEoIISoBs5Kb6tWrc/78+VzHd+3aRZ06dYodlBAV3Z3em7UX1ha/sjt7Th3+GeLy3pjW09Ga/k2qA7LflBBCmJXcPPnkkzz77LPs378fRVEIDw9n0aJFvPDCC0ycOLGkYxSiwhngPwCNouHozaOEJeSdkBSaX0fT+BtjJuz4NN9id/abWnMknPjUzOK1KYQQFZhZU8Ffeukl4uPj6d69O2lpaXTp0gVra2teeOEFnn766ZKOUYgKx9POk561e2KjtUGlBLZG6P4qXNoORxZB5yng6perSGtfVwK9HDlzI5GVh68xtqN/8dsVQogKqMgrFBsMBnbt2kXTpk2xsbHh1KlTGI1GGjVqhIODQ2nFWWJkhWJhKXdW7C4xPw+FC1ug5SMw+Nu8i+y7whurT1DH057NU7qWbPtCCFGGivL7u8iPpbRaLX379iU+Ph47Oztat25N27ZtK0RiI4QllXhicWfH8COLIeZCnkWGtqyJvV7LxZvJ7L0QU7LtCyFEBWHWmJumTZty8eLFko5FiErpdOxpNl7eWPyKfNpAvd6gGmD7x3kWcbC2YmirmgAs2l/MsT5CCFFBmZXcvP/++7zwwgusW7eOiIgI2RVciHwciTrCsN+GMW3PNNIN6cWvsMdrpj+PLYWo03kWGdnONLD4j5ORRCWmFb9NIYSoYMxKbvr168fRo0cZNGgQtWrVwtXVFVdXV1xcXHB1dS3pGIWosJp5NsPb3pvEzES2Xd1W/AprtIQG9wAqbMt7NfCG1Z0I8nUly6iy/GDeC/8JIURlJruCC1GKNIqGe+vcyw/Hf+C3C7/R169v8Svt8Tqc/h1OrYGIo1C9ea4iI9vV5tCVW/yyP4wJXeui1cjAYiFE1VHk2VIVncyWEpZ2Mf4ig1cPxkqx4q9hf+Fu6178Slc8DseXQUBfGLks1+m0TAPtp28mLiWTeWNa06OBV/HbFEKIMlSqs6XukF3BhSicOs51aOrRlCw1i98v/l4ylXZ7BRQtnPsDrh7IddpGp2VYUC0AFu2TgcVCiKpFdgUXwgIG1x0MwOoLqymRzlL3utBypOn9lnfzLPJwW9NmmlvORHHtVkrx2xRCiApCdgUXwgL6+fdDr9ETmxpLdGp0yVTa5SXQ6uHSDri4PdfpOp4OdKznjqrCkgNXS6ZNIYSoAGRXcCEswNnamZ8G/MSmYZvwtPMsmUpdfCBorOn9lnchjx6hO9PCl/x9lUyDsWTaFUKIck52BRfCQhq7N0an0d29YFF0fh6sbOHa33Duz1ynezfywtPRmuikdDadulGybQshRDklu4ILYWEGo4GY1BLaGsHRC9o9YXq/5V0w5uyd0Wk1PNTGB4BF+6+UTJtCCFHOmZXcvPTSSwwZMoTu3buTlJREly5deOyxx3jyySdlV3AhCrAvYh99V/TltV2vlVylHSeD3hEij0Po2lynH2pbG40Cu8/HcPFmUsm1K4QQ5ZTZU8Hff/99oqOjOXDgAPv27ePmzZu8+27eszaEECY17WtyI+UGe8L3EJkcWTKV2rlBh6dM77d+AEZDzjZdbOkeWA2AX2S/KSFEFWB2cgPIruBCFJGPkw9BXkGoqKy9kLuXxWwdJoKNC0SfgePLc51+pL1pYPGvh6+RlmnIdV4IISoTs5Kb5ORk3njjDYKDg6lXrx516tTJ8RJC5G9ovaEArD5fQmveANg4Q6fJpvfbpoMhM8fpLvU9qeliS1xKJuuPR5RMm0IIUU6ZtbfUY489xvbt2xk1ahTVq1dHUWTfGiEKq7dvbz7Y/wFXE69y6MYhWnu3LpmK2z4Be2fCrcsQshBaj80+pdUojGhXm0/+OMPCfVe4r1WtkmlTCCHKIbOSmw0bNvD777/TsWPHko5HiErPTmdHX7++rDq/itXnV5dccqO3N00N3/gy7PgEmj8MOpvs08Na1+KLTWc5HBbHqfAEGtWQvdWEEJWTWY+lXF1dcXNzK+lYhKgyhtQbAsCfV/4kJbMEt0YIGgNONSHhOvw9J8epao429G3iDcAvB2RauBCi8jIruXn33Xd58803SUmR/WqEMEfLai15vOnjzOkzB1sr25KrWGdj2lQTYOdnkBaf4/TIdqb9plYdvk5SelbJtSuEEOVIoR9LtWzZMsfYmvPnz+Pl5YWfn1+O/aUA2V9KiLtQFIVJrSaVTuXNR8CeryH6rOnPHq9nn+pQx506HvZcjE5mzZHr2dszCCFEZVLo5GbIkCGlGIYQosRoraDHG7BsFOz9Fto8blrJGFNSNaJdbd77PZSF+8IY0ba2TAgQQlQ6ilpic1ErhoSEBJydnYmPj8fJSQZUirJ17tY5FoUuopZjLR5r+ljJVayqMKcXXD9oSm4Gfpp9Ki4lg3YfbCY9y8jKicG0qu1acu0KIUQpKcrv72It4nfo0CEWLlzIokWLCAkJKU5VQlRJF+MvsuLcChafXozBWIKL6ykK9Hrb9P7QfIi9mH3KxU7PPc1qALBwnwwsFkJUPmYlN1FRUfTo0YM2bdowadIknn76aYKCgujZsyc3b94s6RiFqLS6+3TH2dqZqJQo9kXsK9nK/TtDvV5gzIIt7+c4NaqDaazNumMRxCZnlGy7QghRxsxKbp555hkSEhI4efIksbGx3Lp1ixMnTpCQkMCkSaU0SFKISkiv1TPAfwBgWrG4xPV80/TniV8h4mj24ea1nGla05mMLCPLD14t+XaFEKIMmZXcbNy4kVmzZtGwYcPsY40aNeLbb79lw4YNJRacEFXBnTVvNodtJj49vuDCRVW9OTR5wPR+87Tsw4qiMOr2flML91/BaKxSQ++EEJWcWcmN0WjMNf0bQKfTYTQaix2UEFVJQ7eGBLoGkmnMZP2l9SXfQI/XQGMF5/+CSzuzD9/bvAbOtjquxqay/Zw8ThZCVB5mJTc9evTg2WefJTw8PPvY9evXee655+jZs2eJBSdEVaAoSnbvTak8mnKrY1q5GOCvt0wzqQBbvZZhQaY9pn7eKwOLhRCVh1nJzTfffENiYiJ+fn7UrVuXevXq4e/vT2JiIl9//XVJxyhEpTewzkBqO9amc83OZBlLYeXgLi+Bzg6uH4LQ37IPj7z9aGrrmSiuxsqK40KIyqFY69xs2rSJ06dPo6oqjRo1olevXiUZW6mQdW5EeaWqaukuqLflPdOGmh714X97TYv9AaPm7mfnuWgmdK3LK/0blF77QghRDKW2zs2WLVto1KgRCQkJAPTu3ZtnnnmGSZMm0aZNGxo3bszOnTvvUktOM2fOxN/fHxsbG4KCggp9/e7du7GysqJFixZFak+I8qrUVwoOfgZs3UzbMhz9JfvwnYHFyw5eJS2zBNfaEUKIMlKk5GbGjBk8/vjjeWZMzs7OPPnkk3z++eeFrm/p0qVMnjyZ1157jZCQEDp37kz//v0JCwsr8Lr4+HhGjx4t43tEpZNlzGJz2Gb2hO8p+cptnKHz86b3W6dDZioAPRpUo4azDbHJGWw4EVHy7QohhIUVKbk5evQo/fr1y/d8nz59OHToUKHr+/zzzxk/fjyPPfYYDRs2ZMaMGfj4+DBr1qwCr3vyyScZMWIEHTp0KHRbQlQEv4T+wuStk/km5JvSaaDNY+BUCxLD4cBsAKy0Gkbc3i1cBhYLISqDIiU3N27cyHMK+B1WVlaFXqE4IyODQ4cO0adPnxzH+/Tpw549+f+vdf78+Vy4cIG33nqrUO2kp6eTkJCQ4yVEeTWgzgCsFCuORx/ndOzpkm9AZwPdXzW93/k5pMYBMLxNbXRahcNhcZy4XsJr7QghhIUVKbmpWbMmx48fz/f8sWPHqF69eqHqio6OxmAw4OXlleO4l5cXkZGReV5z7tw5XnnlFRYtWoSVVeE2NJ8+fTrOzs7ZLx8fn0JdJ0RZ8LD1oKev6XHr8jPLS6eR5g+BZ0NIi4OdnwHg6WhN/yamz+6i/dJ7I4So2IqU3AwYMIA333yTtLS0XOdSU1N56623uOeee4oUwH8HUeY3Y8RgMDBixAjeeecd6tevX+j6p06dSnx8fPbr6lVZal6Ub8PqDwPg90u/k5JZCtOzNVro/Y7p/f7v4NZl4J/9plaHhBOfmlny7QohhIUUKbl5/fXXiY2NpX79+nz88cesWbOGtWvX8tFHHxEYGEhsbCyvvfZaoery8PBAq9Xm6qWJiorK1ZsDkJiYyMGDB3n66aexsrLCysqKadOmcfToUaysrNiyZUue7VhbW+Pk5JTjJUR51ta7Lb5OviRnJpfOisUAAX3AvysYMuAvU6LT2teVBt6OpGYaWHHoWum0K4QQFlCk5MbLy4s9e/bQpEkTpk6dytChQxkyZAivvvoqTZo0Yffu3XkmJnnR6/UEBQWxadOmHMc3bdpEcHBwrvJOTk4cP36cI0eOZL8mTJhAYGAgR44coV27dkW5FSHKLUVReCDAtB/U8rOl9GhKUaDv+4ACJ1fC1QMoisIjd/ab2neFYiyBJYQQZapwA1f+xdfXl/Xr13Pr1i3Onz+PqqoEBATg6upa5ManTJnCqFGjaN26NR06dGD27NmEhYUxYcIEwPRI6fr16/z0009oNBqaNGmS4/pq1aphY2OT67gQFd3geoP5KuQrFBQSMhJw0pdCj6N3U2g5EkIWwh+vwvhNDGlZkw83nOZidDJ7LsTQsZ5HybcrhBClrMjJzR2urq60adOmWI0PHz6cmJgYpk2bRkREBE2aNGH9+vX4+pr+9xgREXHXNW+EqIxcbVxZN3QdNRxqlG5D3V+HEyvh2t9wchUOTe7jvlY1+WnvFX7ee0WSGyFEhVSs7RcqItl+QYj/2PYhbJsOLr7w9N+cjcmgzxc70GoUdr3cnerOtmUdoRBClN72C0IIy0vISODcrXOl10DwM+BYHeKuwP7vqe/lSPs6bhiMKosPyOxCIUTFI8mNEOXY7uu76bmsJ6/uerX0Bvjq7aHHG6b3Oz6F5BhGtfcDYPGBMDINxtJpVwghSokkN0KUY43dG2NUjZyOPc3JmJOl11Dzh00DjNPjYfuH9GnshaejNTcT09lwIu9FNYUQoryS5EaIcszFxoXefr2BUpwWDqDRQJ/3Te8PzkN36wIjb+83NX/3pdJrVwghSoEkN0KUc3dWLN5waQOJGYml11CdrlC/HxizYNObjGzni16rISQsjsNht0qvXSGEKGGS3AhRzrWq1oo6znVIzUrl94u/l25jvd8FRQtn1uMZfYBBLUxT0efvvly67QohRAmS5EaIck5RlOzem+Vnl5fuysGe9aH1ONP7P19jbLDp0dT64xFExKeWXrtCCFGCJLkRogK4t+69WGutuRB3gSsJpbxrd7dXwNoJIo7S+ObG7GnhP+2V3cKFEBWDJDdCVADO1s583OVj/nzgT/yc/Uq3MXsP6Py86f3maTzWzhuAX/aHkZphKN22hRCiBEhyI0QF0aN2D6rZVbNMY+0mgHNtSAynR/QiarvZEZ+aycoQ2S1cCFH+SXIjRAWUkplSug3obKDvewBo9nzFMy1M29DN23UJo7FK7dgihKiAJLkRogIJTwrnyU1PMuy3YRiMpfyIqOEgqNMNDOkMjfoGR2srLtxMZuf56NJtVwghikmSGyEqEBdrF05EnyAsMYxd13eVbmOKAv0/Bo0VVuc38mpAGABzd8mifkKI8k2SGyEqEDudHfcH3A/AwtCFpd+gZyC0nwjAAze/wUbJYMfZm5y7UYqLCQohRDFJciNEBfNQg4fQKBr2RezjQtyF0m+w60vg4I0u/jIfVt8BwPw9l0u/XSGEMJMkN0JUMDUcatDDpwcAv4T+UvoNWjtCH9Pg4nvjf6EG0aw8fI1byRml37YQQphBkhshKqCRDUcC8NvF34hPjy/9Bps+AL4d0RrS+NhxKWmZRhb/HVb67QohhBkkuRGiAgryCiLQNZDUrFTWnF9T+g3eGVysaOmUuZuOmuP8tOcKmQZj6bcthBBFJMmNEBWQoig83fJppgVPY3iD4ZZp1LsJtH0cgPf0PxGTkMSGE5GWaVsIIYpAkhshKqhuPt0YGjAUa621BRudCvae+HOdMdqNzJNp4UKIckiSGyEqgVLdKfzfbF2g1zsAPGu1kvCrlzgcdssybQshRCFJciNEBbfszDIGrR7Emdgzlmmw+cNQqw0OShqv6hZJ740QotyR5EaICu5A5AEuJ1zml9MWmBYOoNHAgE9QURii3UP0ya1cj0u1TNtCCFEIktwIUcHdmRb++8XfuZVmoUdENVqitB4LwFva+czfec4y7QohRCFIciNEBdfCswUN3RqSbkhnxbkVlmu4xxtk6l1oqLmK7u/viZVF/YQQ5YQkN0JUcIqi8EijRwBYcnoJmcZMyzRs54ZV32kATFKWsWpzKW/kKYQQhSTJjRCVQD+/frjZuHEj5QZbwrZYrF2l1WhiPNthq2TQ6PCbJKZK740QouxJciNEJaDX6hlWfxgAi0IXWa5hRcH1wZmko6cDxzm4Zqbl2hZCiHxIciNEJTE8cDi9fXvzTMtnLNquxrMeZxo8BUCr05+QdivCou0LIcR/KarFVv8qHxISEnB2diY+Ph4nJ6eyDkeISiEzM4OLH7QjUL3IZe9++E1YWtYhCSEqmaL8/paeGyFEsel0es62+4AsVYNf5EYMob+XdUhCiCpMkhshKpnI5Eg+OvARc47PsWi7vXr0YZHmXgDS10yBtASLti+EEHdIciNEJXMi+gQLQxcy78Q8kjKSLNaurV5LWscXuWz0wi4tEvWvdyzWthBC/JskN0JUMj1q98DPyY/EjER+PfurRdt+uFMD3lWeBEA5OAfC9lm0fSGEAEluhKh0NIqGcU3GAfDTqZ/IMFhu7RknGx2BHQayNKsbAOraZyAzzWLtCyEESHIjRKV0T5178LLz4mbqTdZcWGPRtsd18udTHuGm6owSfRZ2fmbR9oUQQpIbISohnVbHmMZjAJh/Yj5ZxiyLte3hYM3Ato14M9PUPrs+hxsnLda+EEJIciNEJXVfwH24WLtwNfEqm65ssmjbj3epwyba8achCIxZsHYSGA0WjUEIUXVJciNEJWWns+Oxpo/xaKNHCfIKsmjbNV1sGdKyFm9kjiVFsYfrB2HvNxaNQQhRdckKxUKIUnE+KoneX2znQc1WPtL9ABodPLEVvJuWdWhCiApIVigWQpS5etUc6N/Em6WGbhx36AjGTFjxuMyeEkKUOkluhKgCDt84zIS/JhASFWLRdid2qwcojIsdjcHOE26GwmZZ3E8IUbokuRGiClh7YS27r+9m7vG5Fm23SU1nutT35KbRkXnuL5gO7psJF7ZYNA4hRNUiyY0QVcDYJmNRUNh+bTtnb521aNtTetcHYPp5H+IajzYdXD0RUmItGocQouqQ5EaIKsDXyZfevr0BLN5708LHhb6NvTCq8HrKcHAPgMQIWDcZqtZ8BiGEhUhyI0QVMb7peAA2Xt7I1cSrFm37hT6BaBRYFxrPmY6fg8YKTq2Bo0ssGocQomqQ5EaIKqKReyOCawRjVI38ePJHi7Yd4OXIfa1qATDtsB66vWI6sf5FuHXForEIISo/SW6EqEIea/oYAKvOrSI6NdqibU/uFYBeq2H3+Rh2eY0Gn/aQkQirnpTVi4UQJUqSGyGqkNZerRngP4BX2r2Ck96yi1jWcrVjRLvaAHyy6Rzq0O9A7wBhe2H3lxaNRQhRuUlyI0QVoigKH3X5iGH1h6HX6i3e/tM96mGn13L0Wjx/hNtA/49NJ7a+D+FHLB6PEKJykuRGiCrM0ruveDhYM76TPwCf/nkWQ7OHoeEg0+aaKx+HjBSLxiOEqJwkuRGiCjKqRtacX8PQNUMtPvbm8S51cLHTcT4qiZUh1+HeL8HBG6LPmgYYy/RwIUQxSXIjRBWkoLDk9BIuxF9g9rHZFm3byUbH/7rWBWDGX+dI1zvDfbNB0cCRhXDYsjO5hBCVjyQ3QlRBiqIwOWgyAMvPLrf4ujePBvvh5WTN9bhUftkfBnW6Qo83TCfXvwjXD1s0HiFE5SLJjRBVVLvq7QiuEUyWMYtvQr6xaNs2Oi3P9jRty/DNlvMkp2dBp+cgcCAYMmDZaEiOsWhMQojKQ5IbIaqwZ1s9C8D6S+s5HXvaom0Pa10LP3c7YpIzmLfrEigKDJ0FbnUg/iqsfEzWvxFCmKXMk5uZM2fi7++PjY0NQUFB7Ny5M9+yK1eupHfv3nh6euLk5ESHDh34448/LBitEJVLI/dG9PPrB8CXhy271oxOq2FKn0AAZu+4yK3kDLBxhuELwcrWtHP4tg8tGpMQonIo0+Rm6dKlTJ48mddee42QkBA6d+5M//79CQsLy7P8jh076N27N+vXr+fQoUN0796de++9l5CQEAtHLkTl8UzLZ7BSrNh1fRdnYs9YtO17mlanUXUnEtOzmLX9gumgV2MY9JXp/Y6P4cxGi8YkhKj4FNXSC138S7t27WjVqhWzZs3KPtawYUOGDBnC9OnTC1VH48aNGT58OG+++WahyickJODs7Ex8fDxOTpZdoVWI8mpR6CIaujWklVcri7e99XQUYxf8jbWVhm0vdqO6s63pxO8vwN8/mHpzntgObv4Wj00IUX4U5fd3mfXcZGRkcOjQIfr06ZPjeJ8+fdizZ0+h6jAajSQmJuLm5pZvmfT0dBISEnK8hBA5jWw4skwSG4BugZ609XMjPcvIRxv+Ne6n7wdQqw2kxcOyUZCZWibxCSEqnjJLbqKjozEYDHh5eeU47uXlRWRkZKHq+Oyzz0hOTubBBx/Mt8z06dNxdnbOfvn4+BQrbiEqu5spN8kyZlmsPUVReOOeRigKrD4Szt4Lt2dJWelh2I9g5wGRx+H352WBPyFEoZT5gGJFUXJ8rapqrmN5Wbx4MW+//TZLly6lWrVq+ZabOnUq8fHx2a+rVy27nocQFcm8E/MYsHIAv134zaLtNq3lzMjbm2q+ueYEmQaj6YRzTXhg3u0F/hbBoQUWjUsIUTGVWXLj4eGBVqvN1UsTFRWVqzfnv5YuXcr48eNZtmwZvXr1KrCstbU1Tk5OOV5CiLxpFS1phjS+PfItaVlpFm37xT4NcLPXcy4qifm7L/1zok5X6Hl7TN2Gl+Dq3xaNSwhR8ZRZcqPX6wkKCmLTpk05jm/atIng4OB8r1u8eDFjxozhl19+YeDAgaUdphBVykMNHsLb3psbKTdYcnqJRdt2ttPxSv8GgGlbhoj4f42x6TgZGtxjWuBv8XCIvWjR2IQQFUuZPpaaMmUKc+bMYd68eYSGhvLcc88RFhbGhAkTANMjpdGjR2eXX7x4MaNHj+azzz6jffv2REZGEhkZSXx8fFndghCVirXWmonNJwLww/EfSMiw7AD8B1rVIsjXlZQMA+//HvrPCUWBod9D9RaQEgMLH5AVjIUQ+SrT5Gb48OHMmDGDadOm0aJFC3bs2MH69evx9fUFICIiIseaN99//z1ZWVk89dRTVK9ePfv17LPPltUtCFHpDKo7iLrOdUnISGD+ifkWbVujUZg2uDEaBdYdi2D3+X/tWG7tACOWgXNtiL0Aix+SGVRCiDyV6To3ZUHWuRHi7jaHbWby1snYaG34/b7fqWaX/6D90vD22pMs2HOZOp72bHy2C3qrf/0/7OYZmNvbNEW84b2mGVUarUXjE0JYXoVY50YIUX718OlBc8/mqKicjD5p8faf610fDwdrLt5MZs6u/4yv8QyEhxaDVg+hv8Gfr1s8PiFE+SbJjRAiF0VReCf4HVYNWkX32t0t3r6zrY5XB5gGF3+9+TzX4/7z+MmvIwy5vbL5vpmwd6aFIxRClGeS3Agh8lTXpS4+TmW36OXQljVp6+dGaqaB99adyl2g6QPQ6x3T+z9ehVNrLBugEKLckuRGCHFXIVEhLApdZNE2FUVh2pDGaDUKG05Esv3szdyFOj4LrccDKqx8AsL2WzRGIUT5JMmNEKJA526dY/SG0Xzy9yecjj199wtKUANvJ8YG+wHw1poTpGcZchZQFOj/MdTvB1lpphlUMRcsGqMQovyR5EYIUaAA1wD6+PbBoBp4e8/bGIyGu19Ugp7tFUA1R2sux6Tww448Fu/TWpm2aKjRElJjYeH9kBRl0RiFEOWLJDdCiLt6pe0rOOocORlzkiVnLLtysaONjtfvaQTAN1vPczU2JXchvb1pDRyX2nDrEvw4CJLyeIwlhKgSJLkRQtyVp50nk4MmA/DV4a+ITI4s+IISdm+z6gTXdSct08jzy45iMOaxPJdDNRi1Ghyrw81Q+PFeSXCEqKIkuRFCFMoD9R+ghWcLUrJS+GD/BxZtW1EUPrq/GQ7WVhy4HMvMrefzLuheFx5dJwmOEFWcJDdCiELRKBre7PAmVooVW69u5e9Iy+7O7eNmx7tDGgMwY/M5Dl25lXdBj3qS4AhRxUlyI4QotADXACa1msR7Hd+jtVdri7c/tGUtBreogcGoMnlpCIlpmXkXlARHiCpNkhshRJGMbTKWwfUGoyhKmbT/7pAm1HK15WpsKm+tKWBrCElwhKiyJLkRQpgtMSORKwlXLNqmk42OGcNboFFgZch11hy5nn9hSXCEqJIkuRFCmOXYzWMMXj2YKdumkGnM5/FQKWnt58akngEAvL7qRN7Tw++QBEeIKkeSGyGEWWo51iLTmMnZW2eZdWSWxdt/uns9Wvu6kpiexeSlR8gyGPMvnGeCIwv9CVFZSXIjhDCLm40br7V/DYAfjv/Ajms7LNq+lVbDF8Nb4GhtxaErt/gmv+nhd/w3wZnTE26esUywQgiLkuRGCGG2fn79eCjwIQCm7pxKeFK4Rdv3cbPjvaFNAPhq8zkOXYkt+AKPejDmd3CrA3FhMLc3XLJsUiaEKH2S3AghiuXFNi/SxL0JCRkJvLD9BTIMGRZtf3CLmgxtWROjCs8uOUJCftPD73CvC+P/Ap92kBYPP98HRy27pYQQonRJciOEKBa9Vs9n3T7DSe/E8ejjzD0x1+IxTBvcGB83W67dSuXN1SfufoG9O4xeC42HgjETVj0J2z4CNY9tHYQQFY4kN0KIYqvhUIPpnafT168voxqOsnj7jjY6ZgxviVajsPpIOL/sD7v7RTobuH8edJxs+nrbB7DmKciybM+TEKLkKapatf6rkpCQgLOzM/Hx8Tg5OZV1OEKIEvTt1vN88scZtBqFOaNb071BtcJdeHAe/P4CqAbw7wIP/gy2LqUaqxCiaIry+1t6boQQJU5VVVafX01KZgHrz5SCid3qcn+rWhiMKk/9cpjj1+ILd2HrcTBiKegdTAOM5/U1DTgWQlRIktwIIUrcB/s/4I3db/DevvewZOewoih8eH9TOgd4kJJhYOyCvwte4O/fAnrD2A23p4qfhjm9IGx/6QYshCgVktwIIUpcH78+aBQNv138jRXnVli0bZ1Ww8yRrWhY3YnopHQenX+AuJRCjqOp3gwe2wxeTSDpBszvD7tmgLGABQKFEOWOJDdCiBLXxrsNz7R8BoDp+6cTGhNq0fYdbXQsGNuGGs42XLyZzGM/HiQt01C4i51rwriN0OQB0xicv96CXx6E5JjSDVoIUWIkuRFClIpxTcbRtVZXMowZTNk2hfj0Qo5/KSFeTjYsGNcWRxsrDl65xZRlRzAaC/mIzNoR7p8D934JVjZwfhN81wmu7CndoIUQJUKSGyFEqdAoGt7v9D417GtwLekaEzdPtPgA4/pejswe1Rq9VsP645G8v74IPUiKAkFjTI+p3AMgMRwWDIQdn8hjKiHKOUluhBClxtnama97fo2T3olTMac4EV2IBfZKWIe67nwyrBkAc3ddYt6uS0WrwLsJPLENmj0EqhG2vAcL75OdxYUox2SdGyFEqTsRfYK49Dg61exUZjHM2naBjzaeRlFg5ohW9G9aveiVhCyC35+HrFRw8DY9uvLvXPLBCiFykXVuhBDlShOPJjkSm5spN8kyZlk0hgld6zCqvS+qCs8uPcKuc9FFr6TlSFMvjmcDSIqEnwbBn69DhmUftwkhCibJjRDCoq4kXGHE+hG8sfsNjKrlxq4oisLbgxrTu5EXGVlGxi44wNqjZuxiXq0BPL4VWj5ieky152uY2R4ubCn5oIUQZpHkRghhUVcSrhCdEs26i+t4d9+7Fl3kT6tR+Prhlgxo6k2mQWXS4pCij8EB0NvB4G9hxDJwqgVxV+DnobDySZkyLkQ5IMmNEMKiutTqwvTO09EoGn49+yufHPzEogmOjU7L1w+34tEOvgBMW3eK6RtCzYuhfl94ah+0mwAocGwJfNsGji6VHcaFKEOS3AghLK6ffz/e7vA2AD+f+plvjnxj0fa1GtMjqhf7BgLw/faLPL/8KJkGMx6TWTtC/4/gsb+gWmNIiYFVT8DC++HW5ZINXAhRKJLcCCHKxNCAobza7lUAZh+bzZzjcyzavqIoPNW9Hp880AytRmHl4euM//EgyelmDnSu1Rqe3A493gCtNVzYDDM7mMbkGCw7eFqIqk6SGyFEmXm4wcM8F/QcADuu7bDoAOM7hrX24YfRQdjoNOw4e5MRP+wjJindvMq0OujyAvxvD/h1hswU02yqme3h1Fp5VCWEhcg6N0KIMvfzqZ/p798fD1uPMoshJOwW4xb8za2UTPzc7fhpXDtqu9uZX6GqQsjP8NfbpkdVADVbQ+93wK/s1vsRoqIqyu9vSW6EEOXOTyd/olOtTtRxrmPRdi/cTGL03ANcj0vFw8Ga7x5pRWs/t+JVmpZgejS191vITDYdC+gDPd8yrX4shCgUWcRPCFFh/XH5Dz45+Akjfx/J9qvbLdp2XU8HVk4MpoG3I9FJ6Tz4/V4++eM0GVnFeFxm4wQ9XoNJIdDmMdBYwbk/TRtxrpoAcWEldwNCCECSGyFEORPkFUSraq1IykzimS3P8MOxHyw6VdzLyYblEzpwX6uaGFX4dusF7pu1m/NRicWr2NELBn4GTx2AxkMBFY4uhq+DYOOrkHijROIXQshjqbIORwiRh0xDJh8e+JBlZ5cB0Nu3N+91fA87XTHGwJhh/fEIXl11nLiUTKytNLw2sCGj2vuiKErxK79+CDa9BZd3mr7W6qHZg9DhaajWsPj1C1HJyJibAkhyI0TFsfzscj7Y/wFZxizqu9bny+5fUsuxlkVjuJGQxgvLj7Lz9l5UXet78skDzajmZFP8ylXVNGV820dw7cA/x+v1MiU5dbpBSSRSQlQCktwUQJIbISqWwzcOM2XbFGLSYvis62f08etj8RiMRpWf9l5m+obTpGcZcbXTMf2+pvRrYsbO4vm5esA08Pj0OtOeVQBeTaHDU9DkfrDSl1xbQlRAktwUQJIbISqeyORIdlzbwYOBD2YfS8xIxFHvaNE4zt1IZPLSI5wMTwDggaBavD6wIS52JZh4xF6CfbMgZOE/s6scq0PbJ6DlKHDwLLm2hKhAJLkpgCQ3QlR8EUkR3L/2fu4LuI//tfgf9jp7i7WdkWVkxl9nmbX9AqoKjjZWTOhal7Ed/bDTW5VcQ6m34OB82P89JEWajmmsoF5vaPEw1O8HVtYl154Q5ZwkNwWQ5EaIim/BiQV8dugzAKrZVuOFNi/Qz69fyQz0LaS/L8fyxuoTnI40zaLycLDm2Z71GN6mNnqrEpyImpUBJ1bA3z+YBiHfYetqelzVfATUbCVjc0SlJ8lNASS5EaJy2HFtBx8e+JCriVcBaOvdllfbvUpdl7oWi8FoVPntWDif/XmWsNgUAGq72fF8n/rc26wGGk0JJxw3z8CRX+DYUkiM+Oe4R31o/jA0Gw7ONUu2TSHKCUluCiDJjRCVR7ohnfkn5jPn+BzSDelYKVb8r8X/eKLZExaNIyPLyNK/w/hqy3luJpr2pWrg7chL/QLpHlit5HuUjAa4uM20Tk7oOshKvX1CAZ+2pkdWgf3Bs4H06IhKQ5KbAkhyI0Tlcy3xGh/9/RHbrm5jYvOJ/K/F/wBQVdWij6pSMrKYv/sy322/QGKaaSfwNn6ujO9Uh54Nq6HTlsK6qWkJcGo1HFkMYXtynnP1g/r9TYmOb7BpY08hKihJbgogyY0QldeOazto5N4oewPOLWFb+PHkj4xsOJIetXtgpSnBAb8FiEvJYNb2CyzYfZn021s3eDhYc39QTYa39qGOp0PpNBx/Dc5uhDMb4dJ2MGT8c87aGQJ6mZKdOl3BoVrpxCBEKZHkpgCS3AhRdTzx5xPsjdgLgLe9N8MDh3N/wP242rhapP3I+DTm77nEikPXiU5Kzz7e1s+N4W18GNC0OrZ6bek0np4EF7fCmQ1w9g9Iic553j3A1Jvj18n0p7NlF0cUoqgkuSmAJDdCVB03km+w7Owyfj37K7FpsQBYa60Z4D+A/v796VCjg0XiyDQY2XI6imV/X2XrmSiMt//VdbS2YnDLGgxvXZsmNZ1K7xGa0QDXDsLZDXBuE9w4Cfznn36X2uB7O9HxDQa3OjJeR5QrktwUQJIbIaqedEM6f1z+g4WnFhIaGwpAO+92zOk7J7tMfHo8ztbOpR5LZHwavx66ytKDV7kam5p9vIazDZ0DPOlS35OO9dxLdmHA/0qJhbB9cGW36RVx9J9Vke+wcYHqzaB6c/BubnrvXg80pdTTJMRdSHJTAEluhKi6VFXl6M2jrD6/muaezRkaMBSA6NRoei7vSSO3RnTx6UK3Wt1o4NagVAcjG40q+y7GsOTvq2w8GUlG1j/JhUaBZrVc6BLgQZf6nrTwccGqNAYj35GeCFf3w+XdcGWPaT0dY2bucjo78GpiSniqNzPNxnKvB3ZupRebELdJclMASW6EEP+1+cpmJm+bnOOYp60nzT2b08i9Eb18e+Hv7F9q7admGDhwOZYdZ2+y89xNzt5IynHe0dqKDnXdaVHbhUbVnWhUw4lqjiWwcWd+stIhKhQij5l6dSKOQeTxf005/w9bN/AIMCU67vX+ee9WR1ZRFiVGkpsCSHIjhMjLzZSb7Ly+k+1Xt7M3Yi+p//pF/mnXT+nr1xeA07Gn2XV9F43cGtHQvWGpDE6OiE9l57lodpy9ya7z0cSl5O5F8XCwplENp+xkp1F1J/w97NGW9MKBdxgNEHP+drJz1JTsxJyHhOsFXKSY9sVyrpXPy8e00rKM7RGFUKGSm5kzZ/LJJ58QERFB48aNmTFjBp07d863/Pbt25kyZQonT56kRo0avPTSS0yYMKHQ7UlyI4S4m3RDOsduHuNUzClOxpxkcqvJ1HCoAcCc43P48vCX2WWd9E7UdKhJLcda1HKsxUOBD2WXLQkGo8qJ6/HsvRjDqfAETkUkcPFmUvag5H+z0WnwcbWjlqstNV1tqXXnvYvpvYeDvuQftWUkQ8wFiDkH0edNCc+d9xmJd79eZ2ealm7vCfbVTBuD5nh/+5ytK9i6yFo9VViFSW6WLl3KqFGjmDlzJh07duT7779nzpw5nDp1itq1a+cqf+nSJZo0acLjjz/Ok08+ye7du5k4cSKLFy/m/vvvL1SbktwIIYpja9hWNlzawMmYk4QlhuU6v3LQSgJcAwDTHlgLQxfibuuOq40rbtZuuNm4md7buNHdpzsuNi4ApGWloSgKes3dE5DUDAOnI02Jzp2E53REIqmZhgKvs9FpqOFii4eDNW52elzt9bjZ63Czt8bNXoernR43ez2udnocbayw01uZv0+WqkLyTYi/alp/J/v1r6+Tbxa9Xp29Kcmxccn9p7UjWDuA3h70jqY///u1ztb0srKRHqMKpsIkN+3ataNVq1bMmjUr+1jDhg0ZMmQI06dPz1X+5ZdfZu3atYSGhmYfmzBhAkePHmXv3r2FalOSGyFESUnJTOF60nWuJ13nWuI1riVdY1LLSdjp7AB4Z+87/Hr213yvXzN4DXVc6gDwdcjXzD42GyvFCludLfY6e+yt7LHT2WGns+PN9m9S28n0n76d13ay6/ourLXW6LQ6rLXWWCk6ktMgMU3FS9uaW4nWXLuVwsW4S0SkXuZWShaqqgVVAyigKoCCIa0GGG0BULSJKLo403kUUEGn0WKj12Krs8Je44mDtQMO1lp0ugwUqwT0Wg16Kw06jRadlYK11gqdlQZnvRsOOnt0WgUD6aQb49FoFKw0ClqNBiuNgk7NxCYjGg9DJi4ZSejTYzCmRpKSFolVxi106XFYpcVilR6HVaapF8jZYMT+9q+tdAVitPnP3nIyGHG4XTYDiLbKWdaosUa1skbVWuOgtcZBa4uqtSHTSsdNrRWqxgpVowOt/vZ7PWh12FnZ4qi1A42WTEVDNJmoihYUK1SNFjSmPxVFg63WFicre9BoyELhpiEFFOV2eS0oGlRFAUWDrdYGZ50jKBoMqNzMSABFY/pxoJiuw/S1tcYaF72T6XoVojJMSx2od352kJ286TV6XK2ds49HpUWbfvzwT9nbf1pp9bjpXbKPRaVFo6L+JxG8XVaxwt36n8eyN9NjMKpGtIqGanbV8PIJyPdnY46i/P62zHKdecjIyODQoUO88sorOY736dOHPXv25HnN3r176dOnT45jffv2Ze7cuWRmZqLT5e6uTE9PJz39n8WzEhISSiB6IYQAO50dAa4B2T01/zWp5STuq3cfsWmxxKbFciv9FrGpt/9Mi8Xd1j27bEqmaePNLDWLxIxEEv/zSCdLzcp+f/TmUX45/Uu+cS0ZuITGHoEAzD1+mBmH52Obz4Smzg5voM2owa3kDK5kHibOblmuMkYgGbgZNhZDlKleK+e/sa2xArKA9FyXkHptJFmJTU1lHY9iW2txvvGmhg8jKz4IcEHrkI6dz2+gAWxvv3C6/QK7yO54xdXDWUnGaHeZYz55/74AeDhaw7D4TOyUNK5YZ/JkTZd8y068dYv/xVwB4IpOx/21qv9z0nD7ddvYuASm3IoD4KqVlgd88t+s9OH4RF6NvQVAtEbDQ775L5Y4JDGJd6NNSUqyotDbzyffsn2Tkvn0Zgxg+vn09s/9tOOOrimpfHPjn16ye3xrka7Ju0eubWoacyOjsr8eVrsm8fkkkE3T0vkl4kb216N9ahBpZUXtzEwWXEuDty/nG1NpK7PkJjo6GoPBgJeXV47jXl5eREZG5nlNZGRknuWzsrKIjo6mevXqua6ZPn0677zzTskFLoQQheRq41roAccvtH6BiS0mkpKZQnJWsunPzOTsr73s/vm3L8griMebPk66IZ0MQwaZxszs9wbVgJP1P/+rrWZXjVbVWmFQDRiMBrLULIyqMfv1VKfGNPZoDMDq8xHMOrILg2pARcWo3n4ZjRhVlecGNMHHtjkpGVnsj7rJnzccABVVxfS/e0zT7QHa1/HERa2Owahyw3CVC+o/s6ZyPi5QqensgK2tI6oKKVb2xKl6chb+54osa08iHRoRARj0rmD8O9/v6WJlAIs1d8ZwhqEaZ/LfB1HK7bpXGnty0NgIazJJV2PQGDf/q0zOeI+rAcw1emGFgSRjKlbGSznO//u6SNWTnWptNBhJJgudMTnfeJNUB86o9mhQSceIPq+BVbdlqnoiVTcUVFTUAsuqqpY41T77a50Kaj7lFVVDkmqb/bVeBet8ympVhWTVOvv7Y3W7rJWqkKGU4jpNhVBmj6XCw8OpWbMme/bsoUOHf1YJff/99/n55585ffp0rmvq16/P2LFjmTp1avax3bt306lTJyIiIvD29s51TV49Nz4+PvJYSgghhKhAKsRjKQ8PD7Raba5emqioqFy9M3d4e3vnWd7Kygp3d/c8r7G2tsbaWtZZEEIIIaqKUlzysmB6vZ6goCA2bdqU4/imTZsIDg7O85oOHTrkKv/nn3/SunXrPMfbCCGEEKLqKbPkBmDKlCnMmTOHefPmERoaynPPPUdYWFj2ujVTp05l9OjR2eUnTJjAlStXmDJlCqGhocybN4+5c+fywgsvlNUtCCGEEKKcKbPHUgDDhw8nJiaGadOmERERQZMmTVi/fj2+vr4AREREEBb2zzoS/v7+rF+/nueee45vv/2WGjVq8NVXXxV6jRshhBBCVH5lvkKxpck6N0IIIUTFU5Tf32X6WEoIIYQQoqRJciOEEEKISkWSGyGEEEJUKpLcCCGEEKJSkeRGCCGEEJWKJDdCCCGEqFQkuRFCCCFEpSLJjRBCCCEqFUluhBBCCFGplOn2C2XhzoLMCQkJZRyJEEIIIQrrzu/twmysUOWSm8TERAB8fHzKOBIhhBBCFFViYiLOzs4Flqlye0sZjUbCw8NxdHREUZQSrTshIQEfHx+uXr1aKfetquz3B5X/HuX+Kr7Kfo9yfxVfad2jqqokJiZSo0YNNJqCR9VUuZ4bjUZDrVq1SrUNJyenSvuXFir//UHlv0e5v4qvst+j3F/FVxr3eLcemztkQLEQQgghKhVJboQQQghRqUhyU4Ksra156623sLa2LutQSkVlvz+o/Pco91fxVfZ7lPur+MrDPVa5AcVCCCGEqNyk50YIIYQQlYokN0IIIYSoVCS5EUIIIUSlIsmNEEIIISoVSW4KMHPmTPz9/bGxsSEoKIidO3cWWH779u0EBQVhY2NDnTp1+O6773KVWbFiBY0aNcLa2ppGjRqxatWq0gq/UIpyjytXrqR37954enri5OREhw4d+OOPP3KUWbBgAYqi5HqlpaWV9q3kqSj3t23btjxjP336dI5y5elnWJT7GzNmTJ7317hx4+wy5ennt2PHDu69915q1KiBoiisXr36rtdUtM9gUe+xon0Gi3p/FfEzWNR7rEifw+nTp9OmTRscHR2pVq0aQ4YM4cyZM3e9rjx8DiW5ycfSpUuZPHkyr732GiEhIXTu3Jn+/fsTFhaWZ/lLly4xYMAAOnfuTEhICK+++iqTJk1ixYoV2WX27t3L8OHDGTVqFEePHmXUqFE8+OCD7N+/31K3lUNR73HHjh307t2b9evXc+jQIbp37869995LSEhIjnJOTk5ERETkeNnY2FjilnIo6v3dcebMmRyxBwQEZJ8rTz/Dot7fl19+meO+rl69ipubG8OGDctRrrz8/JKTk2nevDnffPNNocpXxM9gUe+xon0Gi3p/d1SUzyAU/R4r0udw+/btPPXUU+zbt49NmzaRlZVFnz59SE5OzveacvM5VEWe2rZtq06YMCHHsQYNGqivvPJKnuVfeukltUGDBjmOPfnkk2r79u2zv37wwQfVfv365SjTt29f9aGHHiqhqIumqPeYl0aNGqnvvPNO9tfz589XnZ2dSyrEYinq/W3dulUF1Fu3buVbZ3n6GRb357dq1SpVURT18uXL2cfK08/v3wB11apVBZapiJ/BfyvMPealPH8G/60w91fRPoP/Zc7PsCJ9DqOiolRA3b59e75lysvnUHpu8pCRkcGhQ4fo06dPjuN9+vRhz549eV6zd+/eXOX79u3LwYMHyczMLLBMfnWWJnPu8b+MRiOJiYm4ubnlOJ6UlISvry+1atXinnvuyfW/Sksozv21bNmS6tWr07NnT7Zu3ZrjXHn5GZbEz2/u3Ln06tULX1/fHMfLw8/PHBXtM1gSyvNnsDgqwmewpFSkz2F8fDxArr9v/1ZePoeS3OQhOjoag8GAl5dXjuNeXl5ERkbmeU1kZGSe5bOysoiOji6wTH51liZz7vG/PvvsM5KTk3nwwQezjzVo0IAFCxawdu1aFi9ejI2NDR07duTcuXMlGv/dmHN/1atXZ/bs2axYsYKVK1cSGBhIz5492bFjR3aZ8vIzLO7PLyIigg0bNvDYY4/lOF5efn7mqGifwZJQnj+D5qhIn8GSUJE+h6qqMmXKFDp16kSTJk3yLVdePodVblfwolAUJcfXqqrmOna38v89XtQ6S5u58SxevJi3336bNWvWUK1atezj7du3p3379tlfd+zYkVatWvH111/z1VdflVzghVSU+wsMDCQwMDD76w4dOnD16lU+/fRTunTpYladpc3cWBYsWICLiwtDhgzJcby8/fyKqiJ+Bs1VUT6DRVERP4PFUZE+h08//TTHjh1j165ddy1bHj6H0nOTBw8PD7Raba4sMioqKle2eYe3t3ee5a2srHB3dy+wTH51liZz7vGOpUuXMn78eJYtW0avXr0KLKvRaGjTpo3F/8dRnPv7t/bt2+eIvbz8DItzf6qqMm/ePEaNGoVery+wbFn9/MxR0T6DxVERPoMlpbx+BourIn0On3nmGdauXcvWrVupVatWgWXLy+dQkps86PV6goKC2LRpU47jmzZtIjg4OM9rOnTokKv8n3/+SevWrdHpdAWWya/O0mTOPYLpf4tjxozhl19+YeDAgXdtR1VVjhw5QvXq1Ysdc1GYe3//FRISkiP28vIzLM79bd++nfPnzzN+/Pi7tlNWPz9zVLTPoLkqymewpJTXz2BxVYTPoaqqPP3006xcuZItW7bg7+9/12vKzeewxIYmVzJLlixRdTqdOnfuXPXUqVPq5MmTVXt7++wR7a+88oo6atSo7PIXL15U7ezs1Oeee049deqUOnfuXFWn06m//vprdpndu3erWq1W/fDDD9XQ0FD1ww8/VK2srNR9+/ZZ/P5Utej3+Msvv6hWVlbqt99+q0ZERGS/4uLissu8/fbb6saNG9ULFy6oISEh6tixY1UrKyt1//795f7+vvjiC3XVqlXq2bNn1RMnTqivvPKKCqgrVqzILlOefoZFvb87HnnkEbVdu3Z51lmefn6JiYlqSEiIGhISogLq559/roaEhKhXrlxRVbVyfAaLeo8V7TNY1PuraJ9BVS36Pd5RET6H//vf/1RnZ2d127ZtOf6+paSkZJcpr59DSW4K8O2336q+vr6qXq9XW7VqlWP626OPPqp27do1R/lt27apLVu2VPV6vern56fOmjUrV53Lly9XAwMDVZ1OpzZo0CDHh7YsFOUeu3btqgK5Xo8++mh2mcmTJ6u1a9dW9Xq96unpqfbp00fds2ePBe8op6Lc30cffaTWrVtXtbGxUV1dXdVOnTqpv//+e646y9PPsKh/R+Pi4lRbW1t19uzZedZXnn5+d6YF5/f3rTJ8Bot6jxXtM1jU+6uIn0Fz/p5WlM9hXvcFqPPnz88uU14/h8rtGxBCCCGEqBRkzI0QQgghKhVJboQQQghRqUhyI4QQQohKRZIbIYQQQlQqktwIIYQQolKR5EYIIYQQlYokN0IIIYSoVCS5EUIIIUSlIsmNEMIi3n77bVq0aFFm7b/xxhs88cQTpVZ/VFQUnp6eXL9+vdTaEEIUjqxQLIQoNkVRCjz/6KOP8s0335Cenp69M7Al3bhxg4CAAI4dO4afn1+ptTNlyhQSEhKYM2dOqbUhhLg7SW6EEMUWGRmZ/X7p0qW8+eabnDlzJvuYra0tzs7OZREaAB988AHbt2/njz/+KNV2jh8/Ttu2bQkPD8fV1bVU2xJC5E8eSwkhis3b2zv75ezsjKIouY7997HUmDFjGDJkCB988AFeXl64uLjwzjvvkJWVxYsvvoibmxu1atVi3rx5Odq6fv06w4cPx9XVFXd3dwYPHszly5cLjG/JkiUMGjQox7Fu3brxzDPPMHnyZFxdXfHy8mL27NkkJyczduxYHB0dqVu3Lhs2bMi+5tatW4wcORJPT09sbW0JCAhg/vz52eebNm2Kt7c3q1atMv+bKYQoNkluhBBlZsuWLYSHh7Njxw4+//xz3n77be655x5cXV3Zv38/EyZMYMKECVy9ehWAlJQUunfvjoODAzt27GDXrl04ODjQr18/MjIy8mzj1q1bnDhxgtatW+c69+OPP+Lh4cGBAwd45pln+N///sewYcMIDg7m8OHD9O3bl1GjRpGSkgKYxu2cOnWKDRs2EBoayqxZs/Dw8MhRZ9u2bdm5c2cJf6eEEEUhyY0Qosy4ubnx1VdfERgYyLhx4wgMDCQlJYVXX32VgIAApk6dil6vZ/fu3YCpB0aj0TBnzhyaNm1Kw4YNmT9/PmFhYWzbti3PNq5cuYKqqtSoUSPXuebNm/P6669nt2Vra4uHhwePP/44AQEBvPnmm8TExHDs2DEAwsLCaNmyJa1bt8bPz49evXpx77335qizZs2ad+1JEkKULquyDkAIUXU1btwYjeaf/2N5eXnRpEmT7K+1Wi3u7u5ERUUBcOjQIc6fP4+jo2OOetLS0rhw4UKebaSmpgJgY2OT61yzZs1ytdW0adMc8QDZ7f/vf//j/vvv5/Dhw/Tp04chQ4YQHByco05bW9vsnh4hRNmQ5EYIUWZ0Ol2OrxVFyfOY0WgEwGg0EhQUxKJFi3LV5enpmWcbdx4b3bp1K1eZu7V/ZxbYnfb79+/PlStX+P333/nrr7/o2bMnTz31FJ9++mn2NbGxsfnGIoSwDHksJYSoMFq1asW5c+eoVq0a9erVy/HKbzZW3bp1cXJy4tSpUyUSg6enJ2PGjGHhwoXMmDGD2bNn5zh/4sQJWrZsWSJtCSHMI8mNEKLCGDlyJB4eHgwePJidO3dy6dIltm/fzrPPPsu1a9fyvEaj0dCrVy927dpV7PbffPNN1qxZw/nz5zl58iTr1q2jYcOG2edTUlI4dOgQffr0KXZbQgjzSXIjhKgw7Ozs2LFjB7Vr1+a+++6jYcOGjBs3jtTUVJycnPK97oknnmDJkiXZj5fMpdfrmTp1Ks2aNaNLly5otVqWLFmSfX7NmjXUrl2bzp07F6sdIUTxyCJ+QohKT1VV2rdvz+TJk3n44YdLrZ22bdsyefJkRowYUWptCCHuTnpuhBCVnqIozJ49m6ysrFJrIyoqigceeKBUkychROFIz40QQgghKhXpuRFCCCFEpSLJjRBCCCEqFUluhBBCCFGpSHIjhBBCiEpFkhshhBBCVCqS3AghhBCiUpHkRgghhBCViiQ3QgghhKhUJLkRQgghRKXyf5ZE3S7yclroAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.plot(ts, lcce.real, label='CCE')\n", "plt.plot(ts, lmecce.real, label='ME-CCE')\n", "plt.plot(ts, (lmecce_1 * lcce).real, ls='--', label='ME-CCE1 * CCE2')\n", "\n", "plt.xlabel('Time (ms)')\n", "plt.ylabel('Coherence')\n", "\n", "plt.legend();" ] }, { "cell_type": "markdown", "id": "4b3fd3bd-aafa-4e1c-8d78-e9a42a9d2e53", "metadata": {}, "source": [ "We find that the exact accounting for both coherent and incoherent processes leads to a singificantly different qualitative behaviour of the coherence." ] }, { "cell_type": "markdown", "id": "feba43a4-94b0-4d6b-84bb-b3a4e43c8ab0", "metadata": {}, "source": [ "## Dissipation in the central spin\n", "To account for the central spin coupling to its own Markovian bath we need to use gCCE flavor of the ME approach.\n", "\n", "We set up the dissipation on the central spin in the similar way as the bath spins, by calling `.add_single_jump` method of the `calc.center` object. As an example, consider pure dephasing of the central spin due to the Markovian environment." ] }, { "cell_type": "code", "execution_count": 10, "id": "88e22290", "metadata": { "tags": [] }, "outputs": [], "source": [ "t2 = 1 # in ms \n", "decay_rate = 2 / t2 # in rad / ms\n", "calc.bath['e'].so.clear() # Remove bath spin dissipators\n", "calc.center.add_single_jump('z', rate=decay_rate) # z for operator Sz \n", "calc.center[0].detuning = 1e6 # Detune central spin from the spin bath" ] }, { "cell_type": "code", "execution_count": 11, "id": "a5406ea5-d3cd-47b1-910c-32d2228b3b95", "metadata": { "tags": [] }, "outputs": [], "source": [ "calc.order = 2\n", "lmegcce = calc.compute(ts, method='megcce', pulses=1)" ] }, { "cell_type": "markdown", "id": "d345d9a3-fd8f-4edc-be5b-54742047c51f", "metadata": {}, "source": [ "And now we show that central spin dissipators are trivial to deal with and can be factorized out. We can divide the obtained coherence curve by the expected markovian dephasing, and recover the coherence curve limited purely by the interactions with the bath." ] }, { "cell_type": "code", "execution_count": 12, "id": "6f4ab85f-bf03-44b2-964c-f2356f069ae7", "metadata": { "tags": [] }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAGwCAYAAABVdURTAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABymklEQVR4nO3dd3QUVRvH8e/sbipphIRAIEDovYUaukBoglgQRemIiIiACCIviChiQ7EAilIsCKhUkWJAelMg9NBLKAkQSA8pm533j4VoTALps7t5PufMYTJ7d+Y3SYZ9MuVeRVVVFSGEEEIIG6HTOoAQQgghREGS4kYIIYQQNkWKGyGEEELYFCluhBBCCGFTpLgRQgghhE2R4kYIIYQQNkWKGyGEEELYFIPWAYqayWTi+vXruLq6oiiK1nGEEEIIkQOqqhIXF4evry863YPPzRS74ub69ev4+flpHUMIIYQQeXDlyhXKly//wDbFrrhxdXUFzN8cNzc3jdMIIYQQIidiY2Px8/NL/xx/kGJX3Ny/FOXm5ibFjRBCCGFlcnJLidxQLIQQQgibIsWNEEIIIWyKFDdCCCGEsCnF7p4bIYSwZCaTiZSUFK1jCKEJe3v7hz7mnRNS3AghhIVISUnh4sWLmEwmraMIoQmdToe/vz/29vb5Wo8UN0IIYQFUVSU8PBy9Xo+fn1+B/PUqhDW538lueHg4FSpUyFdHu1LcCCGEBTAajSQmJuLr64uzs7PWcYTQhLe3N9evX8doNGJnZ5fn9cifBkIIYQHS0tIA8n06Xghrdv/3//7xkFdS3AghhAWRMe9EcVZQv/9S3AghhBDCpmha3OzYsYOePXvi6+uLoiisXr36oe/Zvn07AQEBODo6UrlyZb766qvCDyqEEEIIq6FpcZOQkECDBg348ssvc9T+4sWLdO/enTZt2hASEsKbb77J6NGjWbFiRSEnFUIIIYS10PRpqW7dutGtW7cct//qq6+oUKECs2fPBqBWrVocOHCAjz/+mCeffLKQUuZMmkkl/Mo57kSfJs27DoqioFNAh+7evIKD3gFXezfsDTqc7fXY6eWqoBDCNkRERDBjxgx+//13rl27RunSpWnYsCFjxoyhY8eOAISEhPDee++xY8cOYmJiqFChAu3ateP111+nevXqXLp0CX9//yzXv3fvXlq0aFGUuySsmFU9Cr53716CgoIyLOvSpQsLFiwgNTU1y8fGkpOTSU5OTv86Nja2ULLFnNmFbllfnvMvlW2b1Jj6JF3vd+8rFZcaU1FUOxTVET0O6BUHDIojdjpHvO2q0sn3OSp6lqBiKWe83RV8XGQUcyGE5bl06RKtWrXCw8ODDz/8kPr165OamsqmTZt4+eWXOXXqFOvWrePJJ5+kS5cuLFmyhCpVqnDz5k1++eUXpkyZwvLly9PXt3nzZurUqZNhG6VKZf9/qxD/ZVXFTUREBD4+PhmW+fj4YDQaiYyMpGzZspneM3PmTN5+++1Cz5bmVYs4nHEwmTCix4gOFAD1XguVewvMFCOKLhVIRSURI2AE7pdht26ncujo6fT3mgshexwojYddGco4l6NmqWp0rtyUpuWrSIdfQtgYVVW5m5q/x2HzyslOn6unVkaOHImiKPz111+UKFEifXmdOnUYMmQIiYmJDB48mO7du7Nq1ar01/39/WnevDnR0dEZ1leqVCnKlCmT7/0QxZdVFTeQ+TExVVWzXH7fpEmTGDduXPrXsbGx+Pn5FXguby8vvAct5sDiRwEV+v0M1btk2TbFaCI+OYWw2IZEJcYRlRRPTFI8MckJxCbHE5eSAKU8MZYux+XbCVyOvkHKvUIomYvcMF3kRjwciYfll0FJrEeA4xga+HnQ0M+DyqX1VCgpf+UIYc3upqZRe+omTbZ9cnoXnO1z9vFw584dNm7cyIwZMzIUNvd5eHiwatUqIiMjmTBhQpbr8PDwyE9cITKxquKmTJkyREREZFh28+ZNDAZDtqcsHRwccHBwKIp4UKk1tBgJ++bA2ldg5D5w9szUzN6gw9PgiGeJSjle9Y34jhy8dpbjNy5wLvoy1+KvcCPpAkm6qyQnerH18i22nr6Foo/Hpfq76IzeeNtXo2bJ2rSp0Jgu1Rrj4ZT5Px4hhMiPc+fOoaoqNWvWzLbN2bNnAR7Y5t8CAwMznY2OiYlBr9fnPagoVqyquGnZsiW//fZbhmV//PEHTZo0yVc3zQWq4xQ4FwyRZ2D9eHhqYYGs1sfFne41mtC9RpMMy+OS73L0WiRnw40cuRrDXxH7iQVMhlvcMN3ixu09bL8N7xzSYW/ypW6Jx+ju350Gfh5UK+2CQW5qFsIiOdnpOTk967O/RbHtnHrY2fN/t8mp5cuXU6tWrQzLpLARuaFpcRMfH8+5c+fSv7548SKHDx/G09OTChUqMGnSJK5du8b3338PwIgRI/jyyy8ZN24cL7zwAnv37mXBggUsXbpUq13IzM4JHv8Kvu0Mx1dAzUeh7hOFtjlXBydaVfajVeX7SxoSFv0Ev5/5m31XQ7gQG0q06QKKPo5U/VX2XLjBjpBjADi5XMPVdz2+ztWo6lGFhmVq0LpCHSp5li60vEKInFEUJceXhrRUrVo1FEUhNDSU3r17Z9mmevXqAJw6dYqWLVs+dJ1+fn5UrVq1IGOKYkZRc1tSF6Bt27bRoUOHTMsHDhzI4sWLGTRoEJcuXWLbtm3pr23fvp2xY8dy4sQJfH19mThxIiNGjMjxNmNjY3F3dycmJgY3t0J8+ujPGbDjQ3Aqab485ardzXEmk4kjEZcJPvc38bHlOBeu5+jVGJKdt+NY5rfMb0hzpYTiS2O3Z2jp24SqpV2p7uNCKZciurwnRDGUlJTExYsX8ff3x9HRUes4udKtWzeOHTvG6dOnM913Ex0djZ2dHZUqVaJ169YZbij+dxsPD4/0R8FDQkJo2LBhEaUXluRBx0FuPr81LW60UGTFjTEFvu0IEUehWhfotxwsaMwYk0nlrysX2HBuFydvnyY88RKxaddQDXfS2yReepG0u+Y+JwzuB3EuHYyHoRKVXGvQuExdOlVuTO3S5eVJLSEKgDUXNxcvXiQwMBBPT0+mT59O/fr1MRqNBAcHM2/ePEJDQ1mzZg19+vSha9eujB49mqpVqxIZGcnPP/9MWFgYy5YtSy9usnoU3MPDw+q+LyL3pLjJoyIrbgBunIT57SAtBXp9CY37F+72CsCt+Fh2h53k4PVT2CXX4/ItE2duxnHTbjn2nnsyvyHNBVddJR7xeoFHqtahTTVv6ZxQiDyw5uIGIDw8nBkzZrBu3TrCw8Px9vYmICCAsWPH0r59ewAOHDjAzJkz2blzZ/qTq4888givv/46VatWfWAnfkuXLuWZZ54pwj0SWpDiJo+KtLgB2P0ZBE8Fe1d4aTeUrFj42ywEEXFRBJ8LYf+1Y5yJPsWt5Auk6iNQFBMA8Wcmo6a54lnCnlZ1EuhYswy9ajaTszpC5JC1FzdCFAQpbvKoyIsbUxos6g5X9kGlNjBgLdjIB3703QS2XDjC/qsnUOKas+XUTSLjk3GqMB9DiQvojWVo4tWJUc360rBsJa3jCmHRpLgRQoqbPCvy4gbg9nn4qjWkJkLX96HFS0Wz3SJmTDOx/ewN3vt7KjeMB1B0RgBUVcFVrUknv+6MbvEE3jKMhBCZSHEjRMEVN7ZxCsHSlaoCQe+Y5zdPg8izmsYpLAa9jo41y7Kl/zf88eQWepQdjbOpGoqiEq8LZfW1WTzyU1+W/xWW634vhBBCiJyS4qaoNBkKVR4BYxKsehHSjFonKlS+bp68H/QC+wevZGHHVQS4PYPO6EXijSAmrjxGv2/2czEyQeuYQgghbJAUN0VFUcxPTDm6w7WDsOsTrRMVmablq7L48cn8PTCYCW174WinY++F23Rb+BFDV39AYmryw1cihBBC5JAUN0XJvRx0/9g8v/0DuB6ibZ4iZm8w8GK7Kvwxph0tqtlj8FrHXzE/0ur7Xvx6bLfW8YQQQtgIKW6KWr0+ULs3mIywcjik3tU6UZGrUMqZnwZ35MmKr0CaM0bDdaYdfIney17nRnyM1vGEEEJYOSluipqiwKOfgksZ8+Cam9/WOpEmdDodb3ccyJreaymra4WiqJxP3kjnnx/li71rtI4nhBDCiklxowVnT3jsS/P8/nlwYZumcbRU2dOHP/p/xSu1P0RnLIWqj+br01OYuXmr1tGEEEJYKSlutFKtMzQZYp5fPRLuRmsaR2vDm3Zj27PrqGwfRMrNrny9OZH31odiMskj40JYskGDBqEoSpYDGI8cORJFURg0aFCGtv+dunbtWiBZVqxYQfv27XF3d8fFxYX69eszffp07tz5Z8y8lJQUPvzwQxo0aICzszNeXl60atWKRYsWkZqaWiQ5ReGT4kZLQe+CZ2WIvQYbJmidRnMlnV1Y8+wsxjUbDsD8HRd49Zfd8jSVEBbOz8+PZcuWcffuP/cQJiUlsXTpUipUqJChbdeuXQkPD88wLV26NN8ZJk+eTN++fWnatCkbNmzg+PHjzJo1iyNHjvDDDz8A5sKmS5cuvP/++wwfPpw9e/bw119/8fLLL/PFF19w4sSJQs8pioZB6wDFmn0JeHw+LAyCo8uhRneo01vrVJob0b4qpVwceGP132yJfoegJV789vTXlHR20TqaECILjRs35sKFC6xcuZLnnnsOgJUrV+Ln50flypUztHVwcKBMmTK5Wv+ePXsYOXIkp06dom7duvzvf//j8ccfJyQkhIYNG/LXX3/x3nvvMXv2bF599dX091WqVInOnTsTHR0NwOzZs9mxYwcHDhygUaNG6e0qV65Mnz59SElJyVdOYTnkzI3W/JpC63Hm+XVjIC5C0ziWok8TP97o5YbO/hYxylGClj1HWPQtrWMJUXRUFVIStJny0IP44MGDWbRoUfrXCxcuZMiQIfn+NsTFxdGzZ0/q1avHoUOHeOedd5g4cWKGNkuWLMHFxYWRI0dmuQ4PD4/0dp06dcpQ2NxnZ2dHiRIl8p1XWAY5c2MJ2k2Es39AxFFYMwqe+8X8VFUx90KzTrg6zmbGgfEk6S/w2Ip+fNfjG+qXqaR1NCEKX2oivOerzbbfvG4+s5wL/fv3Z9KkSVy6dAlFUdi9ezfLli1j27ZtGdqtW7cOF5eMZ2EnTpzIlClTslzvkiVLUBSFb775BkdHR2rXrs21a9d44YUX0tucPXuWypUrY2dn98CMZ8+epX379jnan9zmFJZFihtLYLCHJ76Br9vCuWA4uOifm42LuWfqt6WU03xe2zkKoyGC53/vz2cd5tKhcj2towkh/sXLy4sePXrw3XffoaoqPXr0wMvLK1O7Dh06MG/evAzLPD09ARgxYgQ//vhj+vL4+HhOnz5N/fr1Mwyi2KxZswzvV1UVJQd/EOa03cNyCssnxY2lKF0TOr0Fm96ETZPBv515wE1B52oNWeT8HUM3DifNcJPR24YxLfFTnqwbqHU0IQqPnbP5DIpW286DIUOGMGrUKADmzJmTZZsSJUpQtWrVLF+bPn0648ePz7Asq4LkvwPvVq9enV27dpGamvrAszfVq1cnNDT0ofvxsJzC8sk9N5ak+UtQqY35dPSqEWBK0zqRxQgoV4WVj/+EQ1olTCYDU1Ze4thV6c1Y2DBFMV8a0mLK42Xxrl27kpKSkv5UUm6VLl2aqlWrpk8ANWvW5OjRoyQn//PU5IEDBzK8r1+/fsTHxzN37tws13v/huJ+/fqxefNmQkIyD31jNBpJSJDBfG2FFDeWRKeD3vPAwQ2u/gV7Ptc6kUWp7OnDhqeXUFN9ncREd4Z+9zfXo4vf8BVCWCq9Xk9oaCihoaHo9fos2yQnJxMREZFhioyMzHad/fr1w2QyMXz4cEJDQ9m0aRMff2weo+/+GZ3mzZszYcIEXnvtNSZMmMDevXu5fPkyW7ZsoU+fPnz33XcAjBkzhlatWtGxY0fmzJnDkSNHuHDhAj///DPNmzfn7Nmzec4pLIsUN5bGww+6vm+e3/oe3Djx4PbFjLeLG4v696C6jws345J5/oefZTwqISyIm5sbbm5u2b6+ceNGypYtm2Fq3br1A9f322+/cfjwYRo2bMjkyZOZOnUqQIb7cD744AN++ukn9u/fT5cuXahTpw7jxo2jfv36DBw4EDA/3h0cHMyECRP4+uuvadGiBU2bNuXzzz9n9OjR1K1bN885hWVR1P9evLRxsbGxuLu7ExMT88ADUFOqCkufhTMboEx9GLbFfNOxSHc1KpGeC+eS6vkj7kpttjy3GEc7+R4J65WUlMTFixfx9/fP8KEtMluyZAmDBw8mJiYGJycnreOIAvSg4yA3n99y5sYSKQr0/AycPM2Ph+/4SOtEFqd8SWemdG0N6IhVjtH310mYTCatYwkhCsH333/Prl27uHjxIqtXr2bixIk8/fTTUtiIbElxY6lcfeDRT8zzO2fBtYPa5rFAT9RpwcCqk1BVhQspfzB6/RdaRxJCFIKIiAief/55atWqxdixY+nTpw/z58/XOpawYFLcWLI6j0PdJ0FNg1UvQarcPPtfr7fpQytP8/X0bZELmL1ntbaBhBAFbsKECVy6dCn9ksWnn36Ks3PeHlcXxYMUN5au+8fg4gORp+HPd7VOY5HmPTqWCnaPoCgq355+h99C/9Y6khBCCA1JcWPpnD2h173LLXvnwKXd2uaxQDqdjl/6fIirqTaKLoWpfy4kPEbOcgkhRHElxY01qN4FGvUHVFj9EiTHa53I4jjbOfDrk/NwTXiSqCvdGbL4APHJRq1jCSGE0IAUN9aiy3vg7gfRlyFYBm7Liq+bJ8ufmYCXiyOh4bGMXhqCMU2eoBJCiOJGihtr4egGj90bq+XAQji3Wds8FsrP05lvBjTBwT6ZPTFzGLX+U60jCSGEKGJS3FiTyu2g2Yvm+TWvwN0obfNYqEYVSvJk6xjsPA6yK/IHgs8e1jqSEEKIIiTFjbXpNA08q0DcddgwUes0FuvdTgNxVxug6NKYuGMSCf8adE8IIYqDKVOmMHz4cK1jpDt27Bjly5cvkgFKpbixNvbO8PhXoOjg6HII/U3rRBZJp9PxbY/3Ic2ZVMNVXvjtPa0jCWGTBg0ahKIojBgxItNrI0eORFEUBg0alKn9f6euXbsWSJ4VK1bQvn173N3dcXFxoX79+kyfPp07d+6kt0lJSeHDDz+kQYMGODs74+XlRatWrVi0aBGpqan5zjlt2jSeeeaZAtmfvLpx4wafffYZb7755gPbVapUidmzZ2f52rZt2yhbtmyW34d/T4MGDeLSpUsMHToUf39/nJycqFKlCm+99RYpKSnp66tXrx7NmjXj008L/3YBKW6skV8zaPWqef63MZAgI9VmpaZ3efpXfQ2Ao/GrWH50p8aJhLBNfn5+LFu2jLt3/+mCISkpiaVLl1KhQoVM7bt27Up4eHiGaenSpfnOMXnyZPr27UvTpk3ZsGEDx48fZ9asWRw5coQffvgBMBc2Xbp04f3332f48OHs2bOHv/76i5dffpkvvviCEyf+Gaw4rznXrl3LY489lu/9yY8FCxbQsmVLKlWqlOd1rF27ll69emXY/9mzZ+Pm5pZh2WeffcapU6cwmUx8/fXXnDhxgk8//ZSvvvoqU3E1ePBg5s2bR1paWj738CHUYiYmJkYF1JiYGK2j5E9qkqrOaaGqb7mp6rLnVNVk0jqRxer0/Ytq3cV11frfdlAjE2K1jiNElu7evauePHlSvXv3rtZRcmXgwIHqY489ptarV0/98ccf05cvWbJErVevnvrYY4+pAwcOzNQ+t3bv3q02aNBAdXBwUAMCAtRVq1apgBoSEqKqqqru379fBdTZs2dn+f6oqChVVVX1gw8+UHU6nXro0KFMbVJSUtT4+Ph85QwLC1Pt7OzSt5eVhQsXqjVr1lQdHBzUGjVqqHPmzEl/bfDgwWq9evXUpKSk9EyNGzdW+/Xrp6qqql68eFEF1KVLl6otW7ZUHRwc1Nq1a6tbt27NsI169eqpX3755QOztmvXTgUyTP9WpUoVdd26dRmWLVq0SHV3d3/Id8Hsww8/VP39/TMsS05OVh0cHNQtW7Zk+Z4HHQe5+fyWMzfWyuBgvjylM5gvTR37RetEFmtRr/dQ0jwwcpe3NmzTOo4QuZKYmpjtlJyWnOO2ScakHLXNq8GDB7No0aL0rxcuXMiQIUPyvL5/i4uLo2fPntSrV49Dhw7xzjvvMHFixnsOlyxZgouLCyNHjsxyHR4eHuntOnXqRKNGjTK1sbOzo0SJEvnKunbtWtq2bZu+vf/65ptvmDx5MjNmzCA0NJT33nuPKVOm8N133wHw+eefk5CQwBtvvAGY75uJjIxk7ty5Gdbz+uuv89prrxESEkJgYCC9evXi9u3bAERFRXH8+HGaNGnywKwrV66kfPnyTJ8+Pf0szH0nTpwgIiKCjh075vVbQUxMDJ6enhmW2dvb06BBA3buLNwz6YZCXbsoXGUbQLs3YOu7sH48VGoNbr5ap7I45d09eaPRB0z+9Qrr0nT0qXeLdtW9tY4lRI40/6l5tq+1KdeGuZ3++dBr/3N77hqz7p27iU8TFnX9p/jouqIrUcmZn7g8NvBYnnL279+fSZMmcenSJRRFYffu3Sxbtoxt27Zlartu3TpcXFwyLJs4cSJTpmTdh9eSJUtQFIVvvvkGR0dHateuzbVr13jhhRfS25w9e5bKlStjZ2f3wJxnz56lffv2Odqn3OYEWLNmzQMvSb3zzjvMmjWLJ554AgB/f39OnjzJ119/zcCBA3FxceHHH3+kXbt2uLq6MmvWLLZs2YK7u3uG9YwaNYonn3wSgHnz5rFx40YWLFjAhAkTuHz5Mqqq4uv74M8DT09P9Ho9rq6ulClTJtN+dOnSBUdHxweuIzvnz5/niy++YNasWZleK1euHJcuXcrTenNKihtr13osnF4P1w/B2lfguV9BUbROZXH6NQrkzJUTLN5zidd/OcIfY9vi4WyvdSwhbIaXlxc9evTgu+++Q1VVevTogZeXV5ZtO3TowLx58zIsu/8X/ogRI/jxxx/Tl8fHx3P69Gnq16+f4YO2WbNmGd6vqipKDv7vy2m7h+XMSmxsLNu3b+ebb77J8vVbt25x5coVhg4dmqEwMxqNGYqXli1bMn78+PQzVG3bts20rpYtW6bPGwwGmjRpQmhoKED6vU///n4tWbKEF198Mf3rDRs20KZNm2z3Zc2aNdmeBXuY69ev07VrV/r06cOwYcMyve7k5ERiYt7PEuaEFDfWTm8wX576qo25Y79D30HAIK1TWaSJXWuy4+wtwpL30HfFb2zq/5nWkYR4qP399mf7ml6nz/D1tqe3ZdtWp2S8C2HjkxvzlSsrQ4YMYdSoUQDMmTMn23YlSpSgatWqWb42ffp0xo8fn2FZVgWJqqoZvq5evTq7du0iNTX1gWdvqlevnl4EPMyDcmZlw4YN1KpVi4oVK2b5uslk7jH9m2++oXnzjGfk9Hp9hna7d+9Gr9dz9uzZHG///vfoflEZFRWFt7f5LHWvXr0ybLNcuXLZriciIoJDhw7Ro0ePHG/7vuvXr9OhQwdatmzJ/Pnzs2xz584dqlSpkut154bcc2MLvGtAx6nm+U2TIeqSpnEslZO9nomPeuHou5zrpj95d9tPWkcS4qGc7ZyznRz0Djlu62hwzFHb/OjatSspKSnpTyTlRenSpalatWr6BFCzZk2OHj1K8r/6qzpw4ECG9/Xr14/4+PhM96bcFx0dnd5u8+bNhISEZGpjNBrz1QfLmjVr6NWrV7av+/j4UK5cOS5cuJBhH6tWrYq/v396u48++ojQ0FC2b9/Opk2bMtzLdN++ffsy5D548CA1a9YEoEqVKri5uXHy5Mn0Nq6urhm25+TkBJjvgfnvk0tr166lZcuW2Z55y861a9do3749jRs3ZtGiReh0WZcYx48fz/Kep4IkxY2taPESVAiElHhY/TKYZEylrHSpUZdGrk8BsPzCbI5HXNE4kRC2Q6/XExoaSmhoaIYzEf+VnJxMREREhikyMvsuLfr164fJZGL48OGEhoayadMmPv74Y+CfsxXNmzdnwoQJvPbaa0yYMIG9e/dy+fJltmzZQp8+fdJv2B0zZgytWrWiY8eOzJkzhyNHjnDhwgV+/vlnmjdvnuFMSW5yGo1GNmzY8NBHwKdNm8bMmTP57LPPOHPmDMeOHWPRokV88sknABw+fJipU6eyYMECWrVqxWeffcarr77KhQsXMqxnzpw5rFq1ilOnTvHyyy8TFRWVfgO3TqejU6dO7Nq164FZwNzPzY4dO7h27Vr6vuXlUfbr16/Tvn17/Pz8+Pjjj7l161b69+zfLl26xLVr1+jUqVOu1p9rOXqey4bYzKPgWbl9QVXfLWt+PHzvXK3TWKyElCS18YLuat3FddWWi55R09LStI4khNU/Cp6drB4F5z+PHwNqjRo1Hrid3bt3q/Xr11ft7e3VgIAA9aefflIB9dSpUxnaLV++XG3btq3q6uqqlihRQq1fv746ffr0DI9mJyUlqTNnzlTr1aunOjo6qp6enmqrVq3UxYsXq6mpqXnKuXnzZrV8+fIP/mbds2TJErVhw4aqvb29WrJkSbVt27bqypUr1bt376q1a9dWhw8fnqH9448/rgYGBqpGozH9UfCffvpJbd68uWpvb6/WqlUr06PVGzduVMuVK/fQ/9/27t2r1q9fX3VwcFABNT4+XnV0dFTPnDmTZfvsHgVftGhRlt+v/5YZ7733ntqlS5ds8xTUo+CKqv7nwqWNi42Nxd3dnZiYGNzc3LSOU/AOLIR1Y8HgCCN2gVc1rRNZpK0XjvHK9gEoOiN9KkxkaofntY4kirmkpCQuXryIv79/np9QKU6WLFnC4MGDiYmJSb/EoqXRo0djNBqzvSxWUC5duoS/vz8hISE0bNgw23aqqtKiRQvGjBnDs88+m+P1r1y5kv/9738ZLmkVlOTkZKpVq8bSpUtp1apVlm0edBzk5vNbLkvZmoDBUOURMCbBqhGQZtQ6kUXqULkeTTzMl6d+vTCf6LuFP9aJECLvvv/+e3bt2sXFixdZvXo1EydO5Omnn7aIwgagbt26vPTSS1rHSKcoCvPnz8dozN1ngIuLCx988EGhZLp8+TKTJ0/OtrApSPK0lK1RFOj1JcxtCdcOwO7Z0Hb8Q99WHH3a5VXaLduIaojinT9XMqtHf60jCSGyERERwdSpU4mIiKBs2bL06dOHGTNmaB0rnSUNUHlfgwYNaNCgQa7eExQUVEhpzE+qVa9evdDW/29yWcpWHVkGq14EnR0M2wy+DbVOZJE+2rGGr7ZexTGtCltfb09pV7kcILQhl6WEkMtS4mHq94VaPcGUCiuGQYpcdsnKa617UbdUAxJS0vjkjzNaxxFCCFEApLixVYoCPT8HV1+4fdbc/43IRKdTmPJobQB+PnyMP88V/E10QuRGMTuZLkQGBfX7L8WNLXP2NPdejAIHF0HoOq0TWaQmlTxpVvcqzpVn8ebOaem9iApRlO73C5OSkqJxEiG0c//3/0H9JOWE3FBs6yq3g1ajYfdn5rGnygWAW1mtU1mcCe07Mij4KxJ0oXy+by1jAntrHUkUMwaDAWdnZ27duoWdnV22vbsKYatMJhO3bt3C2dkZgyF/5YncUFwcGFNgQScIPwKV28Pzq0D+48zkmV+mcCJxNXpjaXY/v54SDg4Pf5MQBSglJYWLFy/K2UNRbOl0Ovz9/bG3zzywcW4+v6W4KS4iz5oH1zTehaB3IfAVrRNZnIi4KDr/0h308bQrNYwvH31V60iiGDKZTHJpShRb9vb22Z61zM3nt1yWKi68qkHXmbBuDGx+G/zbQdn6WqeyKGVcS9K13CA2RnzJ9ptLCIvuRwUPb61jiWJGp9PJo+BC5JNcmyhOAgZBzUfvPR4+FFIStU5kcWZ0HIrB6Av6u4ze+KHWcYQQQuSBFDfFiaJAry/AtSxEnoE//qd1IotjbzAwsv44VJOBU9dSOX8rXutIQgghckmKm+LG2RN6zzPPH1gAp9Zrm8cCvdC0C42Uj0i61ZmZ60O1jiOEECKXpLgpjqp0+OeG4rWjIC5C2zwW6K3uLTDoFDaH3mT3uUit4wghhMgFzYubuXPnpo8hERAQwM6dOx/YfsmSJTRo0ABnZ2fKli3L4MGDuX37dhGltSGPTIEy9SDxNqx+CeTR0wyqlnbh+RYV0TmGMe7PSaTkcmRdIYQQ2tG0uFm+fDljxoxh8uTJhISE0KZNG7p160ZYWFiW7Xft2sWAAQMYOnQoJ06c4JdffuHvv/9m2LBhRZzcBhgc4MmFYHCC83/CrllaJ7I4w9uVp0SFxSQ67OP9ncu0jiOEECKHNC1uPvnkE4YOHcqwYcOoVasWs2fPxs/Pj3nz5mXZft++fVSqVInRo0fj7+9P69atefHFFzlw4EC220hOTiY2NjbDJO7xrg497hU1W9+DC9s0jWNpfN3daerZG4BVF77DmJambSAhhBA5ollxk5KSwsGDBwkKCsqwPCgoiD179mT5nsDAQK5evcr69etRVZUbN27w66+/0qNHj2y3M3PmTNzd3dMnPz+/At0Pq9foOWjUH1STefTw2OtaJ7Io7z7yIqQ5YTRE8PHuX7WOI4QQIgc0K24iIyNJS0vDx8cnw3IfHx8iIrK+wTUwMJAlS5bQt29f7O3tKVOmDB4eHnzxxRfZbmfSpEnExMSkT1euXCnQ/bAJ3T8Cn3qQcAt+HQJpqVonshi+bp40dH8UgJ/PLZZu8YUQwgpofkOxoigZvlZVNdOy+06ePMno0aOZOnUqBw8eZOPGjVy8eJERI0Zku34HBwfc3NwyTOI/7Jzg6e/AwQ3C9sKWt7VOZFFmdHwJ1eRAqv4qn+1do3UcIYQQD6FZcePl5YVer890lubmzZuZzubcN3PmTFq1asXrr79O/fr16dKlC3PnzmXhwoWEh4cXRWzbVaoK9J5rnt/zBYT+pm0eC1LBw5u6Lt0A+PH0Ajl7I4QQFk6z4sbe3p6AgACCg4MzLA8ODiYwMDDL9yQmJmYaUEuv1wPmMz4in2r1hJajzPOrR8KdC9rmsSDvPvIyanIZYm8GsO3MDa3jCCGEeABNL0uNGzeOb7/9loULFxIaGsrYsWMJCwtLv8w0adIkBgwYkN6+Z8+erFy5knnz5nHhwgV2797N6NGjadasGb6+vlrthm3pNA38WkByLPw8AFLvap3IIlQtVYany84mNboZc7ZelGJaCCEsmKbFTd++fZk9ezbTp0+nYcOG7Nixg/Xr11OxYkUAwsPDM/R5M2jQID755BO+/PJL6tatS58+fahRowYrV67Uahdsj94O+iwCZy+IOAYbJmidyGKMaFcFe4OOg5ej2HteOo4UQghLpajF7E/Q2NhY3N3diYmJkZuLH+TCNvi+N6DCY3PNj4wL/rf6MMtPraJM6Qh2Dvla6zhCCFFs5ObzW/OnpYSFqtweOkw2z//+GkQc1zSOpXiquQcOZdYQrd/DT0e2aR1HCCFEFqS4Edlr8xpU7QzGu+b7b5JitE6kuYZlK1HRvi0AX4Z8pXEaIYQQWZHiRmRPp4Mn5oO7H9w5D6tkgE2Aqa1fQVV1xCknWHEi6960hRBCaEeKG/Fgzp7w9Pegd4DTv8OuT7ROpLnmFapRztAKgNkHsh4HTQghhHakuBEPV64x9PjYPP/nu3Bus7Z5LMD/Wo1CVRWiOczvp7MfuFUIIUTRk+JG5EzjAdB4IKCaB9iMuqx1Ik218a+Nj745AB/tm6NxGiGEEP8mxY3Iue4fgW9juBsFy58v9h38TWo5CmN8Na5cCuDMjTit4wghhLhHihuRcwYH6PvDvQ7+jpofES9e3SRl0KlqA9q5TSbtrj9f/nlO6zhCCCHukeJG5I57eXhqISg6OLwEDizUOpGmRj1SFYB1R69zMTJB4zRCCCFAihuRF5XbmcegAtgwEa78pWkcLdUt507rmnbYea9n3OYZWscRQgiBFDcirwJHQ+3HwJRq7uAv/qbWiTTTqZ4O+1I7OHt3IxfvyIjhQgihNSluRN4oCjw2B7xqQFw4/DIY0oxap9LEgEaPYJ9WAUWXyvTt32odRwghij0pbkTeObjCM0vA3hUu74LNb2mdSBM6nY7H/J8F4EDUOmKSEjVOJIQQxZsUNyJ/vKrB4/d66d37JZxcq20ejYxv3QfF6AH6eN7f8ZPWcYQQoliT4kbkX62e5ntwANaMgqhLmsbRgrOdA4GlewOw8cpyTDIGlxBCaEaKG1EwOk6F8s0gOQZ+HQLGFK0TFbmp7QajmhwwGiKY99fvWscRQohiS4obUTD0dvDUAnD0gGsHYcvbWicqcr5untQt8Sgpd1qx7bgcWkIIoRX5H1gUHI8K0HuueX7vl3B6g7Z5NDAraCLGW734+xycuB6jdRwhhCiWpLgRBatmD2gx0jy/agREX9E2TxEr5+FE93plAViw86LGaYQQoniS4kYUvE5vmwfYTIqGFUMhLVXrREXqhTb+6J0usfHWBxyPCNM6jhBCFDtS3IiCZ7CHPovAwR2u7Ic/39U6UZGqX94Drwpb0Lse552d87WOI4QQxY4UN6JwlKwEj31hnt89G84Ga5mmyD1TvT8AJ+M2cSNe7r0RQoiiJMWNKDy1H4OmL5jnV70Isde1zVOEXmnRC72xNOiTeGfbYq3jCCFEsSLFjShcQe9CmfqQeBt+HVpsxp8y6PV0KtcHgB03VpKUWvz6/RFCCK1IcSMKl50j9FlsHn8qbA9sf1/rREXmf+36Q1oJVMMdPt2zUus4QghRbEhxIwpfqSrQc7Z5fsfHcGm3pnGKiodTCRq6dwdgxfklMiSDEEIUESluRNGo9xQ0eh5QYfVLkByndaIiMbXdC5iSfYi5VZ8Dl29rHUcIIYoFKW5E0ekyE9wrQPRl+GOK1mmKRDWvsjzqOYvUqEAW7LqsdRwhhCgWpLgRRcfRDXrPMc8fXARnN2ubp4i80LYyAH+cvMGlyASN0wghhO2T4kYULf+20Pwl8/zaUXA3Sts8RaBqaVfa1fBA73aAyX/O1TqOEELYPCluRNHrOBVKVYW4cFg/Qes0RaJV3VicfH/lSMJyIuJsv6ATQggtSXEjip69Mzz+NSg6OPYznFyjdaJCNyygC3pjGRRdMu9s/07rOEIIYdOkuBHaKN8EWo8zz68bC/E3tc1TyHQ6HZ3KPQnArhurSTEWj84MhRBCC1LcCO20mwg+9cy9F/82BlRV60SFanLb5yHNCZPhNnP2/6Z1HCGEsFlS3AjtGOzh8a9AZwenf4cjy7ROVKhKOrtQy7UzAMvP/KRxGiGEsF1S3AhtlakLHSaZ5zdMgJir2uYpZJNaDUNVdSToTrHpbIjWcYQQwiZJcSO0F/gqlG8KybGw5mWw4WEKGvn6460LwBhflVWHrmgdRwghbJIUN0J7egP0/goMTnBhGxxYoHWiQjWz9fvcvTKMrcfsuJMgo4ULIURBk+JGWAavqtD5bfN88FS4c1HbPIWoub839cq5k2w0sfSvMK3jCCGEzZHiRliOpi9ApTaQmgi/j7PZp6cURWFwq0oohlgWHPuWxNRkrSMJIYRNkeJGWA6dDnp+BnoHOP8nHP1Z60SFpls9H1z855Hq/juf7F6hdRwhhLApUtwIy1KqCrSfaJ7fNAkSbmubp5A42dnRuGQQAKsv2PYj8EIIUdSkuBGWJ3A0lK5j7tzvj8lapyk0b7YZgmrSk6y/yK/HdmsdRwghbIYUN8Ly6O2g1+eAAkeWwvmtWicqFDW8fSlnFwjAvMOLtQ0jhBA2RIobYZnKN4Fmw83z68ZASqKmcQrLqIDBANxI+4vjEfLklBBCFAQpboTl6jgF3MpB1CXY/oHWaQpFz1pNcUqrhqKYeG/XQq3jCCGETZDiRlguB1foMcs8v+cLCD+qbZ5C8lS1Z1BNdhy/Hk1SaprWcYQQwupJcSMsW41uUPsxUNPgt9Fgsr0P/1dbPo7bzbeJvd6ZNYevaR1HCCGsnhQ3wvJ1+xAc3OF6COz/Wus0Bc7BYMeglrUBWLT7EqqNdl4ohBBFRYobYflcy/wzNMOf70K07d1427dJBZzs9JyJPsmvx/7WOo4QQlg1KW6EdWg8ECoEQmoC/P6azQ3N4O5sR706BynhP4cvQ+ZqHUcIIayaFDfCOqQPzWAPZ/+AEyu1TlTghjd5FIDb6iFCrtvuwKFCCFHYpLgR1sO7OrQZb57fMBHuRmmbp4B1qtqAEqaaKIrKB7vlsXAhhMgrKW6EdWk9BrxqQMIt8/03Nuapan0BOBH3BzFJttlxoRBCFDYpboR1MTj80/fN3wvMT1DZkFHNe6EYS4I+kQ93yoCaQgiRF1LcCOvj3wbq9QFU+H08mExaJyowjnb2NPfqCcDGK79isqF9E0KIoqJ5cTN37lz8/f1xdHQkICCAnTt3PrB9cnIykydPpmLFijg4OFClShUWLpT7E4qdoHfB3hWuHYCQH7ROU6AmtRmAajKQZIpl67kLWscRQgiro2lxs3z5csaMGcPkyZMJCQmhTZs2dOvWjbCw7Psxefrpp9myZQsLFizg9OnTLF26lJo1axZhamERXMtAhzfN85unQeIdTeMUpMqePgQ6v0XCuQmsORindRwhhLA6iqphd6jNmzencePGzJs3L31ZrVq16N27NzNnzszUfuPGjTzzzDNcuHABT0/PHG0jOTmZ5OTk9K9jY2Px8/MjJiYGNze3/O+E0E6aEb5uCzdPmPvB6fW51okKzPFrMTz6xS4MOoU9bzxCaTdHrSMJIYSmYmNjcXd3z9Hnd57P3Jw/f57//e9/PPvss9y8eRMwFx8nTpzI0ftTUlI4ePAgQUFBGZYHBQWxZ8+eLN+zdu1amjRpwocffki5cuWoXr0648eP5+7du9luZ+bMmbi7u6dPfn5+OdxDYfH0BujxsXn+0Pdw9YC2eQpQ3XLuBFQsidFkZO7urI8HIYQQWctTcbN9+3bq1avH/v37WblyJfHx8QAcPXqUt956K0friIyMJC0tDR8fnwzLfXx8iIiIyPI9Fy5cYNeuXRw/fpxVq1Yxe/Zsfv31V15++eVstzNp0iRiYmLSpytXruRwL4VVqBgIDZ7FfHPxOJsaWLNbY4USVT9gRfj/SPjX2UchhBAPlqfi5o033uDdd98lODgYe3v79OUdOnRg7969uVqXoigZvlZVNdOy+0wmE4qisGTJEpo1a0b37t355JNPWLx4cbZnbxwcHHBzc8swCRvTebp5YM3wI3BwkdZpCsyzDRuj06mgj+OTPb9qHUcIIaxGnoqbY8eO8fjjj2da7u3tze3bt3O0Di8vL/R6faazNDdv3sx0Nue+smXLUq5cOdzd3dOX1apVC1VVuXr1ai72QNgUl9LwyP/M81umQ/wtbfMUkBIODjTy6A7A2ou/aJxGCCGsR56KGw8PD8LDwzMtDwkJoVy5cjlah729PQEBAQQHB2dYHhwcTGBgYJbvadWqFdevX0+/DAZw5swZdDod5cuXz8UeCJvTdCiUqQ9JMeanp2zEpNaDUFUdSfrz/Bb6l9ZxhBDCKuSpuOnXrx8TJ04kIiICRVEwmUzs3r2b8ePHM2DAgByvZ9y4cXz77bcsXLiQ0NBQxo4dS1hYGCNGjADM98v8e339+vWjVKlSDB48mJMnT7Jjxw5ef/11hgwZgpOTU152RdgKnR56fGKeP/wjhO3TNk8BqVW6PGX0zQCYc/A7jdMIIYR1yFNxM2PGDCpUqEC5cuWIj4+ndu3atG3blsDAQP73v//leD19+/Zl9uzZTJ8+nYYNG7Jjxw7Wr19PxYoVAQgPD8/Q542LiwvBwcFER0fTpEkTnnvuOXr27Mnnn9vOI8AiH/yaQqP+5vnfXzM/Km4DXmhgLvCvpu7hwp0bGqcRQgjLl69+bi5cuMChQ4cwmUw0atSIatWqFWS2QpGb5+SFFUq4DV80hqRo6PoBtBihdaJ8M5lMNF3ckxR9GB1KvcLnjw7XOpIQQhS53Hx+G/KzocqVK1O5cuX8rEKIglWiFHR6C9aNha0zoM7j4Jr1DerWQqfT8VzVV/jyzzAOXK+CsZsJg17zkVOEEMJi5el/yKeeeor3338/0/KPPvqIPn365DuUEPnSeCD4NoLkWPjzHa3TFIiRLTvjoffnekwSm0Nvah1HCCEsWp478evRo0em5V27dmXHjh35DiVEvuj00O1D83zIj3D9sKZxCoKjnZ5nmpp71164J1TjNEIIYdnyVNzEx8dn6LzvPjs7O2JjY/MdSoh882sG9foAKmx8A7QbQq3APNe8Ao5lVnHS7jWCzx7WOo4QQlisPBU3devWZfny5ZmWL1u2jNq1a+c7lBAFotM0MDhB2F44sUrrNPlWrqQzZUoaUXSpfHHge63jCCGExcrTDcVTpkzhySef5Pz58zzyyCMAbNmyhaVLl/LLL9KTqrAQ7uWh9VjY9h4ET4Ua3cDOuvtD6l/nWWafOMiFpB3ciI/Bx8X94W8SQohiJk9nbnr16sXq1as5d+4cI0eO5LXXXuPq1ats3ryZ3r17F3BEIfIh8BVwKw8xV2DPF1qnybfBjTujN5ZG0SXz4a6ftI4jhBAWKV/93Fgj6eemGDq+An4dAnbOMOoAuOdsiBBLNWrdbLbfXoDB6MvBwRvQ6eSxcCGE7cvN53e+/ldMSUnh6tWrhIWFZZiEsCh1noAKLSE10SbGnZrYpj+qyQ6j4TpLj8nTiUII8V95Km7Onj1LmzZtcHJyomLFivj7++Pv70+lSpXw9/cv6IxC5I+iQNeZgALHfoYr1j0ApZ97KfzsWwGw6OhSjdMIIYTlydMNxYMGDcJgMLBu3TrKli2LoigFnUuIguXbCBo9Z+73ZsNEGLYFrPhyzqiAgYz9zY7Lcc2JjE/Gy8VB60hCCGEx8lTcHD58mIMHD1KzZs2CziNE4XlkKpxYA9cPwdHl0PBZrRPlWY+aTfhmSzKHb0ez/O8rvNyhqtaRhBDCYuTpT9fatWsTGRlZ0FmEKFyuPtB2vHl+8zRIjtc0Tn71b1ERgJ/2h5FmKlbPBQghxAPlqbj54IMPmDBhAtu2beP27dvExsZmmISwWC1egpL+EB8Buz7ROk2+9KhfFnfPc9xx+4I5+9ZqHUcIISxGnh4Fv//o6X/vtVFVFUVRSEtLK5h0hUAeBReEroPlz4HeAUb9BSUraZ0oz/r8/Can7v6Gm1qP3YOk3xshhO3Kzed3nu652bp1a56CCWERavYA/3ZwcTv8MQX6/qB1ojx7rcVAXtj6GzEcZ1/YaVpUqKF1JCGE0Jx04ieKpxsn4KvWoJpg0Hqo1ErrRHnWavGzxCrHqeXUi5+fnqF1HCGEKBRF0onfzp07ef755wkMDOTatWsA/PDDD+zatSuvqxSi6PjUgcYDzfPBU6x61PA+1Z4GIDR+M9F3EzROI4QQ2stTcbNixQq6dOmCk5MThw4dIjk5GYC4uDjee++9Ag0oRKFpPwnsSsC1g1Y9avjI5j1RjCVBn8gne2TgWiGEyFNx8+677/LVV1/xzTffYGdnl748MDCQQ4cOFVg4IQqVqw+0etU8v+VtMCZrmyeP7A0GmpbqDsD6yys1TiOEENrLU3Fz+vRp2rZtm2m5m5sb0dHR+c0kRNFp+TK4+EDUJTiwUOs0eTah1QBMdysRfaMJx65Gax1HCCE0lafipmzZspw7dy7T8l27dlG5cuV8hxKiyDi4QIc3zfPbP4C70ZrGyasa3r50dJ+OMSaAJftl8FohRPGWp+LmxRdf5NVXX2X//v0oisL169dZsmQJ48ePZ+TIkQWdUYjC1fB58K4Jd6Ng16dap8mz/i3NPRavOXydmLupGqcRQgjt5Km4mTBhAr1796ZDhw7Ex8fTtm1bhg0bxosvvsioUaMKOqMQhUtvgE5vm+f3zYPoK9rmyaMmFUtSrYweo8t2ZmxfonUcIYTQTK6Lm7S0NLZv385rr71GZGQkf/31F/v27ePWrVu88847hZFRiMJXvQtUagNpyfDnu1qnyRNFUahX4yKOZdbxx7UfMJlMWkcSQghN5Lq40ev1dOnShZiYGJydnWnSpAnNmjXDxcWlMPIJUTQUBTpPN88fXQ7hR7TNk0cT2/RDNTmQZrjJ4kObtY4jhBCayNNlqXr16nHhwoWCziKEtso1hrpPASoET7XKjv18XNzxd2wDwA8nl2mcRgghtJGn4mbGjBmMHz+edevWER4eLqOCC9vRcQro7eHCNji/Res0efJyQH8AbpkOcurWVY3TCCFE0cvXqOCQcWRwGRVc2IRNk2Hvl1C6DozYCTq91olyrenC3iTpz9PM/TkW9H5D6zhCCJFvMiq4EPnR5jUI+QFunoAjy6DRc1onyrXuFR9n5dWP+fv2elKM47E35OlQF0IIqySjgguRld2fmwfUdPWFVw6CvbPWiXIl+m4CbZZ2JjXBn4/av0XPetW1jiSEEPkio4ILkV/NhoN7BYi7DvvnaZ0m1zycSvBU6bkkXevH6oMxWscRQogiJaOCC5EVO0fzzcUAOz+FhEht8+RB/+ZVAfjz9E2uRiVqnEYIIYqOjAouRHbqPgVlG0BKHOz4WOs0uVbZ24VWVUuB3S0+3rlG6zhCCFFkZFRwIbKj00Gnaeb5Awsg6rKmcfKiRe0oXKrMYkvkFySmJmsdRwghioSMCi7Eg1R5BPzbQVoKbJupdZpcG9r0EUhzBX0cn+9drXUcIYQoEjIquBAP0+kt879HlsGNE9pmySVnOwfquwcBsOb8Co3TCCFE0ZBRwYV4mHIBUPsxQIUt07VOk2uvtRyAqirE60LZefGk1nGEEKLQ5aufm8TERE6ePInJZKJ27dpWMXim9HMj8iTyHMxpBmoaDN4IFVtqnShX2ix+nmjlCNUde7Ci7/taxxFCiFwrkn5uABkVXBQfXlWhsXnMJja/ZXWDavat8TQAZxK2En03QeM0QghRuPJU3CQkJDBlyhQCAwOpWrUqlStXzjAJYZPavQEGJ7iyH05v0DpNrgxv2h3FWBJVhR8P/aV1HCGEKFR5GnBm2LBhbN++nf79+1O2bNkMg2cKYbPcykKLEbDrU/O9N9W7WM2gmvYGA33Kv8WCbXH8mWzPqFZaJxJCiMKTp3tuPDw8+P3332nVyvr+h5R7bkS+3I2GzxpAUjQ8NteqBtW8GZdE4Mw/MZpU1o9uQ21f+f0XQliPQr/npmTJknh6euYpnBBWzckD2owzz299D1KTNI2TG6VdHelStwygMn/vPq3jCCFEoclTcfPOO+8wdepUEhNlvBpRDDUbDm7lIPYq/P2t1mlypXtDJ0pUnsXm2De5ES8DagohbFOOL0s1atQow701586dQ1VVKlWqlGF8KcCix5eSy1KiQBz6Hta+Ak4l4dUj4OiudaIcMZlMNF7UmTTDTbqWGcVHXV7UOpIQQuRIbj6/c3xDce/evfObSwjb0aAf7PkCIs/A7s//GUHcwul0OgJLP8rOOwvZcnU1JtML6HT56hFCCCEsTr468bNGcuZGFJiTa+Hn/mDnDKNDwLWM1olyJCz6Ft1XBaHojExpPI+n67XWOpIQQjxUkXXid/DgQX788UeWLFlCSEhIflYlhPWp1RPKNYHURNj+odZpcqyChzfl7Mw9LH97ZInGaYQQouDlqbi5efMmjzzyCE2bNmX06NGMGjWKgIAAOnbsyK1btwo6oxCWSVGg0zTz/KHv4PZ5TePkxvBGzwNwPXUfF+7c0DiNEEIUrDwVN6+88gqxsbGcOHGCO3fuEBUVxfHjx4mNjWX06NEFnVEIy+XfBqp2BpMR/nxH6zQ59nitFtinVUDRGfl0969axxFCiAKVp+Jm48aNzJs3j1q1aqUvq127NnPmzGHDBuvqll6IfOs0DVDgxCq4dlDrNDmi0+l4uvIIEi+/QEhoTUymYnXrnRDCxuWpuDGZTJke/waws7PDZDLlO5QQVqVMXWjwjHk+2HoG1Rwd2B0XtSZX7ySx/axcThZC2I48FTePPPIIr776KtevX09fdu3aNcaOHUvHjh0LLJwQVqPDm6C3h0s74dxmrdPkiJO9nj4B5QH4bq/13C8khBAPk6fi5ssvvyQuLo5KlSpRpUoVqlatir+/P3FxcXzxxRcFnVEIy+dRwdxzMZjP3pjStM2TQ88298Oh9DoOpI3h76vntI4jhBAFIk+jgvv5+XHo0CGCg4M5deoUqqpSu3ZtOnXqVND5hLAebV6DQz/AzRNw9Gdo+KzWiR6qircrpUpGEadL5OO9i1ne512tIwkhRL7l6szNn3/+Se3atYmNjQWgc+fOvPLKK4wePZqmTZtSp04ddu7cmasAc+fOxd/fH0dHRwICAnL8/t27d2MwGGjYsGGutidEoXH2hDZjzfNbZ1jNoJpPVXsagJNxm4lJkvHihBDWL1fFzezZs3nhhRey7BnQ3d2dF198kU8++STH61u+fDljxoxh8uTJhISE0KZNG7p160ZYWNgD3xcTE8OAAQPk/h5heZqPAFdfiLkCf3+jdZocGdn8URSjB+gT+HSPPBYuhLB+uSpujhw5QteuXbN9PSgoiIMHc/4o7CeffMLQoUMZNmwYtWrVYvbs2fj5+TFv3rwHvu/FF1+kX79+tGzZMsfbEqJI2DmZby4G2PEx3I3SNk8OONrZ06RUdwB+v7RC4zRCCJF/uSpubty4keUj4PcZDIYc91CckpLCwYMHCQoKyrA8KCiIPXv2ZPu+RYsWcf78ed56660cbSc5OZnY2NgMkxCFqmE/8K4FSdGw61Ot0+TIxFYDUFU9SfoL/Bb6t9ZxhBAiX3JV3JQrV45jx45l+/rRo0cpW7ZsjtYVGRlJWloaPj4+GZb7+PgQERGR5XvOnj3LG2+8wZIlSzAYcnYv9MyZM3F3d0+f/Pz8cvQ+IfJMp/9nWIZ9X0HMVU3j5EQN73L46JsC8FXIDxqnEUKI/MlVcdO9e3emTp1KUlLmGyXv3r3LW2+9xaOPPpqrAIqiZPhaVdVMywDS0tLo168fb7/9NtWrV8/x+idNmkRMTEz6dOXKlVzlEyJPqneBCoGQlgxbZ2qdJkeG1e9Pyu02XDjXjJi7qVrHEUKIPFNUNefdqd64cYPGjRuj1+sZNWoUNWrUQFEUQkNDmTNnDmlpaRw6dCjT2ZispKSk4OzszC+//MLjjz+evvzVV1/l8OHDbN++PUP76OhoSpYsiV6vT19mMplQVRW9Xs8ff/zBI4888tDt5mbIdCHy5crfsKATKDoYsRt8amud6IFUVaXbZzs5FRHH1EdrM6S1v9aRhBAiXW4+v3N15sbHx4c9e/ZQt25dJk2axOOPP07v3r158803qVu3Lrt3785RYQNgb29PQEAAwcHBGZYHBwcTGBiYqb2bmxvHjh3j8OHD6dOIESOoUaMGhw8fpnnz5rnZFSEKn19TqNULVBNseVvrNA+lKArPt6gIwI/7LpOLv3uEEMKi5LoTv4oVK7J+/XqioqI4d+4cqqpSrVo1SpYsmeuNjxs3jv79+9OkSRNatmzJ/PnzCQsLY8SIEYD5ktK1a9f4/vvv0el01K1bN8P7S5cujaOjY6blQliMjlPh1O9wZiNc2g2VWmmd6IF6NyrH+1t/57rDDhYeTGFok6CHv0kIISxMnnooBihZsiRNmzbN18b79u3L7du3mT59OuHh4dStW5f169dTsaL5r8fw8PCH9nkjhEXzqgYBA+HAQgieCsM2Qxb3lFkKFwcDlSud42LKCX448ZMUN0IIq5Sre25sgdxzI4pc3A34vCGkJsJTi6DuE1oneqDN544wdvfzqKqOn7quoX6ZSlpHEkKIwrvnRgiRB64+0OpV8/zmtyx+WIZOVRvgbKqOopj4aPf3WscRQohck+JGiKIQ+Aq4loXoMNj/4B64LUHvyn0AOByzkcTUZI3TCCFE7khxI0RRsC8BHe/1qr1jFsTf1DbPQ7wa+DikuYI+jk92y5AMQgjrIsWNEEWlfl/wbQQpceZRwy2Ys50DjTzM402tvrBM4zRCCJE7UtwIUVR0Oujynnn+0Pdw44S2eR5icpshmJLLEHOrPgcu39Y6jhBC5JgUN0IUpYqBUPsxc8d+m94EC35YsYa3L109PiI1qiXf7ZEuGYQQ1kOKGyGKWqe3QW8PF7bBmU1ap3mg+0MwrD8WTnjMXY3TCCFEzkhxI0RR8/SHFi+Z5//4H6RZ7iCVdXzdaVbZFcX1LyZvma91HCGEyBEpboTQQpvXwNkLbp81915swVrUicTJdwV/Ry8hKjFe6zhCCPFQUtwIoQVHd+jwpnl+20y4G6VtngcY1bwXOmMp0N9lxo4ftY4jhBAPJcWNEFppPBC8a5kLm+0fap0mW/YGA218Hgdg87VfMaalaZxICCEeTIobIbSiN0CXe/3d/DUfIs9pm+cBprQfCCYH0gw3+ObARq3jCCHEA0lxI4SWqnaEakFgMkLwFK3TZMvHxYNqzh0B+CFULk0JISybFDdCaC3oXVD0cHo9XNiudZpsvRE4DFVViFOO8+f5o1rHEUKIbElxI4TWvGtA06Hm+U2TwWSZ97Q086tGKaURxvjqrAy5onUcIYTIlhQ3QliC9pPMT1DdOAaHvtM6Tbbea/Uhd68MYfMRHVEJKVrHEUKILElxI4QlcPaE9vceDd8yHRLvaJsnG4FVSlPH142kVBNL/5YhGYQQlkmKGyEsRdNhULqO+dHwLdO1TpMlRVEY0sofxRDDt8e+ITE1WetIQgiRiRQ3QlgKvQG6f2SeP7gYrodoGic73ev74OI/j1S39Xyye4XWcYQQIhMpboSwJJVaQb0+gAq/jweTSetEmTjZ2dG4ZFcAVl9YpnEaIYTITIobISxN53fA3gWuHYAjP2mdJkv/azsE1aQnWX+Rn4/t0jqOEEJkIMWNEJbGrSy0f8M8H/yWRY47Vc2rLOXtWgHwVchibcMIIcR/SHEjhCVqPgK8akBiJGydqXWaLI1uOgSAm6a/Cbl+UeM0QgjxDyluhLBEejvofm8wzb+/gYhj2ubJQvcaAZQw1UBRTLy782ut4wghRDopboSwVJXbQ+3eoJpg/eugqlonymRgnSGoJjtOXb/LHenUTwhhIaS4EcKSdZkBds4QtheO/qx1mkxebNIVv4QZJNwIYvFuuTQlhLAMUtwIYcncy0Pb8eb54CmQFKttnv/Q6XSMbtcQgMV7LhGXlKptICGEQIobISxfy1HgWQXib8D2D7ROk0mXOmWo7F2CeOU87/5peWeXhBDFjxQ3Qlg6gwN0u3dz8b55cDNU2zz/odMpdGx8mxKV5rE+Yg4xSYlaRxJCFHNS3AhhDap1gpqPgppmkTcXj239KIqxJOjjeHub5Y5qLoQoHqS4EcJadHkPDI5waScc+1XrNBk42znQoezTAGy+tpykVHlySgihHSluhLAWJSv+c3PxpkmQeEfbPP8xrf1gSHNBNdzm/Z0y5pQQQjtS3AhhTQJfNfdcnHALNk/TOk0GJZ1L0MyzNwBrLv6IMS1N20BCiGJLihshrInBHnrONs8f+g4u79E0zn+988gLYHLEaAjns71rtI4jhCimpLgRwtpUDITGA8zzv40Bo+Xc3+Lr5kldl26YUjzYcPw6qoXd+CyEKB6kuBHCGnV6G0p4Q+Rp2POZ1mky+LDzOIyXJ3Lukj97zt/WOo4QohiS4kYIa+TsCV3ujRa+/SO4fV7bPP/i5+HBs838AZiz9ZzGaYQQxZEUN0JYq3pPQZVHIC0Z1o21qL5vXmhbGYPOxF+Rm1h1Yr/WcYQQxYwUN0JYK0WBHrPMfd9c3A5Hl2udKF05Dydq1tmBk++vfHpgjtZxhBDFjBQ3Qlgzz8rQboJ5ftObFtX3zfiWA1FVhShCCD57WOs4QohiRIobIaxdy1fAuxYk3oY/pmidJl37ynXx1gUA8OG+eRqnEUIUJ1LcCGHtDPbQ894TU4d/hIs7tc3zL+OavgRAeNo+9oed1TiNEKK4kOJGCFtQoTkEDDbPrxsLxmRt89zTs1YzXNW6KIqJqTs+1TqOEKKYkOJGCFvR6S0oURpun4VdllNIvNb0FQCuGXfx5/mjGqcRQhQHUtwIYSucSkK3983zO2fBrTPa5rnnyTqBeNKYtER/Fu2+qHUcIUQxIMWNELakzhNQtTOkpcCakWCyjMEr53T+mOQrL7DzpIHDV6K1jiOEsHFS3AhhSxTFPLCmvStc/Rv2zdU6EQB1fb15orEfAB9tOqVxGiGErZPiRghb414euswwz//5LkRaxlNKYzpVw97uLgdiv+PbA39oHUcIYcOkuBHCFjUeAJU7gDEJ1rxsEZenypd0pm6dA9iX2sVXRz/HZDJpHUkIYaOkuBHCFikK9PrCfHnqyn7Y/5XWiQB4v9NoVJM9yfrLzNq9Qus4QggbJcWNELbKww+C3jHPb5luESOHV/fypYFbTwCWnJlPitGocSIhhC2S4kYIWxYwCCq3N1+eWm0ZT0992Hk0pDmRZojg7a3fax1HCGGDpLgRwpalX55ygSv7YP/XWieinJsngV59APgtbDFxyXc1TiSEsDVS3Ahh6zwqWNzlqQ86v4SS5oZqiOLNzd9oHUcIYWOkuBGiOAgYDP5twXj33tNT2j6p5OHkQrdyg0iJasHuI74kJMu9N0KIgiPFjRDFgaJAry/BrgSE7YW/5mudiHc7DaVsaj/uxDqxcJcMyyCEKDiaFzdz587F398fR0dHAgIC2LlzZ7ZtV65cSefOnfH29sbNzY2WLVuyadOmIkwrhBUrWRGCppvnN0/T/PKUnV7HuKAaAMzfcYHIeLn3RghRMDQtbpYvX86YMWOYPHkyISEhtGnThm7duhEWFpZl+x07dtC5c2fWr1/PwYMH6dChAz179iQkJKSIkwthpQKG/HN5au0rml+eerReWaqVi8fo/S0vrX9X0yxCCNuhqKqqarXx5s2b07hxY+bNm5e+rFatWvTu3ZuZM2fmaB116tShb9++TJ06NUftY2NjcXd3JyYmBjc3tzzlFsKqRV2CuYGQmgCd34FWozWN8+W+3/j69JuoJgNLu62mXpmKmuYRQlim3Hx+a3bmJiUlhYMHDxIUFJRheVBQEHv27MnROkwmE3FxcXh6embbJjk5mdjY2AyTEMVayUr/jD21ZTpc1/bM58hmPXBKq4aiMzI2eIamWYQQtkGz4iYyMpK0tDR8fHwyLPfx8SEiIiJH65g1axYJCQk8/fTT2baZOXMm7u7u6ZOfn1++cgthEwIGQa2eYEqFX4dAcpxmUXQ6HZNbTERVFW6Y9rL44GbNsgghbIPmNxQripLha1VVMy3LytKlS5k2bRrLly+ndOnS2babNGkSMTEx6dOVK1fynVkIq6co0PNzcCsHdy7A+gmaxnmsdnP8HToC8NnhD0hMTdY0jxDCumlW3Hh5eaHX6zOdpbl582amszn/tXz5coYOHcrPP/9Mp06dHtjWwcEBNze3DJMQAnD2hCe+AUUHR36CY79qGmdut/9BWgmMhgjGbfxS0yxCCOumWXFjb29PQEAAwcHBGZYHBwcTGBiY7fuWLl3KoEGD+Omnn+jRo0dhxxTCtlVqBW1fN8+vGwt3tOtvxs+jFD39XgBgV3gw16ITNMsihLBuml6WGjduHN9++y0LFy4kNDSUsWPHEhYWxogRIwDzJaUBAwakt1+6dCkDBgxg1qxZtGjRgoiICCIiIoiJidFqF4Swfm0ngF8LSI6FFcMgLVWzKNMfGYRP8gDiL77EzPWnNcshhLBumhY3ffv2Zfbs2UyfPp2GDRuyY8cO1q9fT8WK5kdBw8PDM/R58/XXX2M0Gnn55ZcpW7Zs+vTqq69qtQtCWD+9AZ78Bhzc4doB2JazbhgKg0Gv57NHh6PDwLqj4ew+F6lZFiGE9dK0nxstSD83QmTjxCr4ZRCgwIA1ULmdZlGmrT3B4j3nKVv+GJtfGE8JBwfNsgghLINV9HMjhLAwdR6HxgMAFVa9CAm3NYsytnN13Py/I951Ga9unK1ZDiGEdZLiRgjxj67vg1d1iAuHtaNAoxO77k52PF61FwD77iwn5LoMrCmEyDkpboQQ/7AvAU8uAL09nF4Pf3+rWZS3OvS/13NxKmOD39YshxDC+khxI4TIqGx96Hxv9PBNkyH8qCYxdDod77SZgqrquM1B5u3/XZMcQgjrI8WNECKz5iOgWhdIS4Zlz0GCNk8tdanWiBpO3QD46vgsYpMTNckhhLAuUtwIITJTFHjia/CsAjFh8PMAzfq/mdN9IqS5YTLc4tX1n2qSQQhhXaS4EUJkzakkPLsU7F3h8m7Y+IYmMcq4luTpyiMxJviz97gvV+7I2RshxINJcSOEyJ53DXjyW0Ax31x8YJEmMSa3fZaGdpNISijNaz8fIc1UrLrnEkLkkhQ3QogHq9EVOk4xz68fD5f3FHkEnU7Hh082wMXBwF+X7vBBcNFnEEJYDyluhBAP13oc1HkCTEZY3h+irxR5BD9PZ95+rCYOpdfzU/hIlh/dWeQZhBDWQYobIcTDKQo8NgfK1IfESFjWD1KK/t6XJxtXxM87FUUx8d7fU4iIiyryDEIIyyfFjRAiZ+yd4ZmfwNkLIo5q1oPx94+9j2L0xGS4zYDVk4p8+0IIyyfFjRAi5zz8oO8PoDPA8RWwq+gfzfZ182RSk+moqkK4aTfTt/5Y5BmEEJZNihshRO5UDITuH5nnt0yHM5uKPMKzDdrR2K0PAL9c+oy/r54r8gxCCMslxY0QIveaDDFPqLBiGNw4WeQR5veciGNaFdAlMfKP10hO1aaTQSGE5ZHiRgiRN10/gIqtIDkWfngc7lwo0s072tkzN+hjMLoTFdGUedtl5HAhhJkUN0KIvDHYwzNLoHQdiI+A73tD7PUijdC0fFWmNfoBY2wjPt9yloOX7xTp9oUQlkmKGyFE3jmVhP6rwLMyRF82n8FJLNoC48nG/jzeqBwmFV5ZvovrsVLgCFHcSXEjhMgfVx8YsAZcfeHWKfjxCUiKLdII0x+rQ1mfa8R4vs+ANdqMgSWEsBxS3Agh8s+jAgxYDc6l4HoILH0WUu8W2eZdHe0Y07EWiiGBG6a9vL7p6yLbthDC8khxI4QoGN414PmV4OAGl3fBL4MgreieYOpbvw3NPZ4FYEP4XL7c91uRbVsIYVmkuBFCFBzfhtBvORgc4cxGWDUCTGlFtvlvek3AV98GRTHxVeg01pzcX2TbFkJYDiluhBAFq2Ig9P3xXi/Gv5pHEi+iYRp0Oh2r+n6Kq6k2ii6F/+0bKx38CVEMSXEjhCh41TrDE98AChxYCFveLrJNO9s5sPKpr7BLKw/6OF76fQbRiSlFtn0hhPakuBFCFI66T0DPz8zzuz6F4KlFdganjGtJfugxH7uElkRe6sWw7w6QlFp0l8eEENqS4kYIUXgCBkKXmeb53Z/BmpchzVgkm67j48eyJ2fh6uDMgctRjPv5MGlppiLZthBCW1LcCCEKV8uR8NgcUPRweAksfw5SEotk09V9XJnfvwn2eoXN4T/y9K9vFsl2hRDakuJGCFH4Gj0Pz/wEBifzU1TfP1ZkPRm3rFKK0d3tcfDezJmk3xm17rMi2a4QQjtS3AghikaNruaejB094OpfsLArxFwtkk2PatWJFh4DANgWuYAPd/xcJNsVQmhDihshRNGp0ByGbDQP1RB5GhYEwa3TRbLpr3u+hr99ZxRF5fvzM/nm701Fsl0hRNGT4kYIUbRK14Khf4BXdYi9Bgu7wJW/Cn2zOp2OX/t8SEkaoeiMfHZiIu9u+6nQtyuEKHpS3Aghip6HHwzZBOWawN0o+K4XnCn8Myn2BgO/9f0KL6UJipLGsksf8OnWvYW+XSFE0ZLiRgihDWdPGLgWqnYG413zYJt75xR6Xzjujs5s6vcNVey7kBzxGJ9tusPMDaGoRdQHjxCi8ElxI4TQjn0JeHYpNHgW1DTY9CYsf958NqcwN2swsOqZjxjT3HyT8dfbLzDq5+0kpiYX6naFEEVDihshhLb0dtB7HnT/GPT2cGodfN0Wrh4s1M0qisLLHary0VP10dslsi3mbTotGcSt+NhC3a4QovBJcSOE0J6iQLMXzDcae1SE6DDzjcb75hX6Zao+TfyY2NMDxS6GOOU43X7ux/nbEYW6TSFE4ZLiRghhOXwbwYs7oFYvMKXCxjfuXaaKLtTNDm/WmWlNvoA0Z5L1l3lidT/2h50t1G0KIQqPFDdCCMvi5AFPfw/dPgKd3T+Xqa4V7mWqp+q1Yk6HhShGT0yGWwzbPJClR7YX6jaFEIVDihshhOVRFGg+/F+XqS7Dgi6w76tCvUzV1r8OSx/9EYOxHOjjmBHyCi+vWUSKUQbcFMKaSHEjhLBc5RqbL1PVfPTeZaqJsLgH3DhZaJus4+PHhqeXUlbfGtXoxu/7XXhi3m7O3YwrtG0KIQqWohazzh1iY2Nxd3cnJiYGNzc3reMIIXJCVWH/17B5mrlPHEUPLV6CdhPBsfCO419CTjHjtzCiE1NxMCj0bRfPtI5PodPJ34VCFLXcfH5LcSOEsB7RYea+cEJ/M3/tUga6zIC6T5ovZRWCG7FJjP/lCPtubsDRdwVual2+6fYhtX38CmV7Qois5ebzW/78EEJYD48K0PdHeG4FeFaG+AhYMRS+6wk3TxXKJn3cHPlucDN6NvBBNRmIVY7T9/en+GjnL4WyPSFE/smZGyGEdUpNgr1fwI5Z5ktVOsM/l6ocXAtlk3+eP8r47W+Qqr8CgK++Dd88+jYVPLwLZXtCiH/IZakHkOJGCBsTddl8qerUOvPXrmWh9Tho9DzYOxf45hKSkxm2dgbHElajKCqkOdKp9MvMCHoeZ3tDgW9PCGEmxc0DSHEjhI068wdsmABRF81fO5eCZi+aez529izwzS05so2PD7yP0XCNhIsv42moyqsdq9K3aQXsDXLFX4iCJsXNA0hxI4QNS02CkB9gzxfmvnEA7EpAwEBo+TK4ly/QzRnT0vhy7x+s3ONI2J1EAEqX38ujdWryZttnMej1Bbo9IYozKW4eQIobIYqBNCOcXA27Z0PEMfMynQHq9YFWr0LpWgW6uRSjieV/h/HZ9r9J8pmBokvDYCzHkNov8XLznvLouBAFQIqbB5DiRohiRFXh/J/mIufijn+WV+8KjfpDtc5gcCiwzd1OjGPCH3P5686voE8CwCmtKn2qPcfLLXribFdw2xKiuJHi5gGkuBGimLp2EHbNvtdHzr3/9hzdzYN01usDlVqDrmAuI4VF32J88GxOxq9H0RnNC9Nc6FRqAq8Edqayt0uBbEeI4kSKmweQ4kaIYi7yHBxcBMdXQFz4P8tdypg7A6z3lHl08gLoFPB4RBjTd35NaNxWVOUu8WffBJMzzSp50qm+nj4Na1PSWQodIXJCipsHkOJGCAGAKQ0u74Fjv8DJNZAU/c9rnlXMRU7VzuZCR5+/R7wTU5P5KWQfe086s/X0TUwqOFeci97xJpUcWzOswTP0rNlU7s0R4gGkuHkAKW6EEJkYk+HcFnOhc3qDuVPA++xdoVIr8G8H/m2hdG3IRxESEZPEkr9P8f3lMZgMt9OXK0YPyjk2pG351jzX4BHpGFCI/5Di5gGkuBFCPFBynLnACV0LF3dmPKMD4OxlLnLuT56V83QJy5iWxvchW1gS+gs3jAf+uTcHMMY0pIZ+BG2redGmmhd1y7vgJDcji2JOipsHkOJGCJFjpjTzo+QXt5uftrq8B1ITM7ZxdIcy9aFMvX8mrxpgsM/xZqIS41l2fAebL+7gQvwh4m60xRjbCACdQzglKn1NSV1tqnvUpnGZOjxSuSE1vQu2zx4hLJ0UNw8gxY0QIs+MKeanru4XO1f/hrSUzO10dlC6JpRpAGXqms/ulKxkHvjTzumhm7kenciuc7fZceYWO26sQPVcm7lRmiuuuoq0KNmX9pWaUbusG/5eJdDrCmd0dCG0ZlXFzdy5c/noo48IDw+nTp06zJ49mzZt2mTbfvv27YwbN44TJ07g6+vLhAkTGDFiRI63J8WNEKLAGFMg8rT57E76dBSSYrJ/j0sZc6FTshKUrHiv6KkILj7g4g0Obhkuc6UYjaw7/Tcbz+/ifMwZbqdcxKi/aR7XCki8PJy0xMoAOJU8iqPXn5TQe1PSoQxlS/hSyb0c1UtVoL5PZaqW8pGbloXVspriZvny5fTv35+5c+fSqlUrvv76a7799ltOnjxJhQoVMrW/ePEidevW5YUXXuDFF19k9+7djBw5kqVLl/Lkk0/maJtS3AghCpWqQsyVf4qdG8ch6pJ5gM/k2Ie/X+8ALqWhhBeUKG0ueEp4m+cd3cHRjWjVwM47N9l/5yp2d5twPNKB4xFJmEquwb7UrmxXbbw6grKOdfBycQCnkyQYjuJu705JR0+8nT0pXcITX1cvyrt5UdGjLO6OTjJOlrAYVlPcNG/enMaNGzNv3rz0ZbVq1aJ3797MnDkzU/uJEyeydu1aQkND05eNGDGCI0eOsHfv3hxtU4obIYQmVBXuRt0rdC6Zx75Knw+D+FuQEpf31esdCHd046ydA1cMBq7odFzTKYTrVG7qTcToTbxyyQ8HozMpGNhV6hYHPW9mu77ml1rhnOyJTqfnmuclLrufwaDaYYc99oo9doo9BsWAnWJHM10rSum9MBj03CCcK+pF7HT22OvsMOj1GBQDBp0evaKntnNtStp7oNfpiE6LIsIYgU7RYVD06BQdOkWHXqegU/T4OZXH1d4FnaIQZ4wnMiUSnaKgKAoKirm9TkGvKHjZl6KEXQkAEo13uZMSjaKYT4IpgKKYizQdCu527pQwOIMCSWlJRKXEoPDP2bJ/X9lztXM1twWS05KJSs3+rJyroQQlDOYMKaZUolOjs21bQu+My728qaZU7qRk39ZZ74Srnbk/pDRTGpEpd7Jt66R3xM3O1dxWTSMyOfu2jjoH3O3Nn4Mm1cSt5NvZtnXQ2eNh7w6AqqrcTI7Mtq29zp5STqXwKV8l2zZ5kZvP7/x13pAPKSkpHDx4kDfeeCPD8qCgIPbs2ZPle/bu3UtQUFCGZV26dGHBggWkpqZiZ2eX6T3JyckkJyenfx0bm4O/nIQQoqApinl0cmdPKNc46zYpiZBw658p/iYk3DQXPgm3zGd+kmIz/psSb159WjK+CbfwzWbzqYCeMHT3/tdvkOzAgSgHonR6ovU6ovQ6onX6e//qmGz4BX/V/ATXFzp35hvcs921t6/+Sd0U871HC91d+dGzZLZtF564QdMk8//JS11d+Mwr+xHb50TcpM5d8zAWq1xKMNW7VLZtP75xiy6J5kf4N5Zw5vXSXtm2fffWbR6LTwBgu5Mjo8qUzrbt5Mg7PBNn/h7/5ejA0LI+2bZ97XYUg2LNBeoxe3v6lyuTbduRUdG8FG3+PDprZ0f/8mWzbTs4OpZxUdEAXDXo6e9XLtu2z8TGMfl2FAC3dTr6V8z+xvPH4uJ5N9Jc/CQqCl0r+WXbNig+gVm3zMWPCnTxz3x15b62iXeZdsMI0y5l26awaVbcREZGkpaWho9Pxl8UHx8fIiIisnxPRERElu2NRiORkZGULZv5l2PmzJm8/fbbBRdcCCEKi70z2Fc034uTU6Y08+Pr9wseY5L5ia7UJHN/Panmyc6YlD6PKZWmaak0TUuBtFTzZEo13xydZv7XVKEGqUYjprRUeqel0OJuCvEmI4lqGomkcVdNIxmVVFSc7MpwR6eCqlLOpNIzPo0UBVIVlTQwTwoYAb1aghidI6gqDiY7aiQbMQGqgvlfIA0FVQFUe+LRoaCiqHZ4G9PgXhsVMN27N0kFFNXAXcyPy5tUA25ppvR23Ft/+rfMpOeuan6azagacDaZ0l/776UMk/pP21TVgKPJhErWN22bVB1Jql16WwdT9hdGVFWf3jbloW11/2qrf2BbRVXS2yaruoe0/We9SSg5bqvCA9vqVIVUJedPCxYGzYqb+5T/9A+hqmqmZQ9rn9Xy+yZNmsS4cePSv46NjcXPL/vqVAghrIpOD04e5qkgV3tvAvC7N+VEl3tTTjxxb8qJ3vemnOh+b8qJjsD+HLZtBfydw7aNgQM5bFsnF20r56Ktby7aOuaiLblsqwXNihsvLy/0en2mszQ3b97MdHbmvjJlymTZ3mAwUKpU1qcrHRwccHCQzq+EEEKI4kKz2+Dt7e0JCAggODg4w/Lg4GACAwOzfE/Lli0ztf/jjz9o0qRJlvfbCCGEEKL40fQZv3HjxvHtt9+ycOFCQkNDGTt2LGFhYen91kyaNIkBAwaktx8xYgSXL19m3LhxhIaGsnDhQhYsWMD48eO12gUhhBBCWBhN77np27cvt2/fZvr06YSHh1O3bl3Wr19PxYrmm+nCw8MJCwtLb+/v78/69esZO3Ysc+bMwdfXl88//zzHfdwIIYQQwvZp3kNxUZN+boQQQgjrk5vPb+l6UgghhBA2RYobIYQQQtgUKW6EEEIIYVOkuBFCCCGETZHiRgghhBA2RYobIYQQQtgUKW6EEEIIYVOkuBFCCCGETZHiRgghhBA2RdPhF7Rwv0Pm2NhYjZMIIYQQIqfuf27nZGCFYlfcxMXFAeDn56dxEiGEEELkVlxcHO7u7g9sU+zGljKZTFy/fh1XV1cURSnQdcfGxuLn58eVK1dsctwqW98/sP19lP2zfra+j7J/1q+w9lFVVeLi4vD19UWne/BdNcXuzI1Op6N8+fKFug03Nzeb/aUF298/sP19lP2zfra+j7J/1q8w9vFhZ2zukxuKhRBCCGFTpLgRQgghhE2R4qYAOTg48NZbb+Hg4KB1lEJh6/sHtr+Psn/Wz9b3UfbP+lnCPha7G4qFEEIIYdvkzI0QQgghbIoUN0IIIYSwKVLcCCGEEMKmSHEjhBBCCJsixc0DzJ07F39/fxwdHQkICGDnzp0PbL99+3YCAgJwdHSkcuXKfPXVV5narFixgtq1a+Pg4EDt2rVZtWpVYcXPkdzs48qVK+ncuTPe3t64ubnRsmVLNm3alKHN4sWLURQl05SUlFTYu5Kl3Ozftm3bssx+6tSpDO0s6WeYm/0bNGhQlvtXp06d9DaW9PPbsWMHPXv2xNfXF0VRWL169UPfY23HYG730dqOwdzunzUeg7ndR2s6DmfOnEnTpk1xdXWldOnS9O7dm9OnTz/0fZZwHEpxk43ly5czZswYJk+eTEhICG3atKFbt26EhYVl2f7ixYt0796dNm3aEBISwptvvsno0aNZsWJFepu9e/fSt29f+vfvz5EjR+jfvz9PP/00+/fvL6rdyiC3+7hjxw46d+7M+vXrOXjwIB06dKBnz56EhIRkaOfm5kZ4eHiGydHRsSh2KYPc7t99p0+fzpC9WrVq6a9Z0s8wt/v32WefZdivK1eu4OnpSZ8+fTK0s5SfX0JCAg0aNODLL7/MUXtrPAZzu4/Wdgzmdv/us5ZjEHK/j9Z0HG7fvp2XX36Zffv2ERwcjNFoJCgoiISEhGzfYzHHoSqy1KxZM3XEiBEZltWsWVN94403smw/YcIEtWbNmhmWvfjii2qLFi3Sv3766afVrl27ZmjTpUsX9Zlnnimg1LmT233MSu3atdW33347/etFixap7u7uBRUxX3K7f1u3blUBNSoqKtt1WtLPML8/v1WrVqmKoqiXLl1KX2ZJP79/A9RVq1Y9sI01HoP/lpN9zIolH4P/lpP9s7Zj8L/y8jO0puPw5s2bKqBu37492zaWchzKmZsspKSkcPDgQYKCgjIsDwoKYs+ePVm+Z+/evZnad+nShQMHDpCamvrANtmtszDlZR//y2QyERcXh6enZ4bl8fHxVKxYkfLly/Poo49m+quyKORn/xo1akTZsmXp2LEjW7duzfCapfwMC+Lnt2DBAjp16kTFihUzLLeEn19eWNsxWBAs+RjMD2s4BguKNR2HMTExAJl+3/7NUo5DKW6yEBkZSVpaGj4+PhmW+/j4EBERkeV7IiIismxvNBqJjIx8YJvs1lmY8rKP/zVr1iwSEhJ4+umn05fVrFmTxYsXs3btWpYuXYqjoyOtWrXi7NmzBZr/YfKyf2XLlmX+/PmsWLGClStXUqNGDTp27MiOHTvS21jKzzC/P7/w8HA2bNjAsGHDMiy3lJ9fXljbMVgQLPkYzAtrOgYLgjUdh6qqMm7cOFq3bk3dunWzbWcpx2GxGxU8NxRFyfC1qqqZlj2s/X+X53adhS2veZYuXcq0adNYs2YNpUuXTl/eokULWrRokf51q1ataNy4MV988QWff/55wQXPodzsX40aNahRo0b61y1btuTKlSt8/PHHtG3bNk/rLGx5zbJ48WI8PDzo3bt3huWW9vPLLWs8BvPKWo7B3LDGYzA/rOk4HDVqFEePHmXXrl0PbWsJx6GcucmCl5cXer0+UxV58+bNTNXmfWXKlMmyvcFgoFSpUg9sk906C1Ne9vG+5cuXM3ToUH7++Wc6der0wLY6nY6mTZsW+V8c+dm/f2vRokWG7JbyM8zP/qmqysKFC+nfvz/29vYPbKvVzy8vrO0YzA9rOAYLiqUeg/llTcfhK6+8wtq1a9m6dSvly5d/YFtLOQ6luMmCvb09AQEBBAcHZ1geHBxMYGBglu9p2bJlpvZ//PEHTZo0wc7O7oFtsltnYcrLPoL5r8VBgwbx008/0aNHj4duR1VVDh8+TNmyZfOdOTfyun//FRISkiG7pfwM87N/27dv59y5cwwdOvSh29Hq55cX1nYM5pW1HIMFxVKPwfyyhuNQVVVGjRrFypUr+fPPP/H393/oeyzmOCywW5NtzLJly1Q7Ozt1wYIF6smTJ9UxY8aoJUqUSL+j/Y033lD79++f3v7ChQuqs7OzOnbsWPXkyZPqggULVDs7O/XXX39Nb7N7925Vr9er77//vhoaGqq+//77qsFgUPft21fk+6equd/Hn376STUYDOqcOXPU8PDw9Ck6Ojq9zbRp09SNGzeq58+fV0NCQtTBgwerBoNB3b9/v8Xv36effqquWrVKPXPmjHr8+HH1jTfeUAF1xYoV6W0s6WeY2/277/nnn1ebN2+e5Tot6ecXFxenhoSEqCEhISqgfvLJJ2pISIh6+fJlVVVt4xjM7T5a2zGY2/2ztmNQVXO/j/dZw3H40ksvqe7u7uq2bdsy/L4lJiamt7HU41CKmweYM2eOWrFiRdXe3l5t3LhxhsffBg4cqLZr1y5D+23btqmNGjVS7e3t1UqVKqnz5s3LtM5ffvlFrVGjhmpnZ6fWrFkzw0GrhdzsY7t27VQg0zRw4MD0NmPGjFErVKig2tvbq97e3mpQUJC6Z8+eItyjjHKzfx988IFapUoV1dHRUS1ZsqTaunVr9ffff8+0Tkv6Geb2dzQ6Olp1cnJS58+fn+X6LOnnd/+x4Ox+32zhGMztPlrbMZjb/bPGYzAvv6fWchxmtV+AumjRovQ2lnocKvd2QAghhBDCJsg9N0IIIYSwKVLcCCGEEMKmSHEjhBBCCJsixY0QQgghbIoUN0IIIYSwKVLcCCGEEMKmSHEjhBBCCJsixY0QQgghbIoUN0KIIjFt2jQaNmyo2fanTJnC8OHDC239N2/exNvbm2vXrhXaNoQQOSM9FAsh8k1RlAe+PnDgQL788kuSk5PTRwYuSjdu3KBatWocPXqUSpUqFdp2xo0bR2xsLN9++22hbUMI8XBS3Agh8i0iIiJ9fvny5UydOpXTp0+nL3NycsLd3V2LaAC89957bN++nU2bNhXqdo4dO0azZs24fv06JUuWLNRtCSGyJ5elhBD5VqZMmfTJ3d0dRVEyLfvvZalBgwbRu3dv3nvvPXx8fPDw8ODtt9/GaDTy+uuv4+npSfny5Vm4cGGGbV27do2+fftSsmRJSpUqxWOPPcalS5cemG/ZsmX06tUrw7L27dvzyiuvMGbMGEqWLImPjw/z588nISGBwYMH4+rqSpUqVdiwYUP6e6Kionjuuefw9vbGycmJatWqsWjRovTX69WrR5kyZVi1alXev5lCiHyT4kYIoZk///yT69evs2PHDj755BOmTZvGo48+SsmSJdm/fz8jRoxgxIgRXLlyBYDExEQ6dOiAi4sLO3bsYNeuXbi4uNC1a1dSUlKy3EZUVBTHjx+nSZMmmV777rvv8PLy4q+//uKVV17hpZdeok+fPgQGBnLo0CG6dOlC//79SUxMBMz37Zw8eZINGzYQGhrKvHnz8PLyyrDOZs2asXPnzgL+TgkhckOKGyGEZjw9Pfn888+pUaMGQ4YMoUaNGiQmJvLmm29SrVo1Jk2ahL29Pbt37wbMZ2B0Oh3ffvst9erVo1atWixatIiwsDC2bduW5TYuX76Mqqr4+vpmeq1Bgwb873//S9+Wk5MTXl5evPDCC1SrVo2pU6dy+/Ztjh49CkBYWBiNGjWiSZMmVKpUiU6dOtGzZ88M6yxXrtxDzyQJIQqXQesAQojiq06dOuh0//yN5ePjQ926ddO/1uv1lCpVips3bwJw8OBBzp07h6ura4b1JCUlcf78+Sy3cffuXQAcHR0zvVa/fv1M26pXr16GPED69l966SWefPJJDh06RFBQEL179yYwMDDDOp2cnNLP9AghtCHFjRBCM3Z2dhm+VhQly2UmkwkAk8lEQEAAS5YsybQub2/vLLdx/7JRVFRUpjYP2/79p8Dub79bt25cvnyZ33//nc2bN9OxY0defvllPv744/T33LlzJ9ssQoiiIZelhBBWo3Hjxpw9e5bSpUtTtWrVDFN2T2NVqVIFNzc3Tp48WSAZvL29GTRoED/++COzZ89m/vz5GV4/fvw4jRo1KpBtCSHyRoobIYTVeO655/Dy8uKxxx5j586dXLx4ke3bt/Pqq69y9erVLN+j0+no1KkTu3btyvf2p06dypo1azh37hwnTpxg3bp11KpVK/31xMREDh48SFBQUL63JYTIOyluhBBWw9nZmR07dlChQgWeeOIJatWqxZAhQ7h79y5ubm7Zvm/48OEsW7Ys/fJSXtnb2zNp0iTq169P27Zt0ev1LFu2LP31NWvWUKFCBdq0aZOv7Qgh8kc68RNC2DxVVWnRogVjxozh2WefLbTtNGvWjDFjxtCvX79C24YQ4uHkzI0QwuYpisL8+fMxGo2Fto2bN2/y1FNPFWrxJITIGTlzI4QQQgibImduhBBCCGFTpLgRQgghhE2R4kYIIYQQNkWKGyGEEELYFCluhBBCCGFTpLgRQgghhE2R4kYIIYQQNkWKGyGEEELYFCluhBBCCGFT/g849LVBLIXVdwAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.plot(ts, np.abs(lcce), label='CCE')\n", "plt.plot(ts, np.abs(lmegcce), label='ME-gCCE')\n", "plt.plot(ts, np.abs(lmegcce) / np.exp(-ts/t2), ls='--', label='ME-gCCE / exp(-t/T2)')\n", "\n", "plt.xlabel('Time (ms)')\n", "plt.ylabel('Coherence')\n", "\n", "plt.legend();" ] }, { "cell_type": "markdown", "id": "857153e6-2be6-4fae-b09e-4b3a0ca73a2a", "metadata": {}, "source": [ "Note that calculations with the ME-based methods are significantly more expensive, as we need to solve the eigenvalue problem for the $N^4$ matrix \n", "where $N$ is the size of the Hilbert space, compared to $N^2$ in the case of the closed system." ] } ], "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.13" } }, "nbformat": 4, "nbformat_minor": 5 }