{ "cells": [ { "cell_type": "markdown", "id": "a2a0c72d-95dd-447e-9cc7-dc7c46a98df9", "metadata": {}, "source": [ "# Dynamics of the Spin-Boson Model" ] }, { "cell_type": "markdown", "id": "348ef9f7-cc2c-41ed-b96b-5f20132ffa4c", "metadata": {}, "source": [ "## Overview\n", "\n", "In this notebook we will simulate the dynamics of a 1-mode spin-boson model.\n", "\n", "$$\n", " \\hat H = \\epsilon \\hat \\sigma_z + \\Delta \\hat \\sigma_x + \\omega \\hat b^\\dagger \\hat b\n", " + g \\sigma_z (\\hat b^\\dagger + \\hat b)\n", "$$\n", "\n", "We will first start with implementing it from scratch, and then show how to use the `TimeEvolution` class for the very same task." ] }, { "cell_type": "markdown", "id": "97c09ef0-9d7f-47bf-81c4-2693b7b83e99", "metadata": {}, "source": [ "## Determine the model and basis\n", "\n", "TenCirChem provides a simple yet versatile interface to define the Hamiltonian and the basis set.\n", "The Hamiltonian operators are specified through a list of `Op` object and simple harmonic oscillator (SHO) basis set is used for the boson in the spin-boson model.\n", "The relevant classes are borrowed from [Renormalizer](https://github.com/shuaigroup/Renormalizer)." ] }, { "cell_type": "code", "execution_count": 1, "id": "c334cf97-2078-4fc8-84a5-b70931543f2a", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "([Op('sigma_z', ['spin'], 0.0),\n", " Op('sigma_x', ['spin'], 1.0),\n", " Op('b^\\\\dagger b', ['boson', 'boson'], 1.0),\n", " Op('sigma_z b^\\\\dagger+b', ['spin', 'boson'], 0.5)],\n", " [BasisHalfSpin(dof: spin, nbas: 2),\n", " BasisSHO(dof: boson, x0: 0.0, omega: 1, nbas: 8)])" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from tencirchem import Op, BasisHalfSpin, BasisSHO\n", "\n", "epsilon = 0\n", "delta = 1\n", "omega = 1\n", "g = 0.5\n", "\n", "ham_terms = [\n", " Op(\"sigma_z\", \"spin\", epsilon),\n", " Op(\"sigma_x\", \"spin\", delta),\n", " Op(r\"b^\\dagger b\", \"boson\", omega),\n", " Op(\"sigma_z\", \"spin\", g) * Op(r\"b^\\dagger+b\", \"boson\"),\n", "]\n", "basis = [BasisHalfSpin(\"spin\"), BasisSHO(\"boson\", omega=omega, nbas=8)]\n", "ham_terms, basis" ] }, { "cell_type": "markdown", "id": "0d6a367d-a97b-4f50-b022-5324213461ee", "metadata": {}, "source": [ "`tencirchem.dynamic.model` offers several shotcuts for common models" ] }, { "cell_type": "markdown", "id": "162ea7fb-39f2-44cd-98d8-fe1841ad3183", "metadata": {}, "source": [ "Quantum circuit is not able to simulate the boson basis, and proper basis transform to qubit basis is required" ] }, { "cell_type": "code", "execution_count": 2, "id": "f3ad77b6-07e9-432d-8a31-e2345b4bf059", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "([Op('X', ['spin'], 1.0),\n", " Op('Z', [('boson', 'TCCQUBIT-0')], -2.0),\n", " Op('Z Z', [('boson', 'TCCQUBIT-0'), ('boson', 'TCCQUBIT-1')], -1.0),\n", " Op('Z Z Z', [('boson', 'TCCQUBIT-0'), ('boson', 'TCCQUBIT-1'), ('boson', 'TCCQUBIT-2')], -0.5),\n", " Op('Z X', ['spin', ('boson', 'TCCQUBIT-2')], 0.9517337620166572),\n", " Op('Z Z X', ['spin', ('boson', 'TCCQUBIT-1'), ('boson', 'TCCQUBIT-2')], -0.040295934250509535),\n", " Op('Z Z X', ['spin', ('boson', 'TCCQUBIT-0'), ('boson', 'TCCQUBIT-2')], -0.26872106012443797),\n", " Op('Z Z Z X', ['spin', ('boson', 'TCCQUBIT-0'), ('boson', 'TCCQUBIT-1'), ('boson', 'TCCQUBIT-2')], -0.14271676764170976),\n", " Op('Z X', ['spin', ('boson', 'TCCQUBIT-1')], 0.4829629131445341),\n", " Op('Z X Z', ['spin', ('boson', 'TCCQUBIT-1'), ('boson', 'TCCQUBIT-2')], -0.4829629131445341),\n", " Op('Z Z X', ['spin', ('boson', 'TCCQUBIT-0'), ('boson', 'TCCQUBIT-1')], -0.12940952255126037),\n", " Op('Z Z X Z', ['spin', ('boson', 'TCCQUBIT-0'), ('boson', 'TCCQUBIT-1'), ('boson', 'TCCQUBIT-2')], 0.12940952255126037),\n", " Op('Z X', ['spin', ('boson', 'TCCQUBIT-0')], 0.25),\n", " Op('Z X Z', ['spin', ('boson', 'TCCQUBIT-0'), ('boson', 'TCCQUBIT-2')], 0.25),\n", " Op('Z X Z', ['spin', ('boson', 'TCCQUBIT-0'), ('boson', 'TCCQUBIT-1')], -0.25),\n", " Op('Z X Z Z', ['spin', ('boson', 'TCCQUBIT-0'), ('boson', 'TCCQUBIT-1'), ('boson', 'TCCQUBIT-2')], -0.25)],\n", " 3.5,\n", " [BasisHalfSpin(dof: spin, nbas: 2),\n", " BasisHalfSpin(dof: ('boson', 'TCCQUBIT-0'), nbas: 2),\n", " BasisHalfSpin(dof: ('boson', 'TCCQUBIT-1'), nbas: 2),\n", " BasisHalfSpin(dof: ('boson', 'TCCQUBIT-2'), nbas: 2)])" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# transform the Hamiltonian and basis\n", "from tencirchem.dynamic import qubit_encode_op, qubit_encode_basis\n", "\n", "boson_encoding = \"gray\"\n", "ham_terms_spin, constant = qubit_encode_op(ham_terms, basis, boson_encoding)\n", "basis_spin = qubit_encode_basis(basis, boson_encoding)\n", "# boson basis is transformed to spin basis\n", "ham_terms_spin, constant, basis_spin" ] }, { "cell_type": "markdown", "id": "a0d2d5ea-2b9e-47ec-b21c-0cda341bcc4d", "metadata": {}, "source": [ "## Construct ansatz\n", "\n", "TenCirChem uses the variational Hamiltonian ansatz. Suppose the Hamiltonian has $N$ terms and is written as \n", "\n", "$$\n", "\\hat H = \\sum_j^N \\hat h_j\n", "$$ \n", "\n", "The the variational Hamiltonian ansatz has the form\n", "$$\n", "| \\psi \\rangle = \\prod_k^M \\prod_j^N e^{-i \\theta_{kj} \\hat h_j} | \\phi \\rangle\n", "$$\n", "\n", "where $| \\phi \\rangle$ is the initial state, $\\theta_{kj}$ is the circuit parameter and $M$ is number of layers for the ansatz." ] }, { "cell_type": "code", "execution_count": 3, "id": "2ffaaa62-c1f3-4741-afbc-fb4d5559a0ab", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "DeviceArray([1.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j,\n", " 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j,\n", " 0.+0.j, 0.+0.j], dtype=complex128)" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import tensorcircuit as tc\n", "from tencirchem import set_backend\n", "from tencirchem.dynamic import get_ansatz, get_jacobian_func, get_deriv\n", "\n", "# dynamics simulation requires auto-differentiation from JAX.\n", "set_backend(\"jax\")\n", "\n", "# the initial state\n", "init_circuit = tc.Circuit(len(basis_spin))\n", "# number of layers\n", "n_layers = 3\n", "# get the ansatz. Note that the spin basis is feed in\n", "ansatz = get_ansatz(ham_terms_spin, basis_spin, n_layers, init_circuit)\n", "\n", "# ansatz accepts parameters and outputs wavefunction\n", "import numpy as np\n", "\n", "ansatz(np.zeros(n_layers * len(ham_terms_spin)))" ] }, { "cell_type": "markdown", "id": "30146ff6-d6b4-4f22-aa8a-f0b2ce681522", "metadata": {}, "source": [ "The jacobian matrix for the wavefunction $\\partial | \\psi \\rangle / \\partial \\theta_{kj}$ is required to perform VQA dynamics based on TDVP. " ] }, { "cell_type": "code", "execution_count": 4, "id": "5913f357-9f6c-4321-9424-eb5da6c475e8", "metadata": {}, "outputs": [], "source": [ "# the funcction to evaluate the jacobian of the wavefunction\n", "jacobian_func = get_jacobian_func(ansatz)" ] }, { "cell_type": "markdown", "id": "f589ba52-156b-4e20-b531-783e5fec3b88", "metadata": {}, "source": [ "Compute time derivative $\\dot \\theta$ based on McLachlan's variational principle\n", "\n", "$$\n", " \\textrm{Re} M \\theta = \\textrm{Im} V \n", "$$\n", "\n", "Here $M$ is a matrix and $V$ is a vector\n", "\n", "$$\n", "M = \\frac{\\partial \\langle \\psi | }{\\partial \\theta_{kj}} \\frac{\\partial | \\psi \\rangle }{\\partial \\theta_{k'j'}}\n", "$$\n", "$$\n", " V = \\frac{\\partial \\langle \\psi | }{\\partial \\theta_{kj}} | \\hat H | \\psi \\rangle\n", "$$" ] }, { "cell_type": "code", "execution_count": 5, "id": "3b1b1882-18d4-4448-82a1-7847739fe9a2", "metadata": {}, "outputs": [], "source": [ "# the Hamiltonian in dense matrix format\n", "from tencirchem import get_dense_operator\n", "\n", "h = get_dense_operator(basis_spin, ham_terms_spin)\n", "\n", "\n", "# time derivative for $\\theta_{kj}$ in the scipy solve_ivp format\n", "def scipy_deriv(t, _theta):\n", " return get_deriv(ansatz, jacobian_func, _theta, h)" ] }, { "cell_type": "markdown", "id": "0f93589d-e65b-491e-aa30-50dedc501737", "metadata": {}, "source": [ "## Run simulation\n", "\n", "Use SciPy RK45 solver to solve the $\\theta$ initial value problem (IVP). Calculate $\\langle \\hat \\sigma_z \\rangle$ for plotting." ] }, { "cell_type": "code", "execution_count": 6, "id": "d9ce150b-a716-4b63-bbb3-af8defe96228", "metadata": {}, "outputs": [], "source": [ "from scipy.integrate import solve_ivp\n", "from scipy.linalg import expm\n", "\n", "# time step\n", "tau = 0.1\n", "# initial value\n", "theta = np.zeros(n_layers * len(ham_terms_spin))\n", "# operator to measure Z\n", "z_op = get_dense_operator(basis_spin, Op(\"Z\", \"spin\"))\n", "z_list = [1]\n", "# for reference\n", "z_exact_list = [1]\n", "for n in range(100):\n", " # time evolution\n", " scipy_sol = solve_ivp(scipy_deriv, [n * tau, (n + 1) * tau], theta)\n", " # time evolved parameter\n", " theta = scipy_sol.y[:, -1]\n", " # calculate expectation\n", " state = ansatz(theta)\n", " z_list.append((state.conj().T @ z_op @ state).real)\n", " state_exact = expm(-1j * (n + 1) * tau * h) @ init_circuit.state()\n", " z_exact_list.append((state_exact.conj().T @ z_op @ state_exact).real)" ] }, { "cell_type": "markdown", "id": "8bd1f7a3-0f9c-4e36-9ead-5eb7a5224407", "metadata": {}, "source": [ "## Plot" ] }, { "cell_type": "code", "execution_count": 7, "id": "dd7616e7-993b-4308-85e5-b620c150a046", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAgoAAAEhCAYAAAAeQog9AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/NK7nSAAAACXBIWXMAAA9hAAAPYQGoP6dpAABfQ0lEQVR4nO3dd1xc95no/8/MUIY+FAECIWBQtyqqLrJsCVmOi9xAijbV2RiS3Wy5+W3gKrl3b9sbWcT37u7NJmtwstlNcyRwU9zBtmx1C1CviAEhBAJRRnSYYeb3x/EQIYFEmZkz5Xm/Xv4Dcc73+4w8OvPMtzxfjd1utyOEEEIIMQqt2gEIIYQQwnNJoiCEEEKIMUmiIIQQQogxSaIghBBCiDFJoiCEEEKIMQWoHcBEmUwmysvLKSkpoaysbFz3FBYWYjAYADCbzeTn57swQiGEEMJ3eNWIQlVVFeXl5ZjNZtrb28d1jyNJyM3NJTc3l8zMTPLy8lwcqRBCCOEbNN5YR6G0tJQdO3ZQWVl512ujo6Opra0dHlEA0Gg0eOHLFkIIIdzOq0YUJspkMmE2m0ckCQ7l5eXuD0gIIYTwMl63RmEiTCbTqH9uMBgwm81j3qdJDYfEMELbgojotpCcnExycvK4+7169eqErnfWvdK39C19S983X3v16tXhny0WC21tbZPuW/gxuxcqKSmxZ2Zm3vW6srIy+2gv0Wg02ouKisa8T7M0yc4HG+18sNG+8PUf2i1W64Tie/LJJyd0vbPulb6lb+lb+h5LQkLClPoW/sunpx7GcreFkNOahlg/uBHscDrsKP/1yNsTan/btm2Tjm0q9zrjfrX6ltetzv1q9S2vW537hZgMn17MaDKZyMjIuG3hokajoaysjKysrFHvS0xM5Nq1azywZwcH9J8Q1ZWE+bl/d1b4Hmvz5s3s2bNH7TDcTl63f/HX1+14rgkxUT49omA0GjEYDKOuVRgrSQCG5wB/tubrMKTlRkQjxScPuCxOT+Gv31bkdfsXf33dU1kbIfybVyYKY00dmEwmCgsLR/zZ9u3bR+xwKC0tJTc3947tO/5BLYlPZn7fMrAE8ItxbMX0dv76AJXX7V/89XVLoiAmy6sSBZPJREFBAUVFRVRVVVFQUEBxcfHw78vLyykqKhpxT35+PmazmeLiYoqLizl69Oht19xJUWYeYW/kcH5vHBebOp32WoQQQghv4JVrFFzt1jnMLf/4Ke8du8o3H8rgp99arWJkQggxOf66NkNMnVeNKKjlbx6bjx07v75wiLOtshhICCGE/5BEYRzumzONiKxKOje8y84T76sdjhBCCOE2kiiMg0ajYXFsCgAfmX1/UaMQQgjhIInCOH0t/X4AGvUmugf7VY5GCCGEcA9JFEZx9epVNm/ezKuvvjr8Z1+euxztgB57oIV/O3NYxeiEEGL8Xn31VTZv3jzi3AchJkIShVEkJyezZ8+eEfutA3Q6jINzANh15aBaoQkhxIRs27aNPXv2+GQdhdLSUjZu3IhGoyEnJ4eqqqoRvzeZTCxfvpzo6Ojbauz4k6qqKpYvX05GRsak7pdEYQK+FLcSgOPW0ypHIoQQIjs7m5KSEgC2bt1KZmbmiN8bjUa2b99OSUkJ+fn5k+6noKCAnJycCd93c52fqbY1FZmZmezcuXPS90uiMAF/uehhsENveCtV1xrUDkcIIfyewWAgOzt7zEJ6JpPpjiX7x2Pjxo1s3bp1wveNdh7RZNuaqpiYmEnfG+DEOHze3Jh4Mi6v42pNEMeCuslMVDsiIYQQeXl5bNy4EbPZjMFgGPG7W3+ejMkkGsXFxVRUVDilLbVJojBB34p7kh37T7P31HX+/KF5aocjhBCTZrfb6R0cUjuMYaFBOjQazYTvy8rKwmAwUFxcPGKKobi4mC1btgz/7CjnbzQaKSsrIy8vb3i6orS0lIKCAjIzM8nLy6OsrAxQpjReeOEFzGYzNTU142qrvLycsrKyEecP5efnU1VVNWpbAIWFhRiNRkAZBXG8jvLycgoKCgB45ZVXMJlMmEwm2traRkwn3CmeqZJEYYKyFiex483TfHLmGtYhGwE6mb0RQnin3sEhEl/YrXYYw669soWw4Ml9LOXm5lJUVDQiUbh1hGHHjh3k5eVhNBrJzs4mIyODysrK4ekLxzUxMTFs3bqVXbt2Dc/v5+XljejvTm05Rg1u/sAHxmwrJyeHvLy8Efdt3LiRsrIysrKyhu9pb28fjjMjI2PEuow7xTNV8ik3QcuNMYTObKdl3kF+c+ao2uEIIYRA+eZvMpmGdz6YTKbhb+gOJpNpxGnCRqNxxM8GgwGTyURmZuaIBYCjze/fra2x3NpWVVUV5eXlI6YkjEYj7e3tw+3FxMTcttbCaDRiMpmmHM94yIjCBOm0WkKWXaIl4iwlDYd5frEcEiWE8E6hQTquvbLl7he6SWiQbtL3ZmZmYjQaKSoqoqioiPLycnJzc0dc49ghYTabMZlMtLe3097ePuKaW5OLsYynrfGoqKgYtU/HFIIjObj1GoPBMKI/Z8UzGhlRmIRlYUo9hZO91SpHIoRwJcuQlY+unsNqs6odiktoNBrCggM85r/JrE+4WV5eHrt3K1MpZrP5tt9XVVWRk5PD7t27MRgM404KRjPRtm7+9n+z0eJ0RzwTISMKk7ApaTFvXnuTa0EN2Gw2tFrJt4TwNdc7+7nvrR9zKeEwgZVh3Buwkh8teppHUhaoHZoYQ25uLgUFBRQUFNy2BdFsNrNhwwYqKyuHP0QdH9KjTVPcyWTaqqqqGvXPs7Ky2LFjx21/bjKZxr2N0pmvbTTyCTeK0Uo43yx71lKwaRjS91LRfMW9wQkhXO7ghRYe+K/v0fRpGgAWfQ+fBexl07m/Jbv85+oGN0H+VMLZsZCwtLT0thX/JpMJs9k84kPTMTR/a0XHuxlPWzevIXCsexhNZmYmK1asGLGewNGGY+HiaG4eiXDmaxuNJAqjGK2E883iQsMJ650GwGsmOU1SCF/y7/sv8NiOj2js6GNubDz70v+dvwn6DvEdswB4zfoWu6pv3x/vqXy5hPNo8vLybttVAMoHcn5+PgUFBZSWllJaWkpJScnwNsjy8nJ27tyJyWSioKBg+EO+qqpq+GfHNsW7tQVKouAY4SgvL8doNI7aFkBZWRllZWUUFxdTXFzMrl27hos1jXZPYWEhFRUVFBUVDSdFd4pnrH7HS2O32+0TvsvHbd68mT179tzxmkVv/IjTYUdZ07eOQ0/9yE2RCSFcqdp8nXn7/hydycjWoKf51+fvJ1wfCIDNZiO19D/RYDiHvjeaxkf/jWh9mMoRj994nmtCjEZGFCZpTZRSbOnSoEw9COErvvbZy9iC+wmc0cq/5T4wnCQAaLVaytb/CF1fKP0B3Xz/7bI7tCSE75BEYZK+lfEwoW89g/bDhxmy2dQORwgxRUevXeZI4AEAts/4KoG629d6z4uLZ0fC3xL6zlO89sde9p9vcXeYQridJAqTtColiShLLL0DQ5xruKF2OEKIKfrmoX8FnY3oGyn8cOUjY173gzUP8Y3li7Hb4e93HXNjhEKoQxKFSdJptWQaYwGoMLWpHI0QYirerzvL2RDlQ3/H3D+/65bn/5a9hACdhsNdZ3j9wil3hCiEaiRRmIJpszvpf2AvP7u+S+1QhBBT8BdVr4DWzvQbc8hbfN9dr08whJCyqZq+R97nh2d/44YIhVCPJApTkJSkxZpWx1ndGbVDEUJMUlN3F7X68wC8tPD5cd/31ws2AnBBf4pL5usuiU0ITyCJwhRkG5cD0BvaxrWeTpWjEUJMxmcnWwl96zlmns1i27zxH8v7rYWrCeuKB90Qf3dYRhWF75JEYQqWJcwgoC8MtHZKq4+rHY4QYhJeO1KPti+U3JmPTuisAa1Wy5YoZdHje/17sQz55nkQQkiiMEXTrSkAlDXJgiYhvI25Z5Dyk00APLtq5oTvf/HeZ9AMBjMY2slPjkldBeGbJFEYxd3OerjZfH06AGd7a10dlhDCyQo+L6HzofdIXNLK/BmGCd8fHxZBpmUlAD+7/LaTo3MOfzrrQbiGJAqjuNtZDzdbFaPUf2/UNLo6LCGEk73Z/hlD05vImDs06TZ+tPBpAK5prnGjf8BJkTmPv531IJxPjpmeokdSFvAPp7QMDMCgZYigQJ3aIQkhxuFCewstYXUA/H8LH5t0O09lLMT42y0014bwWWoLTy5PcVKE4m5MJhNFRUUUFxcTExMz4jCompqa4cOYbj6sSUycJApTdN/0dBL/1zfp7rNRs7ZrUsOXvqy9u59fn67kjYbPOTV4kaXB8/jo2b+e0KIxIVzhxePvgtZOWFc8G9PmTrodrVZLzpxl/Kz2Am9XNkii4EZGo5GdO3dSXl7OihUryM/PH/F7s9lMTk6OW2MqLi4mNzfXrX26mkw9TJFOp2VBcgwAZxrM6gbjYd44f5r4sq/wnzr/B59FvkdHXA2fRLzDM398Re3QhODtDuVch3Uhq6bc1uOZMwB473gDg1bZ/eBuMTExo/65wWBg48aNbo3FcTy0L5FEwQkWphgAOH3FrGocnsRiHeL5k//IUFgXWmsg6f1zWDi4FPqD+fDTHn67z6R2iMKPNfd20RpWD8D3Fzw65fbunTMNzaoqGh79DS+f3j/l9tytx9o35n/9Q4PjvrZvaGDS1zqT2WzGZFKeMZmZmZjNZpf1dbPi4mIqKirc0pc7ydSDE2hnNNP7pT/yCw7z3/l/aofjEf7y4xJuGK7AkI4jK37Givg07HY7f/faPoqvN/BX/3aElNhQ1i1IVDtU4Yf+7cxB0NoJ7I1gfersKbcXoNMyfbqWi/p+fnP5U/566UNTD9KNwj9+aszfPRa3incy/2H45/i9W+i1jf4hvy56MXtXvjT8c9q+r9NqGf3QvBWRczi65l8mGfGdOZIEgKysLED5EC8qKqKqqoqSkhKys7PJy8sbnirYuXMnBoMBs9lMcXHx8NqGvLw8MjNHFuIqLi4e8XNubi7l5eWUlZVhMpkoLCwEuG0qxFtJouAEs+OjsNnbaO7zvBXPamg29/HHN20EzlrE5sVGVsSnAaDRaPjJs2tpaz7Aa0fq+cbLn3HmJ08TFhykbsDC75y/0o2uM4k54TOctl7mucQH2NFznBOcxGaz3fVgKeFcFRUVFBYW0tbWRmlpKSUlJSN+n5ubS25uLtHR0RiNRgCWL19OZWXliERgx44d5OXlYTQayc7OJiMjg8rKSgwGA8BwHzt37gSgvLyc0tJSsrOzASVJ8ZUEwUESBSf4UtpCaAZrSDf1nR3MjIxWOyRV/Zc/HKO7G+69sZFd60ce16vVanj5hXt5x/4el1NP8r8+1/Pi2s0qRSr8Vf2pcEJqHuGHuWuc1ubfLn2YHZ/+K5aQLt4yneKZWUuc1rarda9/a8zf6TQjd3K1PLR7zGu1mpHJUd3aX4/72qkabTHjaEpKSsjJyaGsrIyYmJjbRgtMJhPl5eXDCxKNRiPl5eVkZ2djNpspKCigo6Nj+PqioiJWrlzp1NfiaSRRcIL0qBgC+sKxhnTzTu1pvrtkrdohqebgxWb+cLAOjQb+8Rsr0Y3yrUofpGNuWgitQQP8e9P7vIgkCsJ9OvssHKttB+DB+QlOazc+LILpvUaaIqv51wvlXpUohAWEqH6tM+Xl5Y1Yl1BVVTWcEGRlZZGVlUVOTs6oCw8dIxGOdQ7t7e20tyvvl4qKCgwGw/Dows3X+zIZG3OSWIvywDnQclHlSNSVe/Jf6H3kXdZl2ck0xo553f9YogzTNUfWcKjhsrvCE4I9Zy9gCe7BGB/OjNgwp7b9aPRqAA4NHHdqu2JijEbjiJGCWxcYLl++HLPZTFVV1W33VlVVkZOTw+7duzEYDMPTFMCEF0XevFbCm0miMIqJlHB2mBWk1Ik/2e2/pZx7BgY5F3QSW3wLDy2Kv+O1G1LmEdszA7R2flRZ6qYIhYB/bniN3ud2E7Ta+eez/M2iDQB0hzVTY251evuT4e8lnKuqqkZ82JtMJoxGI0VFRbfVWDCbzWzYsIGdO3eSm5uL0WgcTg5MJtOYOyjGSiBGS0S8kSQKo5hICWeHzCjljVhv889/jAD/dOwTbPo+tIN6vrfwobte/9VpmwDYZzvIgEX2ngv3ODekjPrdnzD13Q63WhKfTHTzbALPLuRgdbPT258Mfyjh7JgaGE1BQcFwomA2m9m5c+fw9ENWVtaIao4mkwmz2TwisXC07Ug4srOzh3c1OJSXlwPKSIZjFMGRWPgCWaPgJFnJC/iX45EMtodht9v9svLgr698DFGwaGgxwQGBd73+f6x8kp+W/wZrWBc7K8r5+3unvp/d31htVr716a/Y1fU+mdrFvLvx74gOdu5wui+51HGdnvDrAPz5fNesJfqm7Zu8cqyaU9P6wLfXuKnOUcK5qqoKs9k84gPcUcLZMYJQWFhIUVHRiPvb29spLVVGNB3bIPPz8ykoKBheoFhSUkJBQQFbt24d8XNhYeFwQuHY8WA0GsnNzaWgoICMjAyfqdCosdvtdrWD8DSbN29mz549E7rHYrWR8MJuLEM2zvzfp5gZ518Pa3NfHzEfbcEeNMDPk37EdxeuG9d9K/f8Nyr0h0g2z6dhyz+7OErf8k7dKb5y7CfcCLs2/Gf6PgMlC/+eJzIWqhiZ59p+8E1e7P45+u4Y+p79g0v6eOtoPV/96X4WzIjiyI8fd0kfkzGZ55oQIFMPThMYoGVuUiQAp6903OVq3/NS1UfYgwYI6A/lhQX3j/u+/74om8Cz99Cxby7XO/tdGKFv+XHFezxx/u+4EXYNzWAQ93avR9cbxoB9kK+/VEnJoTq1Q/RI77coc8ZzNc6fdnB4YF4Cdp2Fk0NnudDW4rJ+hHAXSRSc6J4UA3bsVF3xjLlJd/pd4ycALGMZAdrxz2g9nr6IzOsb0d4w8NlZ//t7m4whm43/Xfd70NqZZs7gwNKfc/DZ/8zpB4p56OpWLD2B/Kf/OEpHz+DdG/MzF+zK+oRN01w3dxwbEYzm8Y/o31DGy2f3uqwfIdxFEgUn6k25SM+W3/Mry+/VDsWtOvsstJ2ehq5hBn81Z+LrDB5aoGwt3Xv22l2uFACvH6lH8956Is8sp+rRF7l3hrLjZl7cNMq/+2XuSTFwo9fCT947pnKknqWmo5W+MGVh2vPzxz/qNRkLApQRiw9b5f+B8H5emSgUFhZSXFxMcXHxbatPR1NeXk5OTg7FxcWUl5dTUFAwvIDFmYwxsRBkoVnjX9+M3zvWgL12JosvPMtX5078JL57F8RgTalj94DMn96NdcjGj984hcYSzI+MX2FGZNSI32u1Gn7wzDz6Vx3kxZD/yalW/92Fc6uzlzsJPrCWhJqVzIu98/bdqXp8+nIALmmqXdqPEO7gdYlCYWEhBoNhuG53ZmbmiO0to3EU1sjLyyMvL4+MjIzhVarO9OB05VtEj74N69CQ09v3VOUnmwB4akXKpHZ7LM2Iov/BvbRkHOVIY72zw/MpRQeOU32tk5jwYL77yNxRr3kmM43QhG7sgRaeP1g86jX+6GxdD4G1GTzOYy7v688X3A82DYOhnRy66r+1VYRv8LpEYceOHWzZsmX456ysrNtO8hpNZWUldrudmpoal21ZeTA5A4a0EDDE4Sb/qDZot9t5u+szbJFmHpw/uW9pqZExRPQo0w+/PLfPmeH5lO7Bfv6T+R/o2/QuX388noiQ0begarVafpT2VQAqA49Q2SzJF0CFqQ2A5XeoGOosyRFRRPYoJ6P++8UDLu9PCFfyqkTBUQzj5jrbDo6CF2rSBwai7zcA8FmTf5RyPthwmWtLPqb3ibeYnz75LaFLg+YD8EnHCWeF5nP+cv/vsei70YT38rfrl9/x2h+s2IDhRgrobHzrcNEdr/UHNpuNT2yfMTStmSXpBrf0uThQeU9/3C7vaeHdvC5RGI3jDPE72b17N6WlpZSWllJQUOCC6BTTbMq36oo2/xhu/E31IQAiexNIDIucdDtPJys18msDq7HZbE6Jzde8Zf4UgCeDHiE2JPSO12o0Gv678esAnAo8Rntfj8vj82SHmi7Ttng/fVkfsCBl8u/TiXg6WSnYUxdwSd7Twqt5VaIwlpiYmDuW8DQajWRlZZGdnT18vvitNb5v5jjrwfHfRM58MAYpZVIv9F4Z/wvwYnu/GAFYGDj6fPl4fWvBGhjSMhTSw3um884Izaecbm3kRriyFuSHS54Y1z3fW7qWwN4I7AFW/s8J9Ufc1PRGrVI/Ibx3GtF3SbKc5Zvz1xBy9F4CP9pAfav7EzXHGQ+O//z1rAcxdT5RwvlOSQJwW73tLVu2DB9DOto0huOsh8lYHTWf/ZcuYLPe3q6vsdvt1GprAHhi+p2Hwu/GEBxGXN8MWsPr+c2lAzw+a4EzQvQZ/+fkhwBEdCayMillXPfodFqWsIQK9rO78TP+N0+5MkSPtq/tLISAUZvmtj5jQ8NYbbmPCnMbRy61khYf4ba+QTnj4ebzajZvluPcxeR41YjCzQd13OzWQzxudetWSEdy4IojQLfNfICQvVn0nchwetue5khjPYNhN8Cm4Rvz751ye6tClLLDRztqptyWr3m/4wgA9+snVijoe8bHCTq6mp59i7BY/Xf4+6JVmQpcY5jn1n5Xz44D4PBFzzhJUojJ8LpEwWAwjPoBn5WVNeo9ZrOZnJycEfc41jPcKbmYrNnTlfnPtq4BWrt8uyTxf1w8CEBEbzxJYVF3ufruvj97M6FvPEfvJ5l+/aF2q6aeTq6FKh903503+vt8LF9duIzkpmV0tQXx2Tn/qu/h0G+xYA5RinltTl3q1r4XzQrDMvsCu61vurVfIZzJqxIFgO3bt4/Y4VBaWjpiu6PJZBpRhMlgMFBUVDQiKSguLiY7O3vUaYepCgsOYGZcGPbgfiobfLvS4N72k8DU1yc4PGxMJ5ZYuvutnKz3v/MyxnLwbCv6/euIrV3KE8aJTcnotFqeXD4DgLcq/GPdzK3eqT0NAVY0lkAeSXXOe3W8lhijGFh9iMaUo9R3yntaeCevSxTy8/Mxm83DlRmPHj064ujQ8vLy244S3bJlC4WFhcP/tbW1UVJS4rIYzWs+oifnD7zauNdlfXgC++FM9B9n8e0Zm5zSnlarYUWGssf9WO2d1534kw+qmgm4ksq3w7ei1U78n+yXViRimXWRX2p+Sb/V/85/ePuKsuA2pi+JQJ17l2UtnDadoN4o0MDvLnzu1r6FcBavXMyYn58/5u8cFRtvZjAY7niPsyUFxdIInOny3aJLl693c7XZgl6XwtYFy5zWbtSsNvpiyvl5u4lv81+c1q63sg7ZeP+4slr9ieXjW8R4q4fmTcey9BgD+j5+fuozvr9sYtMX3k5/eRYhpx/jyfvTVOk/dSiNak7wwbUTbMc5SbUQ7uR1IwreYEFEKgBXrE0qR+I6+88rx+dmpscSFuy8fDMxIYChGQ1c1MoWSYDikwe5lnGEiORu1nyxMG6i9IGBzLXcA8B/XP7EmeF5hROmTnSt8TyTOrWdOZN1r0H5uz89IOc+CO8kiYILrI5T1kO0B15XORLX+VnjGwwsrWTuArtT230mfQkAPaGtdPT1ObVtb/SLyx9iWXyCqBW1BOgm/8/1KzPWAXBGd5ohm/+cQ9IzYOVsww0AVrihdPNonvsiQWkPuUrfoEWVGISYCkkUXGB9irJgyhLSRWtvt8rRuMaJgGNYFp4iKtm5H+YrE9LQDuhBZ+P1S1L69pxdKQX+TOKaKbXzV0vWgTWAoeA+3q0764zQvMJbl07Sm3mYsHlXSYpxT6GlWz2WvgDNYBD2QCuvXzquSgxCTIUkCi4wLyZB+bADPqq/oHI0ztfe10tvqLIv/Nn0pU5tW6PRED+YBEBZ4ymntu1tasyt9Icpizq/Pm9qdSoi9Xpie5XdD7tMR6Ycm7d4t/EYlnnnYI7za6aMV4BOR/xACtjhwyvnVItDiMmSRGEUjhLOEyndfKvIQWWY80DzJWeF5THeuHQSdDa0A3pWJKQ6vf0FemXq5li3f8/p/vai8oEe3B3NvNjJncx5s+V65ZCi/Z2np9yWtzjerRTvmh3k/PfpRHwz6MuE7d6Gvdr5tVvuxlHKWUo4i8mSRGEUjhLON5c/naiFQ4sJPD+fgTZ1hjtdyfFNP34wGY1G4/T218YptQLq8c99/w5l15Q6Fen2dKe093TKChjScr27F7vduWtLPFW9rQGA1TGzVY3jidnz0FiCOVx93e1/99u2bWPPnj0kJye7tV/hOyRRcJGtkV8iuGI1vVcNaofidFVffNN3fPN3tucyloE1AEtvAB29vl3d8k7OWpS/5/u+WDU/VV+bu5qYN76CtvwhLl3rckqbnqzfYqErRFlQvGnGQlVjyUyPJUCnoamjjyttvarGIsRESaLgIrOnKwfAVPvgA/mKRvmmv26aaw5uWhSXzPzyFwgp+xKnL99wSR+ers8yiDlQWQeSY1zhlDbD9UGsSpsOwL4vtrf6so+uXISAITTWAB5KUffsldDgAGJX1dG78T1+dr5M1ViEmChJFFxkVmIk9sABzg+asA75zna01p5e+gOUnRzPGJe6rJ/MdKVmQJWfVmg8Xd9J6GtbSPzkWaeWHV47T1nr8Ok535+v/rDhDACRffFur8g4mojEPmwJzXzcKrt5hHeRRMFFkmP09GTv4sbGPXx+rV7tcJzm7OUuwkq/zKxPv8qiONfNeS5LiwGgwuT733xHc/jidTR2LWvj50+qbPNYZs3S0fulP/K7Gf+IzebbB28d76gDIFU7uYqWznZfjLKYtHqoVuVIhJgYSRRcRB8YSHC/cqLi/ibfWb1fVduGBg2rE52zwG4sUTN66HnyDUqSfu7SfjzVoYvK3Pqa2dOc2u4jc4zYoswM6Xv5sN63q19GX1hBaOlWvh71lNqhAPBsunJE+I3QZroG/HftjfA+kii4UNyQ8pCvaq9TNxAnchzW5PjG7yrr0tOxR91gILyDOrN/nbpns9l4PeEV+tfsZ36G3qltG/QhxPQqI0G/v3TYqW17ErvdzonLHWj7Q1iXmqZ2OAA8PGM22sFg0A1JMTHhVSRRcKGZQcrCsQs9DSpH4jyvxfwbfWs/IXmma4etZxumEdinLAh9veaYS/vyNB9duchATDPW9FruS5/u9PaXBjvqKfhuQavmG/20dg2g1Wi4J8WgdjgAaLVa4geUolfvXJVEQXgPSRRcaF64Mjd6dahZ5Uico7ajnZ64BoZSL7MiNdHl/SVYlQ/Jg9d9Z+pmPEpMRwEw9EwnSh/i9PafmqGcPVAfaPLZdQq/rj5I38NlRK00ERKk/kJGh0X6WQBUdvlexVbhuyRRcKEVsWkAw9vcvN1rX3yzD+qLYJbBuXPno5kVrCRaZ3vqXN6XJ9nfoazWXxA4yyXtf33+ahjSMRTcx6dXfa9yKMAnLacYSr5K0HTP2jWTlbAYTb8ec4dvJmjCN0mi4EIPzZgDgEXfRXtfj8rRTN0nzcoHWNKQe1aRL49W9r432Bvd0p+nqNPUAbAhfrFL2jfoQ4joUbZJ7qk77pI+1Hauvw6Ae0Jdu+h2ol6Yv5aw17Yy+NlKWm7I6ajCO0iiMApnnPUAMC86nvDqhQRVrKKmxfsLL53qU+rmLw51zTfdWz2cpMyld4e0MjTkH9/A2vq66QtRFm8+4+QDt262mEUE1Bq5fs03HwFNWiW5fCB+nsqRjBQdpmd+kgGAozVtbulTznoQU+WbT4kpcsZZD6AsXlpxfRNBFxbQ2DLopOjU06xTHr4PuunhuyFlNrr2OHRXZnKm2Temb+5mj+kUaO1o+0NYmuC6OhXfTXgW/YEHuX4xymV9qKWx6waDoWYAnkhbpG4wo1g5SykmdqDGPSNlctaDmCpJFFzMUcr5YpN3jyg0dncyGKqUU3bXw1cfEMSaM19Df3Atl5sG3NKn2s40taK5EcW0gSSXHLjlsNyonG56vLYdq4+N1rxddxo0oOsPY35sgtrh3CYkvZmeZ3ZTrHlF7VCEGBdJFFxsZqKeoZhW9refUTuUKTlYV4+2PZbA3gjmxrjv4btghgGAs1fMbutTTZYrCYT98RlesL3g0n7mTI8kXK+jW9/Ggcu+dUrn/mZlR0HsoOclCQD3z0zFHtZLe0gjg0MWtcMR4q4kUXCxG7F19D32Nnsj3lE7lClpvxZE6HtP8kTN99za7z0pBuwaG0ebfKcWxZ2cvKysT1iaGuvSfrRaDfasz+jd/Ca/NH3k0r7c7bL5BgwGYgz0jNLNt3o8fQEaawD2QAtvm86qHY4QdyWJgoutjldWXXcHt3v1nvXTX3yjXzwz2q39WuOb6Pnyb3k34ddu7VcNg1Yrp64oicLiVNf/Pc/RzwTgqI/t6Q+5sIiw3X9GXsxzaocyqiBdIDF9SQD8sf64usEIMQ6SKLjYQzNmA2ALGqDa7L0L8k7WKyu0F6W4N1FYm2wEnY2+kA66B3x7ncJ7dWdpe/bXDGaVk5EQ4fL+Hpp2DwCXuezyvtzFbrdztsGMBg3LUuLVDmdM8wKVrb9Hbvj2eRvCN0ii4GJxoeEE9IUDsLfhosrRTI5lyMpnmT+n9/G3mDbdvaMiKxJmorEEgs5GWb1vffO91YcNZyDQSmioMjXgatkZyiFFfaHtNPV0urw/d2jq6KOjZxCdVsPcJM/d0fFwvJKk1VGnbiBCjIMkCm5gsCjboSpavfN42X0NtdiDB7BF3mBV0gy39q3VaokcUKpAftLo29++jt5QSlWn69wzt744PomA3gjQQOmlKrf06Wq/rtlPz+bXCX7gGPogndrhjGnrrBUA9IW10dBpVjcYIe5CEgU3SNYpq6/Pdnvn6vIPryg7NsL74tAHBrq9/xSNMp97zGxye9/uVGOtByAzyj0FrQCmW5XEr6zJNw6IOth6EXtkJyHRnl23ZGFcEhEtaQSen8/ntS1qhyPEHUmi4AZzwpSH8WVLk8qRTM6RDuWb7gyNOgVbFoYrC0JrLN6ZaI2HzWajI/gaABuT7nFbv0tClTU0J3p84+Ctc711AMwJnqluIOOQ3fkNgitXc6nOs5MaISRRcINHYjMJqlxJWLX7PgCcqXpA+aa7MEyduvn3xSsfZq0BvnEK52gONtZiDxqEIS2PGxe4rd9nk1cTeHIJmnNz3danKzVplGR8Zaz7RmUmy1Gh8WiN9y5yFv5BEoVROOusB4esmfMJOncPrRcNDHnhFsmWAOXhuzZenQ+Tx1IXEnA5DU11BuYe39z58E69MvQf1htHeJDebf0+O3sJ+lPLaLsY7fWHFPVbLPSEKKdFbkyer3I0d7ciIxZ7gIV9HaddunVaznoQUyWJwiicddaDQ0psKMGBWgatNupbe53SprvU3zBjCVVWxD+uUt38DEMcGWceI+jUEs5d9Y3V+bdqbhlC15DCbNsct/YbFRrE7MRIAI7Xdbi1b2f76MpF0A2hsQbwQLJR7XDuat6McHqy/0Dz2j0cuOq6hc5y1oOYKkkU3ECn1TI91Yp1Zh2fNnjXXHDVlWYC6tIIbUsiwxCnWhwLZihb3c41mFWLwZW6TLGE7N1AXtRWt/c9Nz0Ia/IV3qivcHvfzrS3Udk+G94XR4DOc3c8OEQE64noV/5NvWbyjV0nwjdJouAm5gWf0//gXt5t/VztUCbkWiPo9z/E49eeVzWO2dPDsYV1cajlkqpxuMqpejPgnoqMt+pLraH/4Y94e/BDt/ftTA1tfWjbY0glVe1Qxm2WVln3s79dSjkLzyWJgpukBk0H4GKPd51Z4CjdvCjFoGocTQmn6H3mNf4YsEfVOFyhuauHK72t2LFzjwp/z+sSlPn85gDv3JUzrC6F0Hc3893wr6gdybg9EKMsXK0e8u2tv8K7SaLgJvPDlS2SV4e8a+V+5bUr2LGz0M1nPNxq5TRlzrldd13VOFzhjdpj9Gbvxrr5baJCg9ze/9PGxQBYQroweXGZcccJo44TR73Bc8blAHSGtdDZ369yNEKMThIFN1kZp3zQmQO950FsHRriyNJf0bP1d0QmqLvbYP0MZceFJaSL9j7vWhB6NwdblHUrBtQpOZxmiCGwR1nQ+FbtCVVimKrO/gFMrcpCVzVGZSZrbVIG2gE96IYovXRc7XCEGJUkCm6y7ovDoawh3bT2dqsczfgcuFqLPdACWhv3Jak777sgJhHNYBBo7Xxc751nZozldLey4j0t0L3lsW+WOKSsiN/bfE61GKbi9ZpjdG/9LdZHypkW6b7tpVOl1WpJGFT+v7/bcFzdYIQYgyQKbjLbEId2MBiATxu8Y0HeR1eVsxVC+2IICXT/kPjNNBoNkQPKCvEDzd7x9zde9UONACyJVKegFcACvdL36V7vnCvf13wRdDbC9Z6/2+FWj4Y8SNDRVQzWJaodihCjkkTBTbRaLeEDMQAcbvGOh/Hn7coH8nSmqxyJYrpGeZCevFGnbiBO1hGk1Pp/MMG9NRRutnbaPACuar1rsa3DiU7l31RagHqjMpP17Yz1BF1YwNkzNux2u9rhCHGbSScKr7/+Oq+//rozY/F59/asR7/3YULaPeOD924u9iulm+eFeMZ2s7mhyoeAacB3KsydbbuGLbgf7PBImnrVBJ8zZhJ8YC0BHz9E74BVtTgm6/KQ8p5YHOn5hZZutSw9hkCdluYb/dS39qgdjhC3mVCicPz4cb7zne+wcuVKTCYTNTU1rFixgu9+97scP37cRSG6n7NLODtsMCwnoCGVa16yC+3aF3XzV3lI3fwNcUsJPLWYoLoMtUNxmrJ6ZU1AUF8U00LDVYtjbtw0ktvvQWM2DG+J9SYdgcqozNqE2SpHMnEhQQHMmavBklHN784fdXr7UsJZTFXA3S7o7OykqKiIXbt2sXLlSvLy8li2bNnw73/wgx9w7NgxXn75ZSorK9m6dSu5ublERka6NHBXcpRwdjZHqdxL1zy/DHHXQD99oUrd/E0p7juk6E6eTF3KD0/U0xqgxTpkI0Dn/TNn7dc1BJ5bwKzYWFXj0Gg0LEmNpuxkEycvd7BqlnpVOCfqfHszQ3rlnIpNqZ7xXp0oy/yzDIQcorTVzg952Kltb9u2jW3btrF582antiv8x5hP2tdee41HHnmEnJwcoqOjqaio4F//9V9HJAkOy5Yt4+WXX+bo0aNERUWxfv16Nm3aJFMTt0iJD8aacpmq8AMePxd5trGdwDOL0DekszLBM47snRkbhj5Qx6DVxmUfGaK90RhKcOUqcoKfVDsUZqbB4Lwz/L7Vuyo0OkZlAnsjmR7unV9QHopdCED1UI3KkQhxu1EThS1btlBRUUFRUREffPAB3/72t8fd4AsvvEBFRQUvv/wyn3/+OVu3ur92vadKjQ+jf90n3Fh0mIsdnl046HLTIMEnMlnbmI1W6xnf3LVaDakztVgTG/nsim88UM823ABgfrI6NRRuFpLYyeCKo3wedEDtUCak6bqFgLo0UvvUWww6VV+etRKA7vAWmru7VI5GiJFGnXrYvXv3lBtOT0/nxRdfnHI7viQuNJyAvnCsId18erWauTHxaoc0pjNfHL7kOIzJU7QuOkB/+Fneuh7C8yxVO5wpsdlsnOi7gD0ozCOqCT6etpgXz0JPSBvdg/1uPe56KnqvRqHf/xBfedI7px0AViSkEtAfhlXfwx8uVvI3mQ+pHZIQw5zyVbGurs4ZzfgFg0WZ+61odd2xss5wsPUCtpBeFnjAN92bGYOUwkAXe71zG9/NTlxvpGXdHnqy/8DMePU/lO9LSkU7EAxaO2/XnlY7nHFzjMosSDaoG8gUaDQaZliV3UXvNx5XNxghbuGURCE6Onp418OxY8fYvn27M5r1Scm6BADOdl1ROZI7+zSphN7ndjMw7ZraoYywJEp5mDbaPSuuySi7opwYGNxnIEofonI0Sq2P6AGlVkV5o3ecZmiz2TjZXo8du8eNfk3UinClTPnxft+qPCq83113PYxHbW0t+fn5/Of//J9Zv349AJs2beKDDz5wRvO3KSwsxGAwAGA2m8nPz3fJPa4wJ2wGJzhCvcVz90g2dndiCVV2Zjyaqt7e/tHclzCLf6qDruBW7HY7Go1G7ZAm7XBrDWghfihB7VCGGXUzaeMyxzq9o/plZcsVmh/dhWYgiPSEHLXDmZInUzIpvfIGLUENDA3Z0PnArh7hG5zyTtywYQOFhYUsX76curo6KioqqKlxzWIzxwd+bm4uubm5ZGZmkpeX5/R7XGVptPKNuFXjuYsZ369Tvk3q+kOZHe1Z6yjWpyjfumzBA5xva1E5mqk513sZgFn6FJUj+ZPlBqVmRu2QZ494OXx4RdnxEGwJJyxY3TLjU5U9aylR+zYS8senudjk+Vuohf9wSqKwY8cOKioqiIqKoqamhqqqKiorK53R9Kh9bdmyZfjnrKwsiouLnX6Pq6ydrhSE6QvpYHDIokoMd7Ov+QIAMRbPShIAYkPCCOyLAOCjhgsqRzM1V+3KqNJyg+cUkNqYfA8A5uAWbDabytHc3ZFWZeQj0eYd1U7vJDQwmPtDM9EM6Dl00XO/SAj/45REITc3l5wcZdhvw4YNGI1GOjo6nNH0CCaTCbPZPDyFcLPy8nKn3eNKq6fPJOLgw4S89yT11z3zuORTnXWA59bNj7FOA+Dz6967RdJms9GtV44cfyhprsrR/MmjafOI+PAJQku30NDep3Y4d3W+TxmVme1BozJTcf9c5b29/7x3j5YJ3+K0SbCoqD8tJPrBD37gskRhNAaDAbPZ7LR7XClIF8g9liXozNFcuuaZx03XDSk7CpaqeJrhnay1P0DwofsJue6Zicx4HG2uV47wtml4eIbnlB0ODQpioX4WmqFATl52/r9hZ2tEGZVZEeM5ozJTsXhOGAOLj1ES+luPL8om/MeEFjPGxsbS1tY2rmtHq+DoKjExMbS3tzvtHsdZDw6OEqjOMmd6JKevmLnY1MmjS5Od1q4z2O324dMM1yZ6zjfdm30pdg3v12i5HuaUtbiqaGyxEPT5amLjlA9nT7JoZjSn6s2cqu/gieWem4wNDlnoCVGeRxuS56kcjXOsMk7DsvAkFq2dfVdMPDhz8gnQq6++OuK8GjnrQUzWhJ60HR0ddHZ2etw5DhNNEu52j6vOenCITu5ncP5p9pg7+Ws8a1dB841eAqsysUV38KW1nlnAZs50ZY1CtRcv+Gq6NkTQxfk8FOF5H8TTZvYzsPIwv+67wHYWqR3OmD67agLdEFh1rJ3hfadGjiY+NJKo3kRuhDfx2+rDU0oUbv2CI2c9iMma8Fey2tpaysrKaGtrY+XKlWRlZbktcTAaR38YmM3mMX83mXtcLq6dweUVVHY1Ad9RJ4YxnGvoItA0m4yECOJC1DvN8E4yEsIZir9GddQNOvs3EalXv1jRRF1oVJKcuUmet/c/KTEAi/Y8l/rD1A7lji439RN4ehHTogMI0gWqHY7TLAqcy36a+LTjJPAVtcMRYuJrFDIzM/nwww+prKwkPz+f6OhofvjDH7oittsYjUYMBsOo6w6ysrKcdo+r3RuvbEHrDBrfNI47OUo335NiUDWOO0mICqH/4Y8YWH2IT696x37/W+3vPs5QbCvp09UvtHSrJ9KVUQSrvodLZs9dfX+9SUvw8eU8Pvi02qE41ZNJywGo1dXIOgXhESacKLz88st8+OGHfPjhh1y6dInq6mra29vZtGmTK+K7zfbt20fsVigtLSU3N3f4Z5PJRGFh4YTucbcNKcrhNbbgfi60e9bq5g9bqhiKv8asGcFqhzImrVZL2EA0AAeveefOh+Npb9P3pbcJiPW8A4BmRhkI6lFGOt6uO6lyNGM7c8UMwMKZBlXjcLavz7sXbBosYZ0cbLisdjhCTDxRuPU0SKPRyMsvv0xubi6/+MUvnBbYWPLz8zGbzRQXF1NcXMzRo0cpKioa/n15efmIn8dzj7vFhYYP1wL4pMGzyrXuDXufvkfepzfBswvuJKDUeDhh9r4H6YX2FmzB/QBkpXrmgtH4oSQAPms+r3IkYzvSfRabvpd7vLx0860SwyKJ7FWqdf66+pDK0QgxwTUKWVlZfPTRRzzzzDO3/e65557jJz/5idMCu5M7lV92VF+cyD1qiLbE0RLSxeetNXyHB9QOBwDr0NDw3v71SZ61yPJW6cHJ1HCamn7vW8n90RXlwzewN4L4UM9cBzIvOJUGznGqxzMPL2vt66Y68zVYDgmJj6sdjtMtDJzLwYF2Trd6/5kmwvtNaETh5Zdf5gc/+AEnTpwY9ffR0dFOCcofpAQoleTOdtWrHMmf7Ltqwh5ohSEtD6fMUjucO1r8xeFQ1+zNKkcycUdalfUy0V8UjvJE98Yp02NXNZ6ZiL1fdw40oO0P4Z64RLXDcbofpX6FsJJttH+eLusUhOomlCg4phmWLVvGX/zFX4xIGOrq6rh0yTsXlqlhfvhMAC5bG1WO5E/Kryp188P6YtEHeNbe/lutiVd2rHQHt3ndg/RslzJdMjPAc8sOPz5TWdA4YB+gd3BQ5Whut7dJea/GDMZ79cFgY1k3dwaBOi0N7b1cbu1ROxzh5ya8RiErK4v29nauX7/OsmXL0Ol06HQ6cnJy3Lb7wRdkJz5AyLtPEnd0g9qhDDvapiR60/HcDzCHh4cXhA5wsd1zV+aPpv6L5HBBRKrKkYxt5fQUEt/+GqFvZlNzzfM+qI53KqMynlpmfKrCggNYbozFjp2yM54z6ij806RKOBsMBkpKSrDZbFy6dIlLly5x9OhRjyvE5MnWzEhB1x5LfZOFAcuQ2uEAcHFAeSDND0lTN5BxiAsJJ+nkBkI+fJSm6555uNZY2gOVxGZ1nOeWHdZqtSxNVBLGEx5YyvlPZcY99+9wquIWNdP7dCk7W3+ldijCz035rIf09HTS0z3zTIDJcpRwvrn8qbPFR+mJCg3EZrdjavaMLXLNGmXhlCd/gN1spW0lupZELjf3qx3KuHX3Wwg4eC9BFSvZNNOzF4wuSVXWHJ3ysETh5jLj6zy0zLgzPGCcgT28h9qgiwxarZNu59VXX2Xz5s1SwllM2qiJwksvvcTx48en3Pjx48d56aWXptyOuzlKODvzfIdbaTQaohY20b/qEK/XH3VZP+PVPziE7uAagj9fw1OpS9UOZ1xmT1dGsLyplHN1UxcBTclMv5pJRnSc2uHckSGlh771H/KrQNdve56I062NyvZSm4YvpXlmmXFn+NaCe9FYArHp+/jd2cpJt7Nt2zb27NlDcrJnnSsjvMeoicKGDRv48Y9/zOzZs3nppZfo7Bz/g7izs5Of/OQnrFixgh07dqhW/dAbWJMasM65wL6O02qHwsWmTjQt8cRfXcw907xjFXlswhCWjGrKeverHcq4XWi8AcDcJM+fprsnKZahpEZaIy5jGZr8N1pnq2saIPjQ/SSaVhMb6tllpqciNEBPyoAyuvebOu95jwvfM2qisGzZMnbv3k11dTVRUVGsX7+eTZs28frrr4/Z0GuvvcYjjzzChg0biI6OpqKigl27drF06VJXxe71ZumVDP9Sf4PKkfypdPOCGQavWUWuibnBwL0HOGk4oHYo4/ZecyWWVBNJKZ6/U2Nj6hywBmAPsPJxQ7Xa4Qy73DhIYM1s1ls3qh2Kyz0So5Rzrhjw3AqZwvfddY3CCy+8QEVFBS+//DKff/45s2bN4rvf/S7Hjx/n+PHjfOc732H27NlUVFRQVFTE0aNH+fa3v+2O2L3e0ug0AJo16pdxfuvaISzGamakecbCyvFYO302AP0hZvot3rGg8WPbpwys/Yyu+NvPHvE0+sBAInuVCpjvX1F/1MvhtKN0swefR+Isf3nPegC6Iq5xukWKLwl1jHsxY3p6Oi+++CKXLl0iOzubH//4x+zYsYOcnByqq6vZsWOHzy1qdLUHEpUPup7gdmw2m6qx7LXvY+C+A/QnqD+6MV6rEmeCVQc6GweuemYFwVtd1ylJ4YpY7zgWOVWbAsCRjgsqR/Ine/sOMzStmVnJvjvt4LB0WgphPXGgtfOz03vVDkf4qQkfMw3KGoYNGzxn/7+3Wpc8C85rsQdaONR0mfuT1Uu02gKVCodrE7xnFXmATkdofzS94a3sv3aJDWlz1A7pjnoHBxkIMQPwULJnx+qwLGIWp+xHuDToGWdq9FsHuTD7fZhnIy7xCbXDcYt1AfdRdr6ay8E6WK92NMIfTXl7pJi80KAgQvpiAChvOKtaHHU32rGGdAPwaKp3rSKPsytlkI931KkbyDjsbbgEWjsaa4AyGuIFNiQp74e24GseUQHz4yvVoLOhsQSyNtk7RmWm6n/es43gitUcr7Rjsao78ij8kyQKKptuV3YYHG9Tr/raH02nAOWQotTIGNXimIy0IOWUw+o+z98jvu+asiAwrC8WrdY7/uk9aVyIpk+Ppj2aS+3taofDx1+UGY/ojyNAp1M5GvdYlhZDXEQwXf1WDld7VxVS4Ru842nlw7YFPkfYrj8jrn6JajE4jhKOH/L80s23WhipfDNvtHn+Qi/HqEciCeoGMgHRISEsO/QCIeWPYmpQv7DV0Q6lzPgMjf/UBNBqNWxYHI81oYl/Of2J2uEIPySJgsrunZGKxhLE2YYbqsVwskdZgT9Pn6ZaDJP1dNJq9OWPEHZondqh3NWl/isAzNJ71/kEy9JiATjpARUaLwzWAbAk3D+mHRwiFzXRv/ED3gh4DeuQ9+xMEr5BEgWV3fPFFq8LTTewDqkz/3gVZdh+dYx3LLC72ZoZKQRcS8J8PYCOHs875fBmMRdWof9kA0/G3q92KBOy+ItSzpWX1R/2bg1sAuDhRO9aSzNVP1r2GBprAJawGxSfOKx2OMLPSKIwCnec9eAwMzYM+/JjmB98T1mo5WaD1iECyjeg/2gjW9JXur3/qYoICWR6dAjg2aWc7XY79XUQcDWF9SnelZBFJ/XT89RrvDbzp6rGcbq1kSF9H9g0PJWxSNVY3C02JJw5g/cA8DPTexO6V856EFMlicIo3HHWg4NWq0E7o5mh6U18dPWMy/u71YXGToZ6gom7kcbiBO8o3Xyr6DltDCypYk/D5Ovhu9rV9l66+60E6DRkJESoHc6EPJiWhj28C2tIN5fM6o0q1F+1EPLBl5h5bgPxYd71d+gMeemPAnAu+ATmvvGvF5GzHsRUSaLgAZI1ysr9ig73V+s7Va/MOy+cGe01pZtv1ZdUi2XRST7tPKZ2KGN6v/4MgwtPEDe3g8AA7/pnlxoVTVCvAYC3ak+oFse5+m501xN4OOg+1WJQ0/cWriNgIBR78AA7jn6gdjjCj3jXE8tHzQ9NBeDSgPu3SP62uYyBJVXEG3vd3rezzAlVFgdeHmxSOZKxfXj9GINLj9E/65zaoUxKolVJZj9tVq/eh2Mx5ZI079rC6yyBugBWaVYA8LtrH6scjfAnkih4gDVxyglxzdpmt/f9ub0Cy6KTBMart+tiqpbFpAFwXav+mRljOd+j7Hhw1H3wNgtDlffoyd4a1WL4QPceFmM1c1JCVItBbfn3KNUom2iiydytcjTCX0ii4AEemams4B4IMWPu63NbvzabjY5gpf7AhqT5buvX2W4+HGrA6pmHQ10dUkY7FkWmqRvIJG1IWAhAo+6KKv1f7mynOeMoA/cdYG6K/61PcNictoTlJ79KyNtP8fph7zmXRXg3SRQ8wOK4JLSDwaC182H9ebf1W9l8BVvwANg0PJZ2j9v6dbZ7p6fBkBZ0Qxxq9IwzCW51I6gVgHvjM1SOZHKem70U7GAJ6eZ8h/tHvvbUKscsB/ZGMis6zu39ewqNRkPe0jVo0PBy2QUp6SzcQhIFD6DVaokamAb9wRy/1ui2ft+5rJRuDumLIUrvvcO5gboAQvqVvf77mty/xfRuam+0K9v6gA0p3nPo1s1So6KJas4goHoOx1Sop/Bps7K2I97qfdVDne3P1hqZFqmn1tzBi5/tVzsc4QckUfAQX+n5c8JKv0zANfdtYTrcpnyoJtu9c978ZtNsyuFQxzzwcKjyL0aJAvrCSIqIUjmayXum+6voj9xHQ737d8ec7P6iemiwHGUfFhzA008G0/N0Kf/Q8TN6Bz270JjwfpIoeIglyfFo0HC2wey2Ps/11wFwT6j3l8N9RrOZ0DeeY3rjUrVDuc2R68oCQIN1msqRTE1murLboNLU5va+GzTK2oj74ryrWJWr5D/wAFotDEZ08L29u9UOR/g4SRQ8xPwZyjfNc1fdt/ugxa4MIT8wzTuHw292b2IG2p4Iapp61A7lNnHN8wj949M8Nvik2qFMyXJjLHbtEAfbzrv1yOmO/h76QpStkU+mqXd4midJDItic9AmAH7b/RZd/QMqRyR8mSQKo3BnCWeHjOmh9D1UzrkH/oOGLtcnC119FnR7Hif0zWfZkuF9pZtvNTsxEoCLHljG+VJjD9obBtbGz1M7lCm5JzWKnuw/cGVtKUeu1bqt33drz4LWjrY/hBWJ3nWgliv9cu3z6Ab1WMJvkPvx78a8Tko4i6mSRGEU7izh7JAYGY4m1ow9vJs3alxfYfB4XTsaNMwMSmSmwXvnzR2MieEMLD7GlcXvU3/DrHY4IzjOoJibFKlyJFMTpQ8mrF+Zfni91n1VMC3N0YS+toV767K9tnqoK8Tow9kS+jgAf7C/xk8rRl/YuG3bNn75291SwllMmiQKHiTRonxb+uia6898qPhinnlFRqzL+3IHQ2gwtjkXsabV8UnDRbXDGWbu6+NcxvsMLjiFMTFU7XCmLE2jVBHd3+a+CpMn6jrQ9oWybtpCt/XpLX754DdJ7E2DQAt/07ST107f/v/lZ8f2MesP290fnPAZkih4kIUhyqLCUz2ur35X1P17+h/YS1yG71R3i7Io++sPt6hXPfBWH1+5iCXjEpaFp5geFaZ2OFO2KkpZz3LB4r5zSSpMSg2K5Ub/LN18JyEBwZx59P8S3ZtEQM0svvf/TvLHyis0dfQyYLGy8a1/4nvN/0BnmvsPnBO+QxIFD7IuXqnQeFXr+oprNfqzWNPqSE/Su7wvd0nSJgBwpsv9Z2aMZV+zMroRORCHVuv9/9wen6EsJuzQX8MyZHV5f009nRw27mZg8TGWZ0S7vD9vFKMP59Kmn/Fwz2N09w/xZ/+8j5m/KyDk3ecoD3kXtHZm9/rXsdzCubz/yeVDns1YCsBAmJn6zg6X9XOypRFrSDfYIWdWpsv6cTfH4VB1g+4rWnU3x8zKor9kjW8UCnosYx5YArAHWCm74voqoruqK7AmXcU+y8TMWP8t3Xw3MSFhlHz/Ib661ogxPhx7lBl7SB8M6fhuyNe58NRLaocovFiA2gGIP5kTE09gbwSW0C7euHSCv8l8yCX9lNZUARDSG8v0cO9fyOiwIsbI6x1wXYXDtcZiGmiAYJgXOlPtUJwiJDAIQ+90zFFX+GP9cR5Lc+26gQ8aT0AAJFtnykLGu4gMCeRfX1gDQG3XvexvusSq+HTmGhJVjkx4O0kUPMwMSzqXrzdxccjssj72Xj8Deki1p7qsDzWsnzEPOpTDoXoGBwgLClY7JFp0StKyOs47z3gYzb32NXxUFY8l1fULYU/0V0M4LA/z7q2l7pYeMY30CO8u8CU8h0w9eJi/0D9P6AeP02Ny3cE35waVxX4rI32ryt2KhBQ0ViX3PdBQp24wQNdAPwMhZuCLJMZHPJ+ykaCzi7h4zrXf8G02G81BSkXGx2ZIoSUh1CKJgodZlqas7D5e1+6S9q1DQ7TplTn8x1N86+Gr0+pYfeJ5wv7wVXpb1T/k6tMrJtDY0VgCyYz3nUJB982NB+DE5Q46+1x3rPfeq9XYggbAqiN7tm+9V4XwJpIoeJilXyQKl1o7aOlyfjniisYmNF0RaAYDeTLd91ZCL4lLQWPXcr7RfaWwx9LfGkrYH77KqjNf84kdDw7JMaEkz7AzkFrD789WuKyf12orAYjsTfTq002F8Ha+8/TyEbERwbDpY3q2/p7fVR9xevu19RZC33uSh4//FaFBQU5vX23zvqh+eN6NZ2aM5fzVTjRDAayI9b0TD3WZZxl44DP+o7HcZX2cu94Cg4HM0Xn/oWVCeDNJFEahxlkPN4vRh4HWzsfNzi+SUlmjVGRcme6bC50iE/rpv28fb0eqf6Ke44Cvecm+s7PE4cEYZTTqjNV1VTBtp+YRVrKNv07IcVkf/kDOehBTJYnCKNQ46+FmwxUae51fYfDzWmUV/soM1y2WVJMxMQKrsYa22GoGrK6bPx+P96N207/6IDGJri9M5G5fmb0agK6QFlp6nX8QV9+glROXO9DYtayd7TvrO9Swbds29uzZI2c9iEmTRMEDPZxwDwBNOudWaOwc6GPfip/R+6U/kpHqe9MOAPdPT4chHQQMceBqnWpxdA30055wEevsi8xO9L1CQffPSCOgJxK0dn5z0flTZFW1bViH7MRH6UmN8/7S10J4M0kUPNAzxqUADIZ2UtPR6rR2Sy4eg0ArhPewOCHBae16kkBdAKF9SqnfvY0XVIvj4yvVoFV2PCyPT1EtDlfRaDSkWpW1F+80Vjm9/R21f6Bn82sYVtRJoSUhVCaJggfKiI4jqFeZ1/7NBed9Wyut/xyApMFUn1qFf6tElEp0lR3uO7joVp81KXP3Ef2+ccbDaO6PUqoynuh3fkJ2rPcC9sguUqb5zlkkQngr33yC+YBZttkAvHPNedvPjvYpiyPXRi52WpueaFaw8g2+uu+KajEcMytJSpKPnPEwmi3GVQC0hzbSPdjvtHZ7rf006esAyJ65ymntCiEmRxIFD7V52hoCao101zqnTG57Xw9tYcqah6/NesApbXqqpdFpADRxTbUYagaUv2tfOeNhNJvS5xBzOIvQN7K52OC8mh//fu4I9gAr2v4Qtsxb6rR2hRCT43VnPRQWFmIwGAAwm83k5+ff8fry8nKKiorYuHEjRqORsrIyVq5cSXZ2thuinby/XbSRf3mlmwagrWtAqa8wBb86exh0Q+j6w9iUOtc5QXqoddPnUHhBQ79lCJvNjlbr/jnuli8OploV67s1AAJ0Oh4OW8UHfY0cvHidTKNzktrd9QchEFIHZhEU4HWPKCF8jleNKDiShNzcXHJzc8nMzCQvL++O95jNZqqqqsjLyyMvL4+MjAyPTxIAEgwhw/vv959vmXJ7b1xV1iekWzN8ds7cIStlLoaSr6N/9wmutDm/uuXd9FssDNqVrZkbZsx3e//udP88pZzz3jPOG72pHFSmyLJifOcIdCG8mVd9YuzYsYMtW7YM/5yVlUVxcfFd76usrMRut1NTU0Nubq4rQ3SqB+fHMxTVwe9qDk65LXNtBAGX0/hSzGonRObZggICmJNgAFCllLOpuYfQt54j8a1vsCLB93Y83GzDwgQG559hT/yvuNI99fNJqs3NdIcpiXHu/HVTbk8IMXVekyiYTCbMZvPwtMPNystdV0ZWTcEZTfQ9+RZ/DHpjSu2Yewa5UhWLft9DFCx+3EnRebZ5ScpozPmrzi8GdDdnG8wAzI/33R0PDotmxqCZbcKS0MQ/nZz6v8N9F68ReGEeUS1GViRJoSUhPIHXTACaTKNvdTMYDJjN5jveu3v3bmJilMOWjh49ys6dO+94vaOEs8O2bdtUqdL4rfn3sePov9Af3s6plkYWxSdNqp3951uw2e3MSowgOSbUyVF6pqGUenqj3+FXvRf5G/67W/s+VW8GYOHMaLf2qwaNRsMy7WIO8wlvtRzi/7Dl7jfdwcmzAwQfXcMLG33rCHQ1vPrqqyPK0EsJZzFZXpMojCUmJob29rGHPI1GI5mZmRiNyqKy9vZ2cnJyKCkpGfMeRwlntc2KiSOsO56eiBZ+eX4//xQ/uYfw76oPY4s0s27BLCdH6LmmxQRgC7lObaf73+L/bvs1fRvaCEny/UQB4Otp6zjc+gmmoIv0WQcICZj8wtuPTytrHdYv9N1tpe5y6xecm7/8CDERqiUKpaWl7Nq1667Xbd++nczMsRc13SlJAG67d8uWLeTl5Y05jeFpFgTM4SgtlLUeh0l+W3sz6DX6N7cREja5EQlvtD55Pv90CbpCWrEMWQnUue+t3hhSy1BMD+kJvle6eTTPL1rJ994NxRbSyytn9/PXizdMqp0DV+q4ZK8hOCCBtfPjnRylEGKyVEsUsrOzJ7T7wDEicCuz2Tzm70BJSG7ux5EcmEymOyYgnuKJxBUcvbGfS9rqSd1/sqWR/vA2sMO3Ftzn5Og818aZc+GCFnughc+umNiQ5p6h7Or261hDlJ0WT6YvckufatMHBjJrcD4XQyr57eXPJp0o/OP5d+nb9B6RHXMI1/+Zk6MUQkyW16y0MhqNGAyGUdcqZGVljXqP2WwmJydnxD2O9Qx3Si48ybcW3A82DYOhNzjUWDfh+//pVBkAYT3xzI72n29p+oAgwvuUEzI/bHD+cd1jebvuFACBvZEkR/je8dJjeTZBSUJP2E9ht9sn1ca+rmMArIm4x2lxCSGmzmsSBVCmIW7e4VBaWjpiu6PJZKKwsHD4Z4PBQFFR0YikoLi4mOzsbK+YdgCYERFFVI8yX/uPJz+c8P1vdnwKwDq9/5XCnaFRjtU90jG50ZjJ2NeinHuQMORfc+x/vexhNP3BDLUaONk08bofx1uv0BJ2GYDvzFnv7PCEEFPgVYsZ8/PzKSwsHK6dUFNTQ1FR0fDvHVUYb67WuGXLlhHJQ1tb2x0XMnqi5yOf4eXyi9QHJGHfZB/3aXqfXamhI7IB7PA/Mp92bZAeaFFYOuc5RvVAvdv6PNNTCxEwJzjVbX16gukR4Ww681fsP9fKgeQOliRN7HTSH1WWgAaizCk8Osu3K4cK4W28KlEA7liy2VGx8WYGg+GuZZ493X+991F+/WofZ6ydnLzcwZK0mHHd9z9PvAlBENedxopE3z1zYCxrE+ZRagqj1+y+t/lVlC1oa2Jmu61PT/HEspnsP9fKG59f4S82zRv3fVablfKBA6CH7OjJrW8QQriOV009+KuY8GAez1SKz/x2f8247rHb7RzsOw7AltiHXRWaR/varHsJeyMH62eraO8ecHl/AxYrA53B0B/MozP9b5792dUzCdRpOdhSze/Pjf/U05fPfMagvgvNQBD/a41s4RPC00ii4CWevX86A8sq+Ofg/zOuI30PXbyO9p1HMBzawH9b9aQbIvQ8hrBgUuPCADhzxezy/qqbutB/kkXSu1/j/iTvWCzrTNOjQ1mY1Urvk2/y/YtFd7/hC7+s/gyA+f1LmR4R7qrwhBCTJImCl3h04UxsGSYsUe3874oP7nr9qwdq0dh0bE1aR3yo/z58F86Mxo6dqvqpH6x1N6e+SEYWpkT7fOnmsfzX1Y+CTUNzyGXeqT951+t7Bqxce38+Ie89wd/P3eqGCIUQE+WfTzMvpA8MZBUrAfh14513P3QPDPDaEWUF+VceSHd5bJ5s0HiR3uxdvNz9e5f3dby+FYBFKf5RkXE0G2fNIsmsnJiZf+K3d73+zc/r6e4fYpYujZyFC10dnhBiEiRRGIXjrIeb66R7gh/e8xQAjeGXON3aNOZ1f7n/Va49souwZTXcP9d/aieMJiM2Bru+nwYaXN7XrwJ/Qc/TJdhn+HdN/e0ZOQCcDTzBWXPjmNf1WPr550OHAPjag0a02vHt5hET8+qrr7J582Y560FMmiQKo3Cc9aDGQVB38kTGQiK6EkFr5+lPX8Jms912zf4GE7/pL8Ee0cWa+TF+//DdNGMBAN0hrQxYLS7rx2630xZ0DXt4D4sTE13Wjzf4i+X3Ed6RDFo7f3Xk12Ne99jHP+Hoot/AnBq+8oD/relwl23btrFnzx6Sk5PVDkV4KUkUvMyLs74NNg01Ead4/pORD+HBIQubD/8D9kALUV1J7N6YO0Yr/mPdjFlg1UHAEOWXL7qsn9Ot17Dp+8AOj6f7346Hm2m1Gr4Vp+xe+JhP+OnZ26fKfna2nM80+0Bn49trlpDkJ6eaCuGNJFHwMn+x5EEetz0OwG9v/JGK2j8t0sv5sIiOyAY01gDeXPVfCNIFqhWmxwjUBRDZPw2AsoazLuvHUbo5uNdAXIj/Lh512Ln2SRKb74GBIH7y+zpqW7qHf1fb1cLfmn4KwOzWVfzfDf65K0cIb+F1BZcEvLnpL1ny+xvUfZbE8wcOsWFlDO8Fvc3FkBMAfFO/lYdS/OdI6buZqZnBaa5x1HzJZX3sazkHgTDd5j8ndN6JPiiA40/9bzb941vU1NvZvPMjtn01gIM3TrOv6wTW4D6CbsTy8WMF4640KoRQh4woeKEAnY59z+STFhFH3fUefvFRNRfDj4HORlLnbH7x0NfUDtGjLIlQ5r9rLK4r5XyiTzlPYlmYe06p9AYJkaG8/zdPkx4fTt31Hv6+4k3eGfyIzuBWsATws4zvM8PgH0dxC+HNZETBS8WEB1P6/XUUl1ejD9ZywtpFemg8L25+2m/38I9lw/SFvFr1KYOtcdjt4z8rY7xsNhvNQcquii/NWOLUtr1doiGEPQXreWrnx9Q3G4kIjGOWPoXn5zzA1zIXqx2eEGIcJFHwYvNnGPjHb6784qflqsbiybbNXsn/99LjWIZs1LZ0Y0xw7rfYiy1mtJdnQlwrOQ8uc2rbviBtWjgnXtqM3f6kTDMI4YXkq6fwefogHUvSlCJIn19qdXr7J02dBB9dw9qL38Cgl9X7Y5EkQQjvJImC8AurZ8VhDxzg/dpzTm/7aI2SfKzMiHV620IIoTZJFIRf0KQ20LP1VV7T/8Hpbe9tOoddZ2VlRpzT2xZCCLVJojAKTy3hLCbvSeMiALrDWmjr63Vauzf6+6hY9Ht6tv6O5JlOa1YIp5ESzmKqJFEYhaeWcBaTd+/0dHT9IaCzsetCpdPafe3SCdANobUEs2a6ZArC80gJZzFVkigIv6DRaJhuSQHgvcYTTmv33QalrYTBGbItVQjhk+TJJvzGsrC5ABzrveC0Nqu6lbYW6Wc7rU0hhPAkkigIv/FYklIM6VrQlVFP3pyMhgCl2uMj0xc5pT0hhPA0kigIv7F1znKwaRjS93Lgau2U2zt9/RqW0E6ww9bZK5wQoRBCeB5JFITfiNaHkXZ1NcFH7uXC5b4pt7f7UgUAIb0xzIgwTLk9IYTwRJIoCL+yVf8UgdVzOVsz9UTh+mU9QZUrWT14vxMiE0IIzySJgvArq2YpRZGcUcr56PF+gs7dw1+lbZ5yW0II4akkURB+ZfXsadgib1ARup/ajvZJt1Pf2sOFxk60Gg0PL5zuxAiFEMKzSKIg/EpyTCi2rL0MLD/K/zv18aTb+X/HP8GSUc2S+Xqiw4KcGKEQQngWSRSE31mmU7Yy7mk5POk2dnV8yMC9BwhebHJWWEII4ZEkURiFnPXg276W+iAAdcEX6LdYJnx/98AAV0NrAPhG+gNOjU0IZ5OzHsRUSaIwCjnrwbc9v2AN2sEgbMED/PL0xEcVfnH6EPZAC9oBPV+eu9wFEQrhPHLWg5gqSRSE39EHBJE+qJRz/s3lzyZ8/64rBwAwWuYQoA1wamxCCOFpJFEQfumphDUAnLCdwm63T+jeE0NnAHg8fpXT4xJCCE8jiYLwS3+9aAPYNPTrb3C4oWHc9x25epm+iFaww18tfNiFEQohhGeQcVPhl1IjY7i/5iscr4DP6ebelPHd9x/nPgc7RPVMJ8MwzbVBCiGEB5ARBeG3vjp7NZqhAN47Nv7V4OcPRBJW+mW+HfQ1F0YmhBCeQxIF4bceXaasAj90qZmLLW13vf5w9XWO1rQRPBTKDx6Q8x2EEP5BEgXht9KmhTN3TSddT7zOn+3/+V2vf7HsCABb70sjwRDi6vCEEMIjSKIg/Fr26jTsEV1U6g9xsrlpzOsOXKnjLeNP6dvwAd9+JN2NEQohhLokURB+LX/ZI0T2JECAlRcO/nLM6/726G9Aayc6IpDMmQlujFAIIdQliYLwa1qtlvwUpQLn50EHOd1y7bZrLt/ooDJImXb4/szn3BqfEEKoTRKFUchZD/5le+ajRPTEQ6CVFw78222/f37vr7AHDaLvjSY/8xEVIhRi8uSsBzFVkiiMQs568C9arZa/m/FlAA4H7ed7H77OkM2GxTrEfa+9yCch7wPwjegn0Gl1aoYqxITJWQ9iqiRREAL4L8sfI64nBQKt/PKTatb9tw9Y9OudHIr4GIAHrev42do/UzlKIYRwP0kUhEAZVah94l/4tuabxHUYOXG5g4aD0wkwR/O9sG/y6WM/ktEEIYRfkhLOQnwhPDCEVzb+GddX9/M/S09w7uoNdq76F5anSalmIYT/kkRBiFtMi9Tz02+tVjsMIYTwCDL1IIQQQogxeV2iYDKZKC4uZuPGjeO+p7CwkOLiYoqLiyksLLzr9f66jchft4PK6/Yv/vq6/fW5JqbOqxKFqqoqysvLMZvNtLe3j+uewsJCDAYDubm55ObmkpmZSV5e3h3v8dd/UP76AJXX7V/89XX763NNTJ1XJQqZmZnk5uZiNBrHfc+OHTvYsmXL8M9ZWVkUFxe7IrxhU3kQTfUhpuZDUF63e+91xv1q9S2vW537hZgMr0oUJspkMmE2mzEYDLf9rry83GX9+uuDRF63e+91xv1q9S2vW537hZgMjd1ut6sdxESVlpayY8cOKisr73hdeXk5Gzdu5NaXGB0dzSuvvEJ2dvao98XGxhIYGDj8c3Jy8oSqml29enXSVdCmcq/0LX1L39L3zdfePN1gsVhoa2ubdN/Cf/nl9siYmJg7rnGQf0xCCCGEQrVEobS0lF27dt31uu3bt5OZmenUvse7EFIIIYTwd6olCtnZ2WMO/TvLWIsezWbzhBZECiGEEP7KpxczGo1GDAYDJpPptt9lZWWpEJEQQgjhXbwyURhr6sBkMt1WUGn79u0jdjiUlpaSm5vr0viEEEIIX+FVux5MJhNFRUWUl5dTVVVFfn4+GRkZwx/8xcXF7Ny5k5qamhH3OYouAdTU1LBz585R27/5OrPZTH5+vsteiydxJFeOv7eioiI1w1HFxo0bKSsrUzsMtykoKCAjIwNQFve6ehrQExQXFw9vl66pqWH79u2jbp32ZiaTifLyckpKSkZ9P/vrM05MkV3Y7Xa7fefOnfaioqLhn8vKyuy5ubkqRuQe+fn5I37Ozc21Z2VlqRSNOkpKSuz+8k+ho6PDnpmZae/o6LDb7XZ7ZWWlX7z2nTt3Dr9mu135e8jOzlYvIBeorKy0FxUV2Xfu3GnPzMy87ff++owTU+dVIwquFB0dTW1t7YhvGBqN5rYaDL7EbDaTk5NDSUnJ8Ouuqqpi+fLl1NTU+MWCT7PZzO7du8nLy/Pp/9cOBQUFxMbGjvgmWV5e7vNrdpYvX35b3RVfHUUaq86MPz7jhHN45RoFZ1OrgqMnqKioGLHY05EcmM1mlSJyr927d48o8e3rCgsLycrKGh6iBv9Y2Gs0Gtm4cePw+9pkMvlFIuzgz884MXWSKMCouyIADAaDT39gGgwGOjo6RtSpcDw0/OEh6g/fpG/meJ87PjSMRiN5eXl+8UHxyiuv0N7eTnR0NAUFBZSXl/vVWhx/fcYJ55BE4Q7uVsHRF+3YsYOioiKfW+Q1Gn+rp+H4sDAYDGRmZmI0Gtm5cyc5OTkqR+Z6BoOBvLw8srOzKSwspKSkRD4g8c9nnJg4SRTuwN/+ARUUFJCXl+cX20eLi4v9YqX/zWJiYgBYsWLF8J85vlH6+qhCQUEBRqORkpISampqaG9vZ/ny5WqHpTp/e8aJyZFEAangCMoCqJu3mvqyqqqqER+W/sLxXr71w2GsomS+wjHV4phmMhqNVFZWYjQaKS0tVTk695BnnJgKvzwU6lY3V3C89R+NP8xhO75NOpIEs9lMe3u7zz5A2tvbqaqqGn7djvoRhYWFGI1Gnx1pMBgMGI3G297nZrPZpxMnk8k06lSaP0y5OPj7M05MjYwofMFfKzhWVVVRVVVFZmYmJpOJqqoqduzYMTxM7YuysrLIz88f/i8vLw+A/Px8n00SHHbu3DliS2BpaSlZWVlOP3jNk2RlZVFVVXXbmoTKykqf/P891nSCvz7jxNRJHYWbjLeCo68wm82kp6ePuqjLX94WjlNMHQ/NnJwcn/+G5ahQCMqR6r7+Pgflvb5jxw5iY2OH12Xk5ub61KLdu1WuBf97xgnnkERBCCGEEGOSqQchhBBCjEkSBSGEEEKMSRIFIYQQQoxJEgUhhBBCjEkSBSGEEEKMSRIFIYQQQoxJEgUhhBBCjEkSBSGEEEKMSRIFIYQQQoxJEgUhVFRVVUVGRobaYQghxJgkURBCRbt27fLZUzqFEL5BEgUhVFReXu7TJzcKIbyfJApCqKiqqoqNGzeqHYYQQoxJEgUh3Ky8vJy8vLzhBKGkpIS8vDxMJpPKkQkhxO3kmGkhVFJQUEB5eTmVlZVqhyKEEGOSEQUhVFJVVUVWVpbaYQghxB1JoiCESsrLy2V9ghDC48nUgxAqqKqqYvny5cg/PyGEp5MRBSFUcOu2SLPZrF4wQghxB5IoCKGCsrKyEesTiouLVYxGCCHGJomCECowGAzDpZvLy8tlUaMQwmNJoiCECrZv305ZWRmlpaUAUp1RCOGxZDGjEEIIIcYkIwpCCCGEGJMkCkIIIYQYkyQKQgghhBiTJApCCCGEGJMkCkIIIYQYkyQKQgghhBiTJApCCCGEGNP/D/BH7K/+8blgAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "from matplotlib import pyplot as plt\n", "import mpl_config\n", "\n", "t = np.arange(101) * tau\n", "plt.plot(t, z_list, label=\"Variational\")\n", "plt.plot(t, z_exact_list, linestyle=\"--\", label=\"Exact\")\n", "plt.xlim(0, 10)\n", "plt.ylim(-1, 1)\n", "plt.xlabel(\"$t$\")\n", "plt.ylabel(r\"$\\langle \\sigma_z \\rangle$\")\n", "plt.legend(loc=(1, 0.75))" ] }, { "cell_type": "markdown", "id": "9699fb37-1a78-41fd-989a-ad43853849e8", "metadata": {}, "source": [ "## Use ``TimeEvolution`` class\n", "In this subsection we show how to use the ``TimeEvolution`` class to greatly simplify the workflow." ] }, { "cell_type": "code", "execution_count": 8, "id": "ed91ebfd-4a46-441b-ba79-b4a743ec4de6", "metadata": {}, "outputs": [], "source": [ "from tencirchem import TimeEvolution\n", "\n", "te = TimeEvolution(ham_terms_spin, basis_spin, property_op_dict={\"Z\": Op(\"Z\", \"spin\")})\n", "for i in range(100):\n", " te.kernel(tau)" ] }, { "cell_type": "code", "execution_count": 9, "id": "e9ab006c-2379-4104-9c93-c4c6357df5c0", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAgoAAAEhCAYAAAAeQog9AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/NK7nSAAAACXBIWXMAAA9hAAAPYQGoP6dpAABfQ0lEQVR4nO3dd1xc95no/8/MUIY+FAECIWBQtyqqLrJsCVmOi9xAijbV2RiS3Wy5+W3gKrl3b9sbWcT37u7NJmtwstlNcyRwU9zBtmx1C1CviAEhBAJRRnSYYeb3x/EQIYFEmZkz5Xm/Xv4Dcc73+4w8OvPMtzxfjd1utyOEEEIIMQqt2gEIIYQQwnNJoiCEEEKIMUmiIIQQQogxSaIghBBCiDFJoiCEEEKIMQWoHcBEmUwmysvLKSkpoaysbFz3FBYWYjAYADCbzeTn57swQiGEEMJ3eNWIQlVVFeXl5ZjNZtrb28d1jyNJyM3NJTc3l8zMTPLy8lwcqRBCCOEbNN5YR6G0tJQdO3ZQWVl512ujo6Opra0dHlEA0Gg0eOHLFkIIIdzOq0YUJspkMmE2m0ckCQ7l5eXuD0gIIYTwMl63RmEiTCbTqH9uMBgwm81j3qdJDYfEMELbgojotpCcnExycvK4+7169eqErnfWvdK39C19S983X3v16tXhny0WC21tbZPuW/gxuxcqKSmxZ2Zm3vW6srIy+2gv0Wg02ouKisa8T7M0yc4HG+18sNG+8PUf2i1W64Tie/LJJyd0vbPulb6lb+lb+h5LQkLClPoW/sunpx7GcreFkNOahlg/uBHscDrsKP/1yNsTan/btm2Tjm0q9zrjfrX6ltetzv1q9S2vW537hZgMn17MaDKZyMjIuG3hokajoaysjKysrFHvS0xM5Nq1azywZwcH9J8Q1ZWE+bl/d1b4Hmvz5s3s2bNH7TDcTl63f/HX1+14rgkxUT49omA0GjEYDKOuVRgrSQCG5wB/tubrMKTlRkQjxScPuCxOT+Gv31bkdfsXf33dU1kbIfybVyYKY00dmEwmCgsLR/zZ9u3bR+xwKC0tJTc3947tO/5BLYlPZn7fMrAE8ItxbMX0dv76AJXX7V/89XVLoiAmy6sSBZPJREFBAUVFRVRVVVFQUEBxcfHw78vLyykqKhpxT35+PmazmeLiYoqLizl69Oht19xJUWYeYW/kcH5vHBebOp32WoQQQghv4JVrFFzt1jnMLf/4Ke8du8o3H8rgp99arWJkQggxOf66NkNMnVeNKKjlbx6bjx07v75wiLOtshhICCGE/5BEYRzumzONiKxKOje8y84T76sdjhBCCOE2kiiMg0ajYXFsCgAfmX1/UaMQQgjhIInCOH0t/X4AGvUmugf7VY5GCCGEcA9JFEZx9epVNm/ezKuvvjr8Z1+euxztgB57oIV/O3NYxeiEEGL8Xn31VTZv3jzi3AchJkIShVEkJyezZ8+eEfutA3Q6jINzANh15aBaoQkhxIRs27aNPXv2+GQdhdLSUjZu3IhGoyEnJ4eqqqoRvzeZTCxfvpzo6Ojbauz4k6qqKpYvX05GRsak7pdEYQK+FLcSgOPW0ypHIoQQIjs7m5KSEgC2bt1KZmbmiN8bjUa2b99OSUkJ+fn5k+6noKCAnJycCd93c52fqbY1FZmZmezcuXPS90uiMAF/uehhsENveCtV1xrUDkcIIfyewWAgOzt7zEJ6JpPpjiX7x2Pjxo1s3bp1wveNdh7RZNuaqpiYmEnfG+DEOHze3Jh4Mi6v42pNEMeCuslMVDsiIYQQeXl5bNy4EbPZjMFgGPG7W3+ejMkkGsXFxVRUVDilLbVJojBB34p7kh37T7P31HX+/KF5aocjhBCTZrfb6R0cUjuMYaFBOjQazYTvy8rKwmAwUFxcPGKKobi4mC1btgz/7CjnbzQaKSsrIy8vb3i6orS0lIKCAjIzM8nLy6OsrAxQpjReeOEFzGYzNTU142qrvLycsrKyEecP5efnU1VVNWpbAIWFhRiNRkAZBXG8jvLycgoKCgB45ZVXMJlMmEwm2traRkwn3CmeqZJEYYKyFiex483TfHLmGtYhGwE6mb0RQnin3sEhEl/YrXYYw669soWw4Ml9LOXm5lJUVDQiUbh1hGHHjh3k5eVhNBrJzs4mIyODysrK4ekLxzUxMTFs3bqVXbt2Dc/v5+XljejvTm05Rg1u/sAHxmwrJyeHvLy8Efdt3LiRsrIysrKyhu9pb28fjjMjI2PEuow7xTNV8ik3QcuNMYTObKdl3kF+c+ao2uEIIYRA+eZvMpmGdz6YTKbhb+gOJpNpxGnCRqNxxM8GgwGTyURmZuaIBYCjze/fra2x3NpWVVUV5eXlI6YkjEYj7e3tw+3FxMTcttbCaDRiMpmmHM94yIjCBOm0WkKWXaIl4iwlDYd5frEcEiWE8E6hQTquvbLl7he6SWiQbtL3ZmZmYjQaKSoqoqioiPLycnJzc0dc49ghYTabMZlMtLe3097ePuKaW5OLsYynrfGoqKgYtU/HFIIjObj1GoPBMKI/Z8UzGhlRmIRlYUo9hZO91SpHIoRwJcuQlY+unsNqs6odiktoNBrCggM85r/JrE+4WV5eHrt3K1MpZrP5tt9XVVWRk5PD7t27MRgM404KRjPRtm7+9n+z0eJ0RzwTISMKk7ApaTFvXnuTa0EN2Gw2tFrJt4TwNdc7+7nvrR9zKeEwgZVh3Buwkh8teppHUhaoHZoYQ25uLgUFBRQUFNy2BdFsNrNhwwYqKyuHP0QdH9KjTVPcyWTaqqqqGvXPs7Ky2LFjx21/bjKZxr2N0pmvbTTyCTeK0Uo43yx71lKwaRjS91LRfMW9wQkhXO7ghRYe+K/v0fRpGgAWfQ+fBexl07m/Jbv85+oGN0H+VMLZsZCwtLT0thX/JpMJs9k84kPTMTR/a0XHuxlPWzevIXCsexhNZmYmK1asGLGewNGGY+HiaG4eiXDmaxuNJAqjGK2E883iQsMJ650GwGsmOU1SCF/y7/sv8NiOj2js6GNubDz70v+dvwn6DvEdswB4zfoWu6pv3x/vqXy5hPNo8vLybttVAMoHcn5+PgUFBZSWllJaWkpJScnwNsjy8nJ27tyJyWSioKBg+EO+qqpq+GfHNsW7tQVKouAY4SgvL8doNI7aFkBZWRllZWUUFxdTXFzMrl27hos1jXZPYWEhFRUVFBUVDSdFd4pnrH7HS2O32+0TvsvHbd68mT179tzxmkVv/IjTYUdZ07eOQ0/9yE2RCSFcqdp8nXn7/hydycjWoKf51+fvJ1wfCIDNZiO19D/RYDiHvjeaxkf/jWh9mMoRj994nmtCjEZGFCZpTZRSbOnSoEw9COErvvbZy9iC+wmc0cq/5T4wnCQAaLVaytb/CF1fKP0B3Xz/7bI7tCSE75BEYZK+lfEwoW89g/bDhxmy2dQORwgxRUevXeZI4AEAts/4KoG629d6z4uLZ0fC3xL6zlO89sde9p9vcXeYQridJAqTtColiShLLL0DQ5xruKF2OEKIKfrmoX8FnY3oGyn8cOUjY173gzUP8Y3li7Hb4e93HXNjhEKoQxKFSdJptWQaYwGoMLWpHI0QYirerzvL2RDlQ3/H3D+/65bn/5a9hACdhsNdZ3j9wil3hCiEaiRRmIJpszvpf2AvP7u+S+1QhBBT8BdVr4DWzvQbc8hbfN9dr08whJCyqZq+R97nh2d/44YIhVCPJApTkJSkxZpWx1ndGbVDEUJMUlN3F7X68wC8tPD5cd/31ws2AnBBf4pL5usuiU0ITyCJwhRkG5cD0BvaxrWeTpWjEUJMxmcnWwl96zlmns1i27zxH8v7rYWrCeuKB90Qf3dYRhWF75JEYQqWJcwgoC8MtHZKq4+rHY4QYhJeO1KPti+U3JmPTuisAa1Wy5YoZdHje/17sQz55nkQQkiiMEXTrSkAlDXJgiYhvI25Z5Dyk00APLtq5oTvf/HeZ9AMBjMY2slPjkldBeGbJFEYxd3OerjZfH06AGd7a10dlhDCyQo+L6HzofdIXNLK/BmGCd8fHxZBpmUlAD+7/LaTo3MOfzrrQbiGJAqjuNtZDzdbFaPUf2/UNLo6LCGEk73Z/hlD05vImDs06TZ+tPBpAK5prnGjf8BJkTmPv531IJxPjpmeokdSFvAPp7QMDMCgZYigQJ3aIQkhxuFCewstYXUA/H8LH5t0O09lLMT42y0014bwWWoLTy5PcVKE4m5MJhNFRUUUFxcTExMz4jCompqa4cOYbj6sSUycJApTdN/0dBL/1zfp7rNRs7ZrUsOXvqy9u59fn67kjYbPOTV4kaXB8/jo2b+e0KIxIVzhxePvgtZOWFc8G9PmTrodrVZLzpxl/Kz2Am9XNkii4EZGo5GdO3dSXl7OihUryM/PH/F7s9lMTk6OW2MqLi4mNzfXrX26mkw9TJFOp2VBcgwAZxrM6gbjYd44f5r4sq/wnzr/B59FvkdHXA2fRLzDM398Re3QhODtDuVch3Uhq6bc1uOZMwB473gDg1bZ/eBuMTExo/65wWBg48aNbo3FcTy0L5FEwQkWphgAOH3FrGocnsRiHeL5k//IUFgXWmsg6f1zWDi4FPqD+fDTHn67z6R2iMKPNfd20RpWD8D3Fzw65fbunTMNzaoqGh79DS+f3j/l9tytx9o35n/9Q4PjvrZvaGDS1zqT2WzGZFKeMZmZmZjNZpf1dbPi4mIqKirc0pc7ydSDE2hnNNP7pT/yCw7z3/l/aofjEf7y4xJuGK7AkI4jK37Givg07HY7f/faPoqvN/BX/3aElNhQ1i1IVDtU4Yf+7cxB0NoJ7I1gfersKbcXoNMyfbqWi/p+fnP5U/566UNTD9KNwj9+aszfPRa3incy/2H45/i9W+i1jf4hvy56MXtXvjT8c9q+r9NqGf3QvBWRczi65l8mGfGdOZIEgKysLED5EC8qKqKqqoqSkhKys7PJy8sbnirYuXMnBoMBs9lMcXHx8NqGvLw8MjNHFuIqLi4e8XNubi7l5eWUlZVhMpkoLCwEuG0qxFtJouAEs+OjsNnbaO7zvBXPamg29/HHN20EzlrE5sVGVsSnAaDRaPjJs2tpaz7Aa0fq+cbLn3HmJ08TFhykbsDC75y/0o2uM4k54TOctl7mucQH2NFznBOcxGaz3fVgKeFcFRUVFBYW0tbWRmlpKSUlJSN+n5ubS25uLtHR0RiNRgCWL19OZWXliERgx44d5OXlYTQayc7OJiMjg8rKSgwGA8BwHzt37gSgvLyc0tJSsrOzASVJ8ZUEwUESBSf4UtpCaAZrSDf1nR3MjIxWOyRV/Zc/HKO7G+69sZFd60ce16vVanj5hXt5x/4el1NP8r8+1/Pi2s0qRSr8Vf2pcEJqHuGHuWuc1ubfLn2YHZ/+K5aQLt4yneKZWUuc1rarda9/a8zf6TQjd3K1PLR7zGu1mpHJUd3aX4/72qkabTHjaEpKSsjJyaGsrIyYmJjbRgtMJhPl5eXDCxKNRiPl5eVkZ2djNpspKCigo6Nj+PqioiJWrlzp1NfiaSRRcIL0qBgC+sKxhnTzTu1pvrtkrdohqebgxWb+cLAOjQb+8Rsr0Y3yrUofpGNuWgitQQP8e9P7vIgkCsJ9OvssHKttB+DB+QlOazc+LILpvUaaIqv51wvlXpUohAWEqH6tM+Xl5Y1Yl1BVVTWcEGRlZZGVlUVOTs6oCw8dIxGOdQ7t7e20tyvvl4qKCgwGw/Dows3X+zIZG3OSWIvywDnQclHlSNSVe/Jf6H3kXdZl2ck0xo553f9YogzTNUfWcKjhsrvCE4I9Zy9gCe7BGB/OjNgwp7b9aPRqAA4NHHdqu2JijEbjiJGCWxcYLl++HLPZTFVV1W33VlVVkZOTw+7duzEYDMPTFMCEF0XevFbCm0miMIqJlHB2mBWk1Ik/2e2/pZx7BgY5F3QSW3wLDy2Kv+O1G1LmEdszA7R2flRZ6qYIhYB/bniN3ud2E7Ta+eez/M2iDQB0hzVTY251evuT4e8lnKuqqkZ82JtMJoxGI0VFRbfVWDCbzWzYsIGdO3eSm5uL0WgcTg5MJtOYOyjGSiBGS0S8kSQKo5hICWeHzCjljVhv889/jAD/dOwTbPo+tIN6vrfwobte/9VpmwDYZzvIgEX2ngv3ODekjPrdnzD13Q63WhKfTHTzbALPLuRgdbPT258Mfyjh7JgaGE1BQcFwomA2m9m5c+fw9ENWVtaIao4mkwmz2TwisXC07Ug4srOzh3c1OJSXlwPKSIZjFMGRWPgCWaPgJFnJC/iX45EMtodht9v9svLgr698DFGwaGgxwQGBd73+f6x8kp+W/wZrWBc7K8r5+3unvp/d31htVr716a/Y1fU+mdrFvLvx74gOdu5wui+51HGdnvDrAPz5fNesJfqm7Zu8cqyaU9P6wLfXuKnOUcK5qqoKs9k84gPcUcLZMYJQWFhIUVHRiPvb29spLVVGNB3bIPPz8ykoKBheoFhSUkJBQQFbt24d8XNhYeFwQuHY8WA0GsnNzaWgoICMjAyfqdCosdvtdrWD8DSbN29mz549E7rHYrWR8MJuLEM2zvzfp5gZ518Pa3NfHzEfbcEeNMDPk37EdxeuG9d9K/f8Nyr0h0g2z6dhyz+7OErf8k7dKb5y7CfcCLs2/Gf6PgMlC/+eJzIWqhiZ59p+8E1e7P45+u4Y+p79g0v6eOtoPV/96X4WzIjiyI8fd0kfkzGZ55oQIFMPThMYoGVuUiQAp6903OVq3/NS1UfYgwYI6A/lhQX3j/u+/74om8Cz99Cxby7XO/tdGKFv+XHFezxx/u+4EXYNzWAQ93avR9cbxoB9kK+/VEnJoTq1Q/RI77coc8ZzNc6fdnB4YF4Cdp2Fk0NnudDW4rJ+hHAXSRSc6J4UA3bsVF3xjLlJd/pd4ycALGMZAdrxz2g9nr6IzOsb0d4w8NlZ//t7m4whm43/Xfd70NqZZs7gwNKfc/DZ/8zpB4p56OpWLD2B/Kf/OEpHz+DdG/MzF+zK+oRN01w3dxwbEYzm8Y/o31DGy2f3uqwfIdxFEgUn6k25SM+W3/Mry+/VDsWtOvsstJ2ehq5hBn81Z+LrDB5aoGwt3Xv22l2uFACvH6lH8956Is8sp+rRF7l3hrLjZl7cNMq/+2XuSTFwo9fCT947pnKknqWmo5W+MGVh2vPzxz/qNRkLApQRiw9b5f+B8H5emSgUFhZSXFxMcXHxbatPR1NeXk5OTg7FxcWUl5dTUFAwvIDFmYwxsRBkoVnjX9+M3zvWgL12JosvPMtX5078JL57F8RgTalj94DMn96NdcjGj984hcYSzI+MX2FGZNSI32u1Gn7wzDz6Vx3kxZD/yalW/92Fc6uzlzsJPrCWhJqVzIu98/bdqXp8+nIALmmqXdqPEO7gdYlCYWEhBoNhuG53ZmbmiO0to3EU1sjLyyMvL4+MjIzhVarO9OB05VtEj74N69CQ09v3VOUnmwB4akXKpHZ7LM2Iov/BvbRkHOVIY72zw/MpRQeOU32tk5jwYL77yNxRr3kmM43QhG7sgRaeP1g86jX+6GxdD4G1GTzOYy7v688X3A82DYOhnRy66r+1VYRv8LpEYceOHWzZsmX456ysrNtO8hpNZWUldrudmpoal21ZeTA5A4a0EDDE4Sb/qDZot9t5u+szbJFmHpw/uW9pqZExRPQo0w+/PLfPmeH5lO7Bfv6T+R/o2/QuX388noiQ0begarVafpT2VQAqA49Q2SzJF0CFqQ2A5XeoGOosyRFRRPYoJ6P++8UDLu9PCFfyqkTBUQzj5jrbDo6CF2rSBwai7zcA8FmTf5RyPthwmWtLPqb3ibeYnz75LaFLg+YD8EnHCWeF5nP+cv/vsei70YT38rfrl9/x2h+s2IDhRgrobHzrcNEdr/UHNpuNT2yfMTStmSXpBrf0uThQeU9/3C7vaeHdvC5RGI3jDPE72b17N6WlpZSWllJQUOCC6BTTbMq36oo2/xhu/E31IQAiexNIDIucdDtPJys18msDq7HZbE6Jzde8Zf4UgCeDHiE2JPSO12o0Gv678esAnAo8Rntfj8vj82SHmi7Ttng/fVkfsCBl8u/TiXg6WSnYUxdwSd7Twqt5VaIwlpiYmDuW8DQajWRlZZGdnT18vvitNb5v5jjrwfHfRM58MAYpZVIv9F4Z/wvwYnu/GAFYGDj6fPl4fWvBGhjSMhTSw3um884Izaecbm3kRriyFuSHS54Y1z3fW7qWwN4I7AFW/s8J9Ufc1PRGrVI/Ibx3GtF3SbKc5Zvz1xBy9F4CP9pAfav7EzXHGQ+O//z1rAcxdT5RwvlOSQJwW73tLVu2DB9DOto0huOsh8lYHTWf/ZcuYLPe3q6vsdvt1GprAHhi+p2Hwu/GEBxGXN8MWsPr+c2lAzw+a4EzQvQZ/+fkhwBEdCayMillXPfodFqWsIQK9rO78TP+N0+5MkSPtq/tLISAUZvmtj5jQ8NYbbmPCnMbRy61khYf4ba+QTnj4ebzajZvluPcxeR41YjCzQd13OzWQzxudetWSEdy4IojQLfNfICQvVn0nchwetue5khjPYNhN8Cm4Rvz751ye6tClLLDRztqptyWr3m/4wgA9+snVijoe8bHCTq6mp59i7BY/Xf4+6JVmQpcY5jn1n5Xz44D4PBFzzhJUojJ8LpEwWAwjPoBn5WVNeo9ZrOZnJycEfc41jPcKbmYrNnTlfnPtq4BWrt8uyTxf1w8CEBEbzxJYVF3ufruvj97M6FvPEfvJ5l+/aF2q6aeTq6FKh903503+vt8LF9duIzkpmV0tQXx2Tn/qu/h0G+xYA5RinltTl3q1r4XzQrDMvsCu61vurVfIZzJqxIFgO3bt4/Y4VBaWjpiu6PJZBpRhMlgMFBUVDQiKSguLiY7O3vUaYepCgsOYGZcGPbgfiobfLvS4N72k8DU1yc4PGxMJ5ZYuvutnKz3v/MyxnLwbCv6/euIrV3KE8aJTcnotFqeXD4DgLcq/GPdzK3eqT0NAVY0lkAeSXXOe3W8lhijGFh9iMaUo9R3yntaeCevSxTy8/Mxm83DlRmPHj064ujQ8vLy244S3bJlC4WFhcP/tbW1UVJS4rIYzWs+oifnD7zauNdlfXgC++FM9B9n8e0Zm5zSnlarYUWGssf9WO2d1534kw+qmgm4ksq3w7ei1U78n+yXViRimXWRX2p+Sb/V/85/ePuKsuA2pi+JQJ17l2UtnDadoN4o0MDvLnzu1r6FcBavXMyYn58/5u8cFRtvZjAY7niPsyUFxdIInOny3aJLl693c7XZgl6XwtYFy5zWbtSsNvpiyvl5u4lv81+c1q63sg7ZeP+4slr9ieXjW8R4q4fmTcey9BgD+j5+fuozvr9sYtMX3k5/eRYhpx/jyfvTVOk/dSiNak7wwbUTbMc5SbUQ7uR1IwreYEFEKgBXrE0qR+I6+88rx+dmpscSFuy8fDMxIYChGQ1c1MoWSYDikwe5lnGEiORu1nyxMG6i9IGBzLXcA8B/XP7EmeF5hROmTnSt8TyTOrWdOZN1r0H5uz89IOc+CO8kiYILrI5T1kO0B15XORLX+VnjGwwsrWTuArtT230mfQkAPaGtdPT1ObVtb/SLyx9iWXyCqBW1BOgm/8/1KzPWAXBGd5ohm/+cQ9IzYOVsww0AVrihdPNonvsiQWkPuUrfoEWVGISYCkkUXGB9irJgyhLSRWtvt8rRuMaJgGNYFp4iKtm5H+YrE9LQDuhBZ+P1S1L69pxdKQX+TOKaKbXzV0vWgTWAoeA+3q0764zQvMJbl07Sm3mYsHlXSYpxT6GlWz2WvgDNYBD2QCuvXzquSgxCTIUkCi4wLyZB+bADPqq/oHI0ztfe10tvqLIv/Nn0pU5tW6PRED+YBEBZ4ymntu1tasyt9Icpizq/Pm9qdSoi9Xpie5XdD7tMR6Ycm7d4t/EYlnnnYI7za6aMV4BOR/xACtjhwyvnVItDiMmSRGEUjhLOEyndfKvIQWWY80DzJWeF5THeuHQSdDa0A3pWJKQ6vf0FemXq5li3f8/p/vai8oEe3B3NvNjJncx5s+V65ZCi/Z2np9yWtzjerRTvmh3k/PfpRHwz6MuE7d6Gvdr5tVvuxlHKWUo4i8mSRGEUjhLON5c/naiFQ4sJPD+fgTZ1hjtdyfFNP34wGY1G4/T218YptQLq8c99/w5l15Q6Fen2dKe093TKChjScr27F7vduWtLPFW9rQGA1TGzVY3jidnz0FiCOVx93e1/99u2bWPPnj0kJye7tV/hOyRRcJGtkV8iuGI1vVcNaofidFVffNN3fPN3tucyloE1AEtvAB29vl3d8k7OWpS/5/u+WDU/VV+bu5qYN76CtvwhLl3rckqbnqzfYqErRFlQvGnGQlVjyUyPJUCnoamjjyttvarGIsRESaLgIrOnKwfAVPvgA/mKRvmmv26aaw5uWhSXzPzyFwgp+xKnL99wSR+ers8yiDlQWQeSY1zhlDbD9UGsSpsOwL4vtrf6so+uXISAITTWAB5KUffsldDgAGJX1dG78T1+dr5M1ViEmChJFFxkVmIk9sABzg+asA75zna01p5e+gOUnRzPGJe6rJ/MdKVmQJWfVmg8Xd9J6GtbSPzkWaeWHV47T1nr8Ok535+v/rDhDACRffFur8g4mojEPmwJzXzcKrt5hHeRRMFFkmP09GTv4sbGPXx+rV7tcJzm7OUuwkq/zKxPv8qiONfNeS5LiwGgwuT733xHc/jidTR2LWvj50+qbPNYZs3S0fulP/K7Gf+IzebbB28d76gDIFU7uYqWznZfjLKYtHqoVuVIhJgYSRRcRB8YSHC/cqLi/ibfWb1fVduGBg2rE52zwG4sUTN66HnyDUqSfu7SfjzVoYvK3Pqa2dOc2u4jc4zYoswM6Xv5sN63q19GX1hBaOlWvh71lNqhAPBsunJE+I3QZroG/HftjfA+kii4UNyQ8pCvaq9TNxAnchzW5PjG7yrr0tOxR91gILyDOrN/nbpns9l4PeEV+tfsZ36G3qltG/QhxPQqI0G/v3TYqW17ErvdzonLHWj7Q1iXmqZ2OAA8PGM22sFg0A1JMTHhVSRRcKGZQcrCsQs9DSpH4jyvxfwbfWs/IXmma4etZxumEdinLAh9veaYS/vyNB9duchATDPW9FruS5/u9PaXBjvqKfhuQavmG/20dg2g1Wi4J8WgdjgAaLVa4geUolfvXJVEQXgPSRRcaF64Mjd6dahZ5Uico7ajnZ64BoZSL7MiNdHl/SVYlQ/Jg9d9Z+pmPEpMRwEw9EwnSh/i9PafmqGcPVAfaPLZdQq/rj5I38NlRK00ERKk/kJGh0X6WQBUdvlexVbhuyRRcKEVsWkAw9vcvN1rX3yzD+qLYJbBuXPno5kVrCRaZ3vqXN6XJ9nfoazWXxA4yyXtf33+ahjSMRTcx6dXfa9yKMAnLacYSr5K0HTP2jWTlbAYTb8ec4dvJmjCN0mi4EIPzZgDgEXfRXtfj8rRTN0nzcoHWNKQe1aRL49W9r432Bvd0p+nqNPUAbAhfrFL2jfoQ4joUbZJ7qk77pI+1Hauvw6Ae0Jdu+h2ol6Yv5aw17Yy+NlKWm7I6ajCO0iiMApnnPUAMC86nvDqhQRVrKKmxfsLL53qU+rmLw51zTfdWz2cpMyld4e0MjTkH9/A2vq66QtRFm8+4+QDt262mEUE1Bq5fs03HwFNWiW5fCB+nsqRjBQdpmd+kgGAozVtbulTznoQU+WbT4kpcsZZD6AsXlpxfRNBFxbQ2DLopOjU06xTHr4PuunhuyFlNrr2OHRXZnKm2Temb+5mj+kUaO1o+0NYmuC6OhXfTXgW/YEHuX4xymV9qKWx6waDoWYAnkhbpG4wo1g5SykmdqDGPSNlctaDmCpJFFzMUcr5YpN3jyg0dncyGKqUU3bXw1cfEMSaM19Df3Atl5sG3NKn2s40taK5EcW0gSSXHLjlsNyonG56vLYdq4+N1rxddxo0oOsPY35sgtrh3CYkvZmeZ3ZTrHlF7VCEGBdJFFxsZqKeoZhW9refUTuUKTlYV4+2PZbA3gjmxrjv4btghgGAs1fMbutTTZYrCYT98RlesL3g0n7mTI8kXK+jW9/Ggcu+dUrn/mZlR0HsoOclCQD3z0zFHtZLe0gjg0MWtcMR4q4kUXCxG7F19D32Nnsj3lE7lClpvxZE6HtP8kTN99za7z0pBuwaG0ebfKcWxZ2cvKysT1iaGuvSfrRaDfasz+jd/Ca/NH3k0r7c7bL5BgwGYgz0jNLNt3o8fQEaawD2QAtvm86qHY4QdyWJgoutjldWXXcHt3v1nvXTX3yjXzwz2q39WuOb6Pnyb3k34ddu7VcNg1Yrp64oicLiVNf/Pc/RzwTgqI/t6Q+5sIiw3X9GXsxzaocyqiBdIDF9SQD8sf64usEIMQ6SKLjYQzNmA2ALGqDa7L0L8k7WKyu0F6W4N1FYm2wEnY2+kA66B3x7ncJ7dWdpe/bXDGaVk5EQ4fL+Hpp2DwCXuezyvtzFbrdztsGMBg3LUuLVDmdM8wKVrb9Hbvj2eRvCN0ii4GJxoeEE9IUDsLfhosrRTI5lyMpnmT+n9/G3mDbdvaMiKxJmorEEgs5GWb1vffO91YcNZyDQSmioMjXgatkZyiFFfaHtNPV0urw/d2jq6KOjZxCdVsPcJM/d0fFwvJKk1VGnbiBCjIMkCm5gsCjboSpavfN42X0NtdiDB7BF3mBV0gy39q3VaokcUKpAftLo29++jt5QSlWn69wzt744PomA3gjQQOmlKrf06Wq/rtlPz+bXCX7gGPogndrhjGnrrBUA9IW10dBpVjcYIe5CEgU3SNYpq6/Pdnvn6vIPryg7NsL74tAHBrq9/xSNMp97zGxye9/uVGOtByAzyj0FrQCmW5XEr6zJNw6IOth6EXtkJyHRnl23ZGFcEhEtaQSen8/ntS1qhyPEHUmi4AZzwpSH8WVLk8qRTM6RDuWb7gyNOgVbFoYrC0JrLN6ZaI2HzWajI/gaABuT7nFbv0tClTU0J3p84+Ctc711AMwJnqluIOOQ3fkNgitXc6nOs5MaISRRcINHYjMJqlxJWLX7PgCcqXpA+aa7MEyduvn3xSsfZq0BvnEK52gONtZiDxqEIS2PGxe4rd9nk1cTeHIJmnNz3danKzVplGR8Zaz7RmUmy1Gh8WiN9y5yFv5BEoVROOusB4esmfMJOncPrRcNDHnhFsmWAOXhuzZenQ+Tx1IXEnA5DU11BuYe39z58E69MvQf1htHeJDebf0+O3sJ+lPLaLsY7fWHFPVbLPSEKKdFbkyer3I0d7ciIxZ7gIV9HaddunVaznoQUyWJwiicddaDQ0psKMGBWgatNupbe53SprvU3zBjCVVWxD+uUt38DEMcGWceI+jUEs5d9Y3V+bdqbhlC15DCbNsct/YbFRrE7MRIAI7Xdbi1b2f76MpF0A2hsQbwQLJR7XDuat6McHqy/0Dz2j0cuOq6hc5y1oOYKkkU3ECn1TI91Yp1Zh2fNnjXXHDVlWYC6tIIbUsiwxCnWhwLZihb3c41mFWLwZW6TLGE7N1AXtRWt/c9Nz0Ia/IV3qivcHvfzrS3Udk+G94XR4DOc3c8OEQE64noV/5NvWbyjV0nwjdJouAm5gWf0//gXt5t/VztUCbkWiPo9z/E49eeVzWO2dPDsYV1cajlkqpxuMqpejPgnoqMt+pLraH/4Y94e/BDt/ftTA1tfWjbY0glVe1Qxm2WVln3s79dSjkLzyWJgpukBk0H4GKPd51Z4CjdvCjFoGocTQmn6H3mNf4YsEfVOFyhuauHK72t2LFzjwp/z+sSlPn85gDv3JUzrC6F0Hc3893wr6gdybg9EKMsXK0e8u2tv8K7SaLgJvPDlS2SV4e8a+V+5bUr2LGz0M1nPNxq5TRlzrldd13VOFzhjdpj9Gbvxrr5baJCg9ze/9PGxQBYQroweXGZcccJo44TR73Bc8blAHSGtdDZ369yNEKMThIFN1kZp3zQmQO950FsHRriyNJf0bP1d0QmqLvbYP0MZceFJaSL9j7vWhB6NwdblHUrBtQpOZxmiCGwR1nQ+FbtCVVimKrO/gFMrcpCVzVGZSZrbVIG2gE96IYovXRc7XCEGJUkCm6y7ovDoawh3bT2dqsczfgcuFqLPdACWhv3Jak777sgJhHNYBBo7Xxc751nZozldLey4j0t0L3lsW+WOKSsiN/bfE61GKbi9ZpjdG/9LdZHypkW6b7tpVOl1WpJGFT+v7/bcFzdYIQYgyQKbjLbEId2MBiATxu8Y0HeR1eVsxVC+2IICXT/kPjNNBoNkQPKCvEDzd7x9zde9UONACyJVKegFcACvdL36V7vnCvf13wRdDbC9Z6/2+FWj4Y8SNDRVQzWJaodihCjkkTBTbRaLeEDMQAcbvGOh/Hn7coH8nSmqxyJYrpGeZCevFGnbiBO1hGk1Pp/MMG9NRRutnbaPACuar1rsa3DiU7l31RagHqjMpP17Yz1BF1YwNkzNux2u9rhCHGbSScKr7/+Oq+//rozY/F59/asR7/3YULaPeOD924u9iulm+eFeMZ2s7mhyoeAacB3KsydbbuGLbgf7PBImnrVBJ8zZhJ8YC0BHz9E74BVtTgm6/KQ8p5YHOn5hZZutSw9hkCdluYb/dS39qgdjhC3mVCicPz4cb7zne+wcuVKTCYTNTU1rFixgu9+97scP37cRSG6n7NLODtsMCwnoCGVa16yC+3aF3XzV3lI3fwNcUsJPLWYoLoMtUNxmrJ6ZU1AUF8U00LDVYtjbtw0ktvvQWM2DG+J9SYdgcqozNqE2SpHMnEhQQHMmavBklHN784fdXr7UsJZTFXA3S7o7OykqKiIXbt2sXLlSvLy8li2bNnw73/wgx9w7NgxXn75ZSorK9m6dSu5ublERka6NHBXcpRwdjZHqdxL1zy/DHHXQD99oUrd/E0p7juk6E6eTF3KD0/U0xqgxTpkI0Dn/TNn7dc1BJ5bwKzYWFXj0Gg0LEmNpuxkEycvd7BqlnpVOCfqfHszQ3rlnIpNqZ7xXp0oy/yzDIQcorTVzg952Kltb9u2jW3btrF582antiv8x5hP2tdee41HHnmEnJwcoqOjqaio4F//9V9HJAkOy5Yt4+WXX+bo0aNERUWxfv16Nm3aJFMTt0iJD8aacpmq8AMePxd5trGdwDOL0DekszLBM47snRkbhj5Qx6DVxmUfGaK90RhKcOUqcoKfVDsUZqbB4Lwz/L7Vuyo0OkZlAnsjmR7unV9QHopdCED1UI3KkQhxu1EThS1btlBRUUFRUREffPAB3/72t8fd4AsvvEBFRQUvv/wyn3/+OVu3ur92vadKjQ+jf90n3Fh0mIsdnl046HLTIMEnMlnbmI1W6xnf3LVaDakztVgTG/nsim88UM823ABgfrI6NRRuFpLYyeCKo3wedEDtUCak6bqFgLo0UvvUWww6VV+etRKA7vAWmru7VI5GiJFGnXrYvXv3lBtOT0/nxRdfnHI7viQuNJyAvnCsId18erWauTHxaoc0pjNfHL7kOIzJU7QuOkB/+Fneuh7C8yxVO5wpsdlsnOi7gD0ozCOqCT6etpgXz0JPSBvdg/1uPe56KnqvRqHf/xBfedI7px0AViSkEtAfhlXfwx8uVvI3mQ+pHZIQw5zyVbGurs4ZzfgFg0WZ+61odd2xss5wsPUCtpBeFnjAN92bGYOUwkAXe71zG9/NTlxvpGXdHnqy/8DMePU/lO9LSkU7EAxaO2/XnlY7nHFzjMosSDaoG8gUaDQaZliV3UXvNx5XNxghbuGURCE6Onp418OxY8fYvn27M5r1Scm6BADOdl1ROZI7+zSphN7ndjMw7ZraoYywJEp5mDbaPSuuySi7opwYGNxnIEofonI0Sq2P6AGlVkV5o3ecZmiz2TjZXo8du8eNfk3UinClTPnxft+qPCq83113PYxHbW0t+fn5/Of//J9Zv349AJs2beKDDz5wRvO3KSwsxGAwAGA2m8nPz3fJPa4wJ2wGJzhCvcVz90g2dndiCVV2Zjyaqt7e/tHclzCLf6qDruBW7HY7Go1G7ZAm7XBrDWghfihB7VCGGXUzaeMyxzq9o/plZcsVmh/dhWYgiPSEHLXDmZInUzIpvfIGLUENDA3Z0PnArh7hG5zyTtywYQOFhYUsX76curo6KioqqKlxzWIzxwd+bm4uubm5ZGZmkpeX5/R7XGVptPKNuFXjuYsZ369Tvk3q+kOZHe1Z6yjWpyjfumzBA5xva1E5mqk513sZgFn6FJUj+ZPlBqVmRu2QZ494OXx4RdnxEGwJJyxY3TLjU5U9aylR+zYS8senudjk+Vuohf9wSqKwY8cOKioqiIqKoqamhqqqKiorK53R9Kh9bdmyZfjnrKwsiouLnX6Pq6ydrhSE6QvpYHDIokoMd7Ov+QIAMRbPShIAYkPCCOyLAOCjhgsqRzM1V+3KqNJyg+cUkNqYfA8A5uAWbDabytHc3ZFWZeQj0eYd1U7vJDQwmPtDM9EM6Dl00XO/SAj/45REITc3l5wcZdhvw4YNGI1GOjo6nNH0CCaTCbPZPDyFcLPy8nKn3eNKq6fPJOLgw4S89yT11z3zuORTnXWA59bNj7FOA+Dz6967RdJms9GtV44cfyhprsrR/MmjafOI+PAJQku30NDep3Y4d3W+TxmVme1BozJTcf9c5b29/7x3j5YJ3+K0SbCoqD8tJPrBD37gskRhNAaDAbPZ7LR7XClIF8g9liXozNFcuuaZx03XDSk7CpaqeJrhnay1P0DwofsJue6Zicx4HG2uV47wtml4eIbnlB0ODQpioX4WmqFATl52/r9hZ2tEGZVZEeM5ozJTsXhOGAOLj1ES+luPL8om/MeEFjPGxsbS1tY2rmtHq+DoKjExMbS3tzvtHsdZDw6OEqjOMmd6JKevmLnY1MmjS5Od1q4z2O324dMM1yZ6zjfdm30pdg3v12i5HuaUtbiqaGyxEPT5amLjlA9nT7JoZjSn6s2cqu/gieWem4wNDlnoCVGeRxuS56kcjXOsMk7DsvAkFq2dfVdMPDhz8gnQq6++OuK8GjnrQUzWhJ60HR0ddHZ2etw5DhNNEu52j6vOenCITu5ncP5p9pg7+Ws8a1dB841eAqsysUV38KW1nlnAZs50ZY1CtRcv+Gq6NkTQxfk8FOF5H8TTZvYzsPIwv+67wHYWqR3OmD67agLdEFh1rJ3hfadGjiY+NJKo3kRuhDfx2+rDU0oUbv2CI2c9iMma8Fey2tpaysrKaGtrY+XKlWRlZbktcTAaR38YmM3mMX83mXtcLq6dweUVVHY1Ad9RJ4YxnGvoItA0m4yECOJC1DvN8E4yEsIZir9GddQNOvs3EalXv1jRRF1oVJKcuUmet/c/KTEAi/Y8l/rD1A7lji439RN4ehHTogMI0gWqHY7TLAqcy36a+LTjJPAVtcMRYuJrFDIzM/nwww+prKwkPz+f6OhofvjDH7oittsYjUYMBsOo6w6ysrKcdo+r3RuvbEHrDBrfNI47OUo335NiUDWOO0mICqH/4Y8YWH2IT696x37/W+3vPs5QbCvp09UvtHSrJ9KVUQSrvodLZs9dfX+9SUvw8eU8Pvi02qE41ZNJywGo1dXIOgXhESacKLz88st8+OGHfPjhh1y6dInq6mra29vZtGmTK+K7zfbt20fsVigtLSU3N3f4Z5PJRGFh4YTucbcNKcrhNbbgfi60e9bq5g9bqhiKv8asGcFqhzImrVZL2EA0AAeveefOh+Npb9P3pbcJiPW8A4BmRhkI6lFGOt6uO6lyNGM7c8UMwMKZBlXjcLavz7sXbBosYZ0cbLisdjhCTDxRuPU0SKPRyMsvv0xubi6/+MUvnBbYWPLz8zGbzRQXF1NcXMzRo0cpKioa/n15efmIn8dzj7vFhYYP1wL4pMGzyrXuDXufvkfepzfBswvuJKDUeDhh9r4H6YX2FmzB/QBkpXrmgtH4oSQAPms+r3IkYzvSfRabvpd7vLx0860SwyKJ7FWqdf66+pDK0QgxwTUKWVlZfPTRRzzzzDO3/e65557jJz/5idMCu5M7lV92VF+cyD1qiLbE0RLSxeetNXyHB9QOBwDr0NDw3v71SZ61yPJW6cHJ1HCamn7vW8n90RXlwzewN4L4UM9cBzIvOJUGznGqxzMPL2vt66Y68zVYDgmJj6sdjtMtDJzLwYF2Trd6/5kmwvtNaETh5Zdf5gc/+AEnTpwY9ffR0dFOCcofpAQoleTOdtWrHMmf7Ltqwh5ohSEtD6fMUjucO1r8xeFQ1+zNKkcycUdalfUy0V8UjvJE98Yp02NXNZ6ZiL1fdw40oO0P4Z64RLXDcbofpX6FsJJttH+eLusUhOomlCg4phmWLVvGX/zFX4xIGOrq6rh0yTsXlqlhfvhMAC5bG1WO5E/Kryp188P6YtEHeNbe/lutiVd2rHQHt3ndg/RslzJdMjPAc8sOPz5TWdA4YB+gd3BQ5Whut7dJea/GDMZ79cFgY1k3dwaBOi0N7b1cbu1ROxzh5ya8RiErK4v29nauX7/OsmXL0Ol06HQ6cnJy3Lb7wRdkJz5AyLtPEnd0g9qhDDvapiR60/HcDzCHh4cXhA5wsd1zV+aPpv6L5HBBRKrKkYxt5fQUEt/+GqFvZlNzzfM+qI53KqMynlpmfKrCggNYbozFjp2yM54z6ij806RKOBsMBkpKSrDZbFy6dIlLly5x9OhRjyvE5MnWzEhB1x5LfZOFAcuQ2uEAcHFAeSDND0lTN5BxiAsJJ+nkBkI+fJSm6555uNZY2gOVxGZ1nOeWHdZqtSxNVBLGEx5YyvlPZcY99+9wquIWNdP7dCk7W3+ldijCz035rIf09HTS0z3zTIDJcpRwvrn8qbPFR+mJCg3EZrdjavaMLXLNGmXhlCd/gN1spW0lupZELjf3qx3KuHX3Wwg4eC9BFSvZNNOzF4wuSVXWHJ3ysETh5jLj6zy0zLgzPGCcgT28h9qgiwxarZNu59VXX2Xz5s1SwllM2qiJwksvvcTx48en3Pjx48d56aWXptyOuzlKODvzfIdbaTQaohY20b/qEK/XH3VZP+PVPziE7uAagj9fw1OpS9UOZ1xmT1dGsLyplHN1UxcBTclMv5pJRnSc2uHckSGlh771H/KrQNdve56I062NyvZSm4YvpXlmmXFn+NaCe9FYArHp+/jd2cpJt7Nt2zb27NlDcrJnnSsjvMeoicKGDRv48Y9/zOzZs3nppZfo7Bz/g7izs5Of/OQnrFixgh07dqhW/dAbWJMasM65wL6O02qHwsWmTjQt8cRfXcw907xjFXlswhCWjGrKeverHcq4XWi8AcDcJM+fprsnKZahpEZaIy5jGZr8N1pnq2saIPjQ/SSaVhMb6tllpqciNEBPyoAyuvebOu95jwvfM2qisGzZMnbv3k11dTVRUVGsX7+eTZs28frrr4/Z0GuvvcYjjzzChg0biI6OpqKigl27drF06VJXxe71ZumVDP9Sf4PKkfypdPOCGQavWUWuibnBwL0HOGk4oHYo4/ZecyWWVBNJKZ6/U2Nj6hywBmAPsPJxQ7Xa4Qy73DhIYM1s1ls3qh2Kyz0So5Rzrhjw3AqZwvfddY3CCy+8QEVFBS+//DKff/45s2bN4rvf/S7Hjx/n+PHjfOc732H27NlUVFRQVFTE0aNH+fa3v+2O2L3e0ug0AJo16pdxfuvaISzGamakecbCyvFYO302AP0hZvot3rGg8WPbpwys/Yyu+NvPHvE0+sBAInuVCpjvX1F/1MvhtKN0swefR+Isf3nPegC6Iq5xukWKLwl1jHsxY3p6Oi+++CKXLl0iOzubH//4x+zYsYOcnByqq6vZsWOHzy1qdLUHEpUPup7gdmw2m6qx7LXvY+C+A/QnqD+6MV6rEmeCVQc6GweuemYFwVtd1ylJ4YpY7zgWOVWbAsCRjgsqR/Ine/sOMzStmVnJvjvt4LB0WgphPXGgtfOz03vVDkf4qQkfMw3KGoYNGzxn/7+3Wpc8C85rsQdaONR0mfuT1Uu02gKVCodrE7xnFXmATkdofzS94a3sv3aJDWlz1A7pjnoHBxkIMQPwULJnx+qwLGIWp+xHuDToGWdq9FsHuTD7fZhnIy7xCbXDcYt1AfdRdr6ay8E6WK92NMIfTXl7pJi80KAgQvpiAChvOKtaHHU32rGGdAPwaKp3rSKPsytlkI931KkbyDjsbbgEWjsaa4AyGuIFNiQp74e24GseUQHz4yvVoLOhsQSyNtk7RmWm6n/es43gitUcr7Rjsao78ij8kyQKKptuV3YYHG9Tr/raH02nAOWQotTIGNXimIy0IOWUw+o+z98jvu+asiAwrC8WrdY7/uk9aVyIpk+Ppj2aS+3taofDx1+UGY/ojyNAp1M5GvdYlhZDXEQwXf1WDld7VxVS4Ru842nlw7YFPkfYrj8jrn6JajE4jhKOH/L80s23WhipfDNvtHn+Qi/HqEciCeoGMgHRISEsO/QCIeWPYmpQv7DV0Q6lzPgMjf/UBNBqNWxYHI81oYl/Of2J2uEIPySJgsrunZGKxhLE2YYbqsVwskdZgT9Pn6ZaDJP1dNJq9OWPEHZondqh3NWl/isAzNJ71/kEy9JiATjpARUaLwzWAbAk3D+mHRwiFzXRv/ED3gh4DeuQ9+xMEr5BEgWV3fPFFq8LTTewDqkz/3gVZdh+dYx3LLC72ZoZKQRcS8J8PYCOHs875fBmMRdWof9kA0/G3q92KBOy+ItSzpWX1R/2bg1sAuDhRO9aSzNVP1r2GBprAJawGxSfOKx2OMLPSKIwCnec9eAwMzYM+/JjmB98T1mo5WaD1iECyjeg/2gjW9JXur3/qYoICWR6dAjg2aWc7XY79XUQcDWF9SnelZBFJ/XT89RrvDbzp6rGcbq1kSF9H9g0PJWxSNVY3C02JJw5g/cA8DPTexO6V856EFMlicIo3HHWg4NWq0E7o5mh6U18dPWMy/u71YXGToZ6gom7kcbiBO8o3Xyr6DltDCypYk/D5Ovhu9rV9l66+60E6DRkJESoHc6EPJiWhj28C2tIN5fM6o0q1F+1EPLBl5h5bgPxYd71d+gMeemPAnAu+ATmvvGvF5GzHsRUSaLgAZI1ysr9ig73V+s7Va/MOy+cGe01pZtv1ZdUi2XRST7tPKZ2KGN6v/4MgwtPEDe3g8AA7/pnlxoVTVCvAYC3ak+oFse5+m501xN4OOg+1WJQ0/cWriNgIBR78AA7jn6gdjjCj3jXE8tHzQ9NBeDSgPu3SP62uYyBJVXEG3vd3rezzAlVFgdeHmxSOZKxfXj9GINLj9E/65zaoUxKolVJZj9tVq/eh2Mx5ZI079rC6yyBugBWaVYA8LtrH6scjfAnkih4gDVxyglxzdpmt/f9ub0Cy6KTBMart+tiqpbFpAFwXav+mRljOd+j7Hhw1H3wNgtDlffoyd4a1WL4QPceFmM1c1JCVItBbfn3KNUom2iiydytcjTCX0ii4AEemams4B4IMWPu63NbvzabjY5gpf7AhqT5buvX2W4+HGrA6pmHQ10dUkY7FkWmqRvIJG1IWAhAo+6KKv1f7mynOeMoA/cdYG6K/61PcNictoTlJ79KyNtP8fph7zmXRXg3SRQ8wOK4JLSDwaC182H9ebf1W9l8BVvwANg0PJZ2j9v6dbZ7p6fBkBZ0Qxxq9IwzCW51I6gVgHvjM1SOZHKem70U7GAJ6eZ8h/tHvvbUKscsB/ZGMis6zu39ewqNRkPe0jVo0PBy2QUp6SzcQhIFD6DVaokamAb9wRy/1ui2ft+5rJRuDumLIUrvvcO5gboAQvqVvf77mty/xfRuam+0K9v6gA0p3nPo1s1So6KJas4goHoOx1Sop/Bps7K2I97qfdVDne3P1hqZFqmn1tzBi5/tVzsc4QckUfAQX+n5c8JKv0zANfdtYTrcpnyoJtu9c978ZtNsyuFQxzzwcKjyL0aJAvrCSIqIUjmayXum+6voj9xHQ737d8ec7P6iemiwHGUfFhzA008G0/N0Kf/Q8TN6Bz270JjwfpIoeIglyfFo0HC2wey2Ps/11wFwT6j3l8N9RrOZ0DeeY3rjUrVDuc2R68oCQIN1msqRTE1murLboNLU5va+GzTK2oj74ryrWJWr5D/wAFotDEZ08L29u9UOR/g4SRQ8xPwZyjfNc1fdt/ugxa4MIT8wzTuHw292b2IG2p4Iapp61A7lNnHN8wj949M8Nvik2qFMyXJjLHbtEAfbzrv1yOmO/h76QpStkU+mqXd4midJDItic9AmAH7b/RZd/QMqRyR8mSQKo3BnCWeHjOmh9D1UzrkH/oOGLtcnC119FnR7Hif0zWfZkuF9pZtvNTsxEoCLHljG+VJjD9obBtbGz1M7lCm5JzWKnuw/cGVtKUeu1bqt33drz4LWjrY/hBWJ3nWgliv9cu3z6Ab1WMJvkPvx78a8Tko4i6mSRGEU7izh7JAYGY4m1ow9vJs3alxfYfB4XTsaNMwMSmSmwXvnzR2MieEMLD7GlcXvU3/DrHY4IzjOoJibFKlyJFMTpQ8mrF+Zfni91n1VMC3N0YS+toV767K9tnqoK8Tow9kS+jgAf7C/xk8rRl/YuG3bNn75291SwllMmiQKHiTRonxb+uia6898qPhinnlFRqzL+3IHQ2gwtjkXsabV8UnDRbXDGWbu6+NcxvsMLjiFMTFU7XCmLE2jVBHd3+a+CpMn6jrQ9oWybtpCt/XpLX754DdJ7E2DQAt/07ST107f/v/lZ8f2MesP290fnPAZkih4kIUhyqLCUz2ur35X1P17+h/YS1yG71R3i7Io++sPt6hXPfBWH1+5iCXjEpaFp5geFaZ2OFO2KkpZz3LB4r5zSSpMSg2K5Ub/LN18JyEBwZx59P8S3ZtEQM0svvf/TvLHyis0dfQyYLGy8a1/4nvN/0BnmvsPnBO+QxIFD7IuXqnQeFXr+oprNfqzWNPqSE/Su7wvd0nSJgBwpsv9Z2aMZV+zMroRORCHVuv9/9wen6EsJuzQX8MyZHV5f009nRw27mZg8TGWZ0S7vD9vFKMP59Kmn/Fwz2N09w/xZ/+8j5m/KyDk3ecoD3kXtHZm9/rXsdzCubz/yeVDns1YCsBAmJn6zg6X9XOypRFrSDfYIWdWpsv6cTfH4VB1g+4rWnU3x8zKor9kjW8UCnosYx5YArAHWCm74voqoruqK7AmXcU+y8TMWP8t3Xw3MSFhlHz/Ib661ogxPhx7lBl7SB8M6fhuyNe58NRLaocovFiA2gGIP5kTE09gbwSW0C7euHSCv8l8yCX9lNZUARDSG8v0cO9fyOiwIsbI6x1wXYXDtcZiGmiAYJgXOlPtUJwiJDAIQ+90zFFX+GP9cR5Lc+26gQ8aT0AAJFtnykLGu4gMCeRfX1gDQG3XvexvusSq+HTmGhJVjkx4O0kUPMwMSzqXrzdxccjssj72Xj8Deki1p7qsDzWsnzEPOpTDoXoGBwgLClY7JFp0StKyOs47z3gYzb32NXxUFY8l1fULYU/0V0M4LA/z7q2l7pYeMY30CO8u8CU8h0w9eJi/0D9P6AeP02Ny3cE35waVxX4rI32ryt2KhBQ0ViX3PdBQp24wQNdAPwMhZuCLJMZHPJ+ykaCzi7h4zrXf8G02G81BSkXGx2ZIoSUh1CKJgodZlqas7D5e1+6S9q1DQ7TplTn8x1N86+Gr0+pYfeJ5wv7wVXpb1T/k6tMrJtDY0VgCyYz3nUJB982NB+DE5Q46+1x3rPfeq9XYggbAqiN7tm+9V4XwJpIoeJilXyQKl1o7aOlyfjniisYmNF0RaAYDeTLd91ZCL4lLQWPXcr7RfaWwx9LfGkrYH77KqjNf84kdDw7JMaEkz7AzkFrD789WuKyf12orAYjsTfTq002F8Ha+8/TyEbERwbDpY3q2/p7fVR9xevu19RZC33uSh4//FaFBQU5vX23zvqh+eN6NZ2aM5fzVTjRDAayI9b0TD3WZZxl44DP+o7HcZX2cu94Cg4HM0Xn/oWVCeDNJFEahxlkPN4vRh4HWzsfNzi+SUlmjVGRcme6bC50iE/rpv28fb0eqf6Ke44Cvecm+s7PE4cEYZTTqjNV1VTBtp+YRVrKNv07IcVkf/kDOehBTJYnCKNQ46+FmwxUae51fYfDzWmUV/soM1y2WVJMxMQKrsYa22GoGrK6bPx+P96N207/6IDGJri9M5G5fmb0agK6QFlp6nX8QV9+glROXO9DYtayd7TvrO9Swbds29uzZI2c9iEmTRMEDPZxwDwBNOudWaOwc6GPfip/R+6U/kpHqe9MOAPdPT4chHQQMceBqnWpxdA30055wEevsi8xO9L1CQffPSCOgJxK0dn5z0flTZFW1bViH7MRH6UmN8/7S10J4M0kUPNAzxqUADIZ2UtPR6rR2Sy4eg0ArhPewOCHBae16kkBdAKF9SqnfvY0XVIvj4yvVoFV2PCyPT1EtDlfRaDSkWpW1F+80Vjm9/R21f6Bn82sYVtRJoSUhVCaJggfKiI4jqFeZ1/7NBed9Wyut/xyApMFUn1qFf6tElEp0lR3uO7joVp81KXP3Ef2+ccbDaO6PUqoynuh3fkJ2rPcC9sguUqb5zlkkQngr33yC+YBZttkAvHPNedvPjvYpiyPXRi52WpueaFaw8g2+uu+KajEcMytJSpKPnPEwmi3GVQC0hzbSPdjvtHZ7rf006esAyJ65ymntCiEmRxIFD7V52hoCao101zqnTG57Xw9tYcqah6/NesApbXqqpdFpADRxTbUYagaUv2tfOeNhNJvS5xBzOIvQN7K52OC8mh//fu4I9gAr2v4Qtsxb6rR2hRCT43VnPRQWFmIwGAAwm83k5+ff8fry8nKKiorYuHEjRqORsrIyVq5cSXZ2thuinby/XbSRf3mlmwagrWtAqa8wBb86exh0Q+j6w9iUOtc5QXqoddPnUHhBQ79lCJvNjlbr/jnuli8OploV67s1AAJ0Oh4OW8UHfY0cvHidTKNzktrd9QchEFIHZhEU4HWPKCF8jleNKDiShNzcXHJzc8nMzCQvL++O95jNZqqqqsjLyyMvL4+MjAyPTxIAEgwhw/vv959vmXJ7b1xV1iekWzN8ds7cIStlLoaSr6N/9wmutDm/uuXd9FssDNqVrZkbZsx3e//udP88pZzz3jPOG72pHFSmyLJifOcIdCG8mVd9YuzYsYMtW7YM/5yVlUVxcfFd76usrMRut1NTU0Nubq4rQ3SqB+fHMxTVwe9qDk65LXNtBAGX0/hSzGonRObZggICmJNgAFCllLOpuYfQt54j8a1vsCLB93Y83GzDwgQG559hT/yvuNI99fNJqs3NdIcpiXHu/HVTbk8IMXVekyiYTCbMZvPwtMPNystdV0ZWTcEZTfQ9+RZ/DHpjSu2Yewa5UhWLft9DFCx+3EnRebZ5ScpozPmrzi8GdDdnG8wAzI/33R0PDotmxqCZbcKS0MQ/nZz6v8N9F68ReGEeUS1GViRJoSUhPIHXTACaTKNvdTMYDJjN5jveu3v3bmJilMOWjh49ys6dO+94vaOEs8O2bdtUqdL4rfn3sePov9Af3s6plkYWxSdNqp3951uw2e3MSowgOSbUyVF6pqGUenqj3+FXvRf5G/67W/s+VW8GYOHMaLf2qwaNRsMy7WIO8wlvtRzi/7Dl7jfdwcmzAwQfXcMLG33rCHQ1vPrqqyPK0EsJZzFZXpMojCUmJob29rGHPI1GI5mZmRiNyqKy9vZ2cnJyKCkpGfMeRwlntc2KiSOsO56eiBZ+eX4//xQ/uYfw76oPY4s0s27BLCdH6LmmxQRgC7lObaf73+L/bvs1fRvaCEny/UQB4Otp6zjc+gmmoIv0WQcICZj8wtuPTytrHdYv9N1tpe5y6xecm7/8CDERqiUKpaWl7Nq1667Xbd++nczMsRc13SlJAG67d8uWLeTl5Y05jeFpFgTM4SgtlLUeh0l+W3sz6DX6N7cREja5EQlvtD55Pv90CbpCWrEMWQnUue+t3hhSy1BMD+kJvle6eTTPL1rJ994NxRbSyytn9/PXizdMqp0DV+q4ZK8hOCCBtfPjnRylEGKyVEsUsrOzJ7T7wDEicCuz2Tzm70BJSG7ux5EcmEymOyYgnuKJxBUcvbGfS9rqSd1/sqWR/vA2sMO3Ftzn5Og818aZc+GCFnughc+umNiQ5p6h7Or261hDlJ0WT6YvckufatMHBjJrcD4XQyr57eXPJp0o/OP5d+nb9B6RHXMI1/+Zk6MUQkyW16y0MhqNGAyGUdcqZGVljXqP2WwmJydnxD2O9Qx3Si48ybcW3A82DYOhNzjUWDfh+//pVBkAYT3xzI72n29p+oAgwvuUEzI/bHD+cd1jebvuFACBvZEkR/je8dJjeTZBSUJP2E9ht9sn1ca+rmMArIm4x2lxCSGmzmsSBVCmIW7e4VBaWjpiu6PJZKKwsHD4Z4PBQFFR0YikoLi4mOzsbK+YdgCYERFFVI8yX/uPJz+c8P1vdnwKwDq9/5XCnaFRjtU90jG50ZjJ2NeinHuQMORfc+x/vexhNP3BDLUaONk08bofx1uv0BJ2GYDvzFnv7PCEEFPgVYsZ8/PzKSwsHK6dUFNTQ1FR0fDvHVUYb67WuGXLlhHJQ1tb2x0XMnqi5yOf4eXyi9QHJGHfZB/3aXqfXamhI7IB7PA/Mp92bZAeaFFYOuc5RvVAvdv6PNNTCxEwJzjVbX16gukR4Ww681fsP9fKgeQOliRN7HTSH1WWgAaizCk8Osu3K4cK4W28KlEA7liy2VGx8WYGg+GuZZ493X+991F+/WofZ6ydnLzcwZK0mHHd9z9PvAlBENedxopE3z1zYCxrE+ZRagqj1+y+t/lVlC1oa2Jmu61PT/HEspnsP9fKG59f4S82zRv3fVablfKBA6CH7OjJrW8QQriOV009+KuY8GAez1SKz/x2f8247rHb7RzsOw7AltiHXRWaR/varHsJeyMH62eraO8ecHl/AxYrA53B0B/MozP9b5792dUzCdRpOdhSze/Pjf/U05fPfMagvgvNQBD/a41s4RPC00ii4CWevX86A8sq+Ofg/zOuI30PXbyO9p1HMBzawH9b9aQbIvQ8hrBgUuPCADhzxezy/qqbutB/kkXSu1/j/iTvWCzrTNOjQ1mY1Urvk2/y/YtFd7/hC7+s/gyA+f1LmR4R7qrwhBCTJImCl3h04UxsGSYsUe3874oP7nr9qwdq0dh0bE1aR3yo/z58F86Mxo6dqvqpH6x1N6e+SEYWpkT7fOnmsfzX1Y+CTUNzyGXeqT951+t7Bqxce38+Ie89wd/P3eqGCIUQE+WfTzMvpA8MZBUrAfh14513P3QPDPDaEWUF+VceSHd5bJ5s0HiR3uxdvNz9e5f3dby+FYBFKf5RkXE0G2fNIsmsnJiZf+K3d73+zc/r6e4fYpYujZyFC10dnhBiEiRRGIXjrIeb66R7gh/e8xQAjeGXON3aNOZ1f7n/Va49souwZTXcP9d/aieMJiM2Bru+nwYaXN7XrwJ/Qc/TJdhn+HdN/e0ZOQCcDTzBWXPjmNf1WPr550OHAPjag0a02vHt5hET8+qrr7J582Y560FMmiQKo3Cc9aDGQVB38kTGQiK6EkFr5+lPX8Jms912zf4GE7/pL8Ee0cWa+TF+//DdNGMBAN0hrQxYLS7rx2630xZ0DXt4D4sTE13Wjzf4i+X3Ed6RDFo7f3Xk12Ne99jHP+Hoot/AnBq+8oD/relwl23btrFnzx6Sk5PVDkV4KUkUvMyLs74NNg01Ead4/pORD+HBIQubD/8D9kALUV1J7N6YO0Yr/mPdjFlg1UHAEOWXL7qsn9Ot17Dp+8AOj6f7346Hm2m1Gr4Vp+xe+JhP+OnZ26fKfna2nM80+0Bn49trlpDkJ6eaCuGNJFHwMn+x5EEetz0OwG9v/JGK2j8t0sv5sIiOyAY01gDeXPVfCNIFqhWmxwjUBRDZPw2AsoazLuvHUbo5uNdAXIj/Lh512Ln2SRKb74GBIH7y+zpqW7qHf1fb1cLfmn4KwOzWVfzfDf65K0cIb+F1BZcEvLnpL1ny+xvUfZbE8wcOsWFlDO8Fvc3FkBMAfFO/lYdS/OdI6buZqZnBaa5x1HzJZX3sazkHgTDd5j8ndN6JPiiA40/9bzb941vU1NvZvPMjtn01gIM3TrOv6wTW4D6CbsTy8WMF4640KoRQh4woeKEAnY59z+STFhFH3fUefvFRNRfDj4HORlLnbH7x0NfUDtGjLIlQ5r9rLK4r5XyiTzlPYlmYe06p9AYJkaG8/zdPkx4fTt31Hv6+4k3eGfyIzuBWsATws4zvM8PgH0dxC+HNZETBS8WEB1P6/XUUl1ejD9ZywtpFemg8L25+2m/38I9lw/SFvFr1KYOtcdjt4z8rY7xsNhvNQcquii/NWOLUtr1doiGEPQXreWrnx9Q3G4kIjGOWPoXn5zzA1zIXqx2eEGIcJFHwYvNnGPjHb6784qflqsbiybbNXsn/99LjWIZs1LZ0Y0xw7rfYiy1mtJdnQlwrOQ8uc2rbviBtWjgnXtqM3f6kTDMI4YXkq6fwefogHUvSlCJIn19qdXr7J02dBB9dw9qL38Cgl9X7Y5EkQQjvJImC8AurZ8VhDxzg/dpzTm/7aI2SfKzMiHV620IIoTZJFIRf0KQ20LP1VV7T/8Hpbe9tOoddZ2VlRpzT2xZCCLVJojAKTy3hLCbvSeMiALrDWmjr63Vauzf6+6hY9Ht6tv6O5JlOa1YIp5ESzmKqJFEYhaeWcBaTd+/0dHT9IaCzsetCpdPafe3SCdANobUEs2a6ZArC80gJZzFVkigIv6DRaJhuSQHgvcYTTmv33QalrYTBGbItVQjhk+TJJvzGsrC5ABzrveC0Nqu6lbYW6Wc7rU0hhPAkkigIv/FYklIM6VrQlVFP3pyMhgCl2uMj0xc5pT0hhPA0kigIv7F1znKwaRjS93Lgau2U2zt9/RqW0E6ww9bZK5wQoRBCeB5JFITfiNaHkXZ1NcFH7uXC5b4pt7f7UgUAIb0xzIgwTLk9IYTwRJIoCL+yVf8UgdVzOVsz9UTh+mU9QZUrWT14vxMiE0IIzySJgvArq2YpRZGcUcr56PF+gs7dw1+lbZ5yW0II4akkURB+ZfXsadgib1ARup/ajvZJt1Pf2sOFxk60Gg0PL5zuxAiFEMKzSKIg/EpyTCi2rL0MLD/K/zv18aTb+X/HP8GSUc2S+Xqiw4KcGKEQQngWSRSE31mmU7Yy7mk5POk2dnV8yMC9BwhebHJWWEII4ZEkURiFnPXg276W+iAAdcEX6LdYJnx/98AAV0NrAPhG+gNOjU0IZ5OzHsRUSaIwCjnrwbc9v2AN2sEgbMED/PL0xEcVfnH6EPZAC9oBPV+eu9wFEQrhPHLWg5gqSRSE39EHBJE+qJRz/s3lzyZ8/64rBwAwWuYQoA1wamxCCOFpJFEQfumphDUAnLCdwm63T+jeE0NnAHg8fpXT4xJCCE8jiYLwS3+9aAPYNPTrb3C4oWHc9x25epm+iFaww18tfNiFEQohhGeQcVPhl1IjY7i/5iscr4DP6ebelPHd9x/nPgc7RPVMJ8MwzbVBCiGEB5ARBeG3vjp7NZqhAN47Nv7V4OcPRBJW+mW+HfQ1F0YmhBCeQxIF4bceXaasAj90qZmLLW13vf5w9XWO1rQRPBTKDx6Q8x2EEP5BEgXht9KmhTN3TSddT7zOn+3/+V2vf7HsCABb70sjwRDi6vCEEMIjSKIg/Fr26jTsEV1U6g9xsrlpzOsOXKnjLeNP6dvwAd9+JN2NEQohhLokURB+LX/ZI0T2JECAlRcO/nLM6/726G9Aayc6IpDMmQlujFAIIdQliYLwa1qtlvwUpQLn50EHOd1y7bZrLt/ooDJImXb4/szn3BqfEEKoTRKFUchZD/5le+ajRPTEQ6CVFw78222/f37vr7AHDaLvjSY/8xEVIhRi8uSsBzFVkiiMQs568C9arZa/m/FlAA4H7ed7H77OkM2GxTrEfa+9yCch7wPwjegn0Gl1aoYqxITJWQ9iqiRREAL4L8sfI64nBQKt/PKTatb9tw9Y9OudHIr4GIAHrev42do/UzlKIYRwP0kUhEAZVah94l/4tuabxHUYOXG5g4aD0wkwR/O9sG/y6WM/ktEEIYRfkhLOQnwhPDCEVzb+GddX9/M/S09w7uoNdq76F5anSalmIYT/kkRBiFtMi9Tz02+tVjsMIYTwCDL1IIQQQogxeV2iYDKZKC4uZuPGjeO+p7CwkOLiYoqLiyksLLzr9f66jchft4PK6/Yv/vq6/fW5JqbOqxKFqqoqysvLMZvNtLe3j+uewsJCDAYDubm55ObmkpmZSV5e3h3v8dd/UP76AJXX7V/89XX763NNTJ1XJQqZmZnk5uZiNBrHfc+OHTvYsmXL8M9ZWVkUFxe7IrxhU3kQTfUhpuZDUF63e+91xv1q9S2vW537hZgMr0oUJspkMmE2mzEYDLf9rry83GX9+uuDRF63e+91xv1q9S2vW537hZgMjd1ut6sdxESVlpayY8cOKisr73hdeXk5Gzdu5NaXGB0dzSuvvEJ2dvao98XGxhIYGDj8c3Jy8oSqml29enXSVdCmcq/0LX1L39L3zdfePN1gsVhoa2ubdN/Cf/nl9siYmJg7rnGQf0xCCCGEQrVEobS0lF27dt31uu3bt5OZmenUvse7EFIIIYTwd6olCtnZ2WMO/TvLWIsezWbzhBZECiGEEP7KpxczGo1GDAYDJpPptt9lZWWpEJEQQgjhXbwyURhr6sBkMt1WUGn79u0jdjiUlpaSm5vr0viEEEIIX+FVux5MJhNFRUWUl5dTVVVFfn4+GRkZwx/8xcXF7Ny5k5qamhH3OYouAdTU1LBz585R27/5OrPZTH5+vsteiydxJFeOv7eioiI1w1HFxo0bKSsrUzsMtykoKCAjIwNQFve6ehrQExQXFw9vl66pqWH79u2jbp32ZiaTifLyckpKSkZ9P/vrM05MkV3Y7Xa7fefOnfaioqLhn8vKyuy5ubkqRuQe+fn5I37Ozc21Z2VlqRSNOkpKSuz+8k+ho6PDnpmZae/o6LDb7XZ7ZWWlX7z2nTt3Dr9mu135e8jOzlYvIBeorKy0FxUV2Xfu3GnPzMy87ff++owTU+dVIwquFB0dTW1t7YhvGBqN5rYaDL7EbDaTk5NDSUnJ8Ouuqqpi+fLl1NTU+MWCT7PZzO7du8nLy/Pp/9cOBQUFxMbGjvgmWV5e7vNrdpYvX35b3RVfHUUaq86MPz7jhHN45RoFZ1OrgqMnqKioGLHY05EcmM1mlSJyr927d48o8e3rCgsLycrKGh6iBv9Y2Gs0Gtm4cePw+9pkMvlFIuzgz884MXWSKMCouyIADAaDT39gGgwGOjo6RtSpcDw0/OEh6g/fpG/meJ87PjSMRiN5eXl+8UHxyiuv0N7eTnR0NAUFBZSXl/vVWhx/fcYJ55BE4Q7uVsHRF+3YsYOioiKfW+Q1Gn+rp+H4sDAYDGRmZmI0Gtm5cyc5OTkqR+Z6BoOBvLw8srOzKSwspKSkRD4g8c9nnJg4SRTuwN/+ARUUFJCXl+cX20eLi4v9YqX/zWJiYgBYsWLF8J85vlH6+qhCQUEBRqORkpISampqaG9vZ/ny5WqHpTp/e8aJyZFEAangCMoCqJu3mvqyqqqqER+W/sLxXr71w2GsomS+wjHV4phmMhqNVFZWYjQaKS0tVTk695BnnJgKvzwU6lY3V3C89R+NP8xhO75NOpIEs9lMe3u7zz5A2tvbqaqqGn7djvoRhYWFGI1Gnx1pMBgMGI3G297nZrPZpxMnk8k06lSaP0y5OPj7M05MjYwofMFfKzhWVVVRVVVFZmYmJpOJqqoqduzYMTxM7YuysrLIz88f/i8vLw+A/Px8n00SHHbu3DliS2BpaSlZWVlOP3jNk2RlZVFVVXXbmoTKykqf/P891nSCvz7jxNRJHYWbjLeCo68wm82kp6ePuqjLX94WjlNMHQ/NnJwcn/+G5ahQCMqR6r7+Pgflvb5jxw5iY2OH12Xk5ub61KLdu1WuBf97xgnnkERBCCGEEGOSqQchhBBCjEkSBSGEEEKMSRIFIYQQQoxJEgUhhBBCjEkSBSGEEEKMSRIFIYQQQoxJEgUhhBBCjEkSBSGEEEKMSRIFIYQQQoxJEgUhVFRVVUVGRobaYQghxJgkURBCRbt27fLZUzqFEL5BEgUhVFReXu7TJzcKIbyfJApCqKiqqoqNGzeqHYYQQoxJEgUh3Ky8vJy8vLzhBKGkpIS8vDxMJpPKkQkhxO3kmGkhVFJQUEB5eTmVlZVqhyKEEGOSEQUhVFJVVUVWVpbaYQghxB1JoiCESsrLy2V9ghDC48nUgxAqqKqqYvny5cg/PyGEp5MRBSFUcOu2SLPZrF4wQghxB5IoCKGCsrKyEesTiouLVYxGCCHGJomCECowGAzDpZvLy8tlUaMQwmNJoiCECrZv305ZWRmlpaUAUp1RCOGxZDGjEEIIIcYkIwpCCCGEGJMkCkIIIYQYkyQKQgghhBiTJApCCCGEGJMkCkIIIYQYkyQKQgghhBiTJApCCCGEGNP/D/BH7K/+8blgAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.plot(te.t_list, te.properties[\"Z\"][:, 0].real, label=\"Variational\")\n", "plt.plot(te.t_list, te.properties[\"Z\"][:, 1].real, linestyle=\"--\", label=\"Exact\")\n", "plt.xlim(0, 10)\n", "plt.ylim(-1, 1)\n", "plt.xlabel(\"$t$\")\n", "plt.ylabel(r\"$\\langle \\sigma_z \\rangle$\")\n", "plt.legend(loc=(1, 0.75))\n", "plt.savefig(\"sbm.pdf\")" ] } ], "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.7.13" } }, "nbformat": 4, "nbformat_minor": 5 }