{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# VQE with gradients, active spaces, and gate fusion\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This tutorial will explore the Variational Quantum Eigensolver, a hybrid quantum classical algorithm for determining the ground state energy of molecules. The first part of this tutorial will walk through the key aspects of the VQE algorithm and how to implement it with CUDA-Q.  The following sections explore advanced topics: parallel gradients, active spaces, and gate fusion"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The Basics of VQE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The VQE algorithm is hybrid quantum-classical algorithm, meaning some subroutines run on a quantum computer or quantum simulator and others run on a traditional (super)computer.  \n",
    "\n",
    "The goal is to take a parameterized quantum circuit and a qubit form of the molecular Hamiltonian, measure an expectation value that corresponds to the ground state energy of the molecule, and then repeat the process to variationally minimize the energy with respect to the parameters in the quantum circuit. The optimization is performed on a classical device while the expectation values are determined on a quantum device.  See the figure below.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![VQE.png](./images/VQE.png)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The next few cells will elaborate on each part of the VQE procedure and show you how to build a VQE simulation to compute the ground state energy of the water molecule."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Installing/Loading Relevant Packages"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Install the relevant packages.\n",
    "!pip install pyscf==2.6.2 openfermionpyscf==0.5 matplotlib==3.8.4 openfermion==1.6.1 -q"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "import openfermion\n",
    "import openfermionpyscf\n",
    "from openfermion.transforms import jordan_wigner, get_fermion_operator\n",
    "\n",
    "import os\n",
    "import timeit\n",
    "\n",
    "import cudaq\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy.optimize import minimize\n",
    "import numpy as np"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Implementing VQE in CUDA-Q"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Like most quantum chemistry programs, the first step is to specify a molecular geometry, basis set, charge, and multiplicity. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "geometry = [('O', (0.1173, 0.0, 0.0)), ('H', (-0.4691, 0.7570, 0.0)),\n",
    "            ('H', (-0.4691, -0.7570, 0.0))]\n",
    "basis = 'sto3g'\n",
    "multiplicity = 1\n",
    "charge = 0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The VQE procedure requires some classical preprocessing.  The code below uses the PySCF package and OpenFermion to compute the Hartree Fock reference state and compute the integrals required for the Hamiltonian."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "molecule = openfermionpyscf.run_pyscf(\n",
    "    openfermion.MolecularData(geometry, basis, multiplicity, charge))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, the Hamiltonian is built using `get_molecular_hamiltonian`. The Hamiltonian must then be converted to a qubit Hamiltonian consisting of qubit operators. The standard Jordan-Wigner transformation is used to perform this mapping. \n",
    "\n",
    "Finally, the Jordan-Wigner qubit Hamiltonian is converted into a CUDA-Q spin operator which can be used to evaluate an expectation value given a quantum circuit."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/tmp/ipykernel_225414/3022299596.py:7: ComplexWarning: Casting complex values to real discards the imaginary part\n",
      "  spin_ham = cudaq.SpinOperator(qubit_hamiltonian)\n"
     ]
    }
   ],
   "source": [
    "molecular_hamiltonian = molecule.get_molecular_hamiltonian()\n",
    "\n",
    "fermion_hamiltonian = get_fermion_operator(molecular_hamiltonian)\n",
    "\n",
    "qubit_hamiltonian = jordan_wigner(fermion_hamiltonian)\n",
    "\n",
    "spin_ham = cudaq.SpinOperator(qubit_hamiltonian)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, the quantum circuit needs to be defined, which models the wavefunction. This is done in CUDA-Q by specifying a CUDA-Q kernel. The kernel takes as an input the number of qubits, the number of electrons, and the parameters of the circuit ansatz (form of the wavefunction) yet to be defined. \n",
    "\n",
    "The number of qubits corresponds to the potential positions of electrons and is therefore twice the number of spatial orbitals constructed with the chosen basis set, as each can be occupied by two electrons. \n",
    "\n",
    "The Hartree-Fock reference is constructed by applying $X$ bitflip operations to each of the first $N$ qubits where $N$ is the number of electrons. Next, a parameterized ansatz is chosen. Theoretically, any set of operations can work as an ansatz, however, it is good practice to use an ansatz that captures the underlying physics of the problem. The most common choice for chemistry is the Unitary Coupled Cluster Ansatz with Single and Double excitations (UCCSD).  This UCCSD hate operations are automatically added to the kernel with the `cudaq.kernels.uccsd(qubits, thetas, electron_num, qubit_num)` function. \n",
    "\n",
    "The STO-3G water molecuule UCCSD ansatz requires optimization of 140 parameters."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "140\n"
     ]
    }
   ],
   "source": [
    "electron_count = 10\n",
    "qubit_count = 2 * 7\n",
    "\n",
    "\n",
    "@cudaq.kernel\n",
    "def kernel(qubit_num: int, electron_num: int, thetas: list[float]):\n",
    "    qubits = cudaq.qvector(qubit_num)\n",
    "\n",
    "    for i in range(electron_num):\n",
    "        x(qubits[i])\n",
    "\n",
    "    cudaq.kernels.uccsd(qubits, thetas, electron_num, qubit_num)\n",
    "\n",
    "\n",
    "parameter_count = cudaq.kernels.uccsd_num_parameters(electron_count,\n",
    "                                                     qubit_count)\n",
    "\n",
    "print(parameter_count)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The classical optimizer requires a custom cost function which is defined below. The `cudaq.observe()` function computes an expectation given the Hamiltonian and the kernel defined above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "def cost(theta):\n",
    "\n",
    "    exp_val = cudaq.observe(kernel, spin_ham, qubit_count, electron_count,\n",
    "                            theta).expectation()\n",
    "\n",
    "    return exp_val\n",
    "\n",
    "\n",
    "exp_vals = []\n",
    "\n",
    "\n",
    "def callback(xk):\n",
    "    exp_vals.append(cost(xk))\n",
    "\n",
    "\n",
    "# Initial variational parameters.\n",
    "np.random.seed(42)\n",
    "x0 = np.random.normal(0, 1, parameter_count)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The final step is to run the optimization using the scipy minimize function and a selected optimizer, in this case COBYLA. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "UCCSD-VQE energy =   -70.21329161891076\n",
      "Total number of qubits =  14\n",
      "Total number of parameters =  140\n",
      "Total number of terms in the spin hamiltonian =  1086\n",
      "Total elapsed time (s) =  21.523745890706778\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAj4AAAHHCAYAAAC/R1LgAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjMsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvZiW1igAAAAlwSFlzAAAPYQAAD2EBqD+naQAAWmBJREFUeJzt3Xd8VfX9x/HXvRk3g0wIGWRA2HsTUZApiFahoDgAwYVQUOssCBa1VagoWtv+qLYUq1ate6AyBEQEZG9I2AQyICRkkJ3c8/sj5NYwQ7jJvTf3/Xw87uPBvefccz/3lHrffKfJMAwDERERETdgdnQBIiIiInVFwUdERETchoKPiIiIuA0FHxEREXEbCj4iIiLiNhR8RERExG0o+IiIiIjbUPARERERt6HgIyIiIm5DwUdERETchoKPiLiEW2+9FT8/P/Ly8i56zpgxY/D29iYzMxOA/Px8/vCHP9CpUyf8/PwICgqib9++vPvuu1xotx6TyXTRx6RJk2rtu4lI3fF0dAEiItUxZswYvv76az7//HPuueee844XFBTw5ZdfcuONN9KwYUNOnDjBoEGD2Lt3L3feeSdTp06lqKiITz/9lHvuuYfFixfz7rvvYjZX/fffDTfccMHrt2rVqta+m4jUHQUfEXEJt956KwEBAbz//vsXDCZffvkl+fn5jBkzBoDx48ezd+9ePv/8c2699VbbeY888ghPPfUUr7zyCl26dOGpp56qcp1WrVoxduzY2v0yIuIw6uoSEZfg6+vLyJEjWb58OSdPnjzv+Pvvv09AQAC33norP//8M0uWLGHChAlVQk+l2bNn07JlS+bMmUNhYWFdlC8iTkLBR0RcxpgxYygrK+Ojjz6q8npWVhZLlizh17/+Nb6+vnz99dcAF2wZAvD09OTuu+8mKyuLtWvXVjlWVFTEqVOnznuUlJTUzpcSkTql4CMiLmPgwIFERkby/vvvV3n9448/prS01NbNtWfPHgA6d+580WtVHqs8t9KCBQsICws77/HZZ5/Z86uIiINojI+IuAwPDw/uvPNOXnvtNY4cOULTpk2Bim6u8PBwBg0aBGCb+RUQEHDRa1UeO3eW2PDhw5k6dep553fs2NEeX0FEHEzBR0RcypgxY3jttdd4//33eeaZZzh+/DirV6/mkUcewcPDA6gaaoKDgy94ncrA07hx4yqvR0dHM3jw4Nr7AiLiUOrqEhGX0r17d9q0acMHH3wAwAcffIBhGLZuLoB27doBsGPHjotep/JYfHx8LVYrIs5GwUdEXM6YMWPYtWsXO3bs4P3336dly5b07NnTdvyWW24B4J133rng+8vLy23dY9dff32d1CwizkHBR0RcTmXrzu9//3u2bdtWpbUH4JprrmHIkCEsXLiQRYsWnff+GTNmsG/fPp5++mk8PdXjL+JOTMaF1m0XEXFy1113nW0q+v79+2nRokWV4ydOnGDgwIEkJiZy991307dvX4qLi/nss8/44YcfGDt2LO+88w4mk8n2HpPJdNGVm8PDw7nhhhtq90uJSK1T8BERl/R///d/TJkyhV69erF+/foLnnPmzBnmzZvHRx99xMGDBykqKgLg2Wef5YUXXjjv/F+GoHP169ePH374wS61i4jjKPiIiNtISUnh2muvpaysjHXr1hEbG+vokkSkjmmMj4i4jSZNmrB48WKKiooYNmwYp0+fdnRJIlLH1OIjIiIibkMtPiIiIuI2FHxERETEbSj4iIiIiNtQ8BERERG3oSVLz2G1WklNTSUgIOCSa3qIiIiI8zAMg7y8PKKiojCbL96uo+BzjtTUVGJiYhxdhoiIiNTAsWPHiI6OvuhxBZ9zBAQEABU3LjAw0MHViIiISHXk5uYSExNj+x2/GAWfc1R2bwUGBir4iIiIuJjLDVPR4GYRERFxGwo+IiIi4jYUfERERMRtKPiIiIiI21DwEREREbeh4CMiIiJuQ8FHRERE3IaCj4iIiLgNBR8RERFxGwo+IiIi4jYUfERERMRtuEzwadq0KSaTqcpjzpw5Vc5ZsmQJ11xzDQEBAYSFhTFq1CiOHDnimIJFRETE6bhM8AF44YUXSEtLsz0efvhh27HDhw8zfPhwBg4cyLZt21iyZAmnTp1i5MiRDqz4fwpKykhKz6O4rNzRpYiIiLgtl9qdPSAggIiIiAse27x5M+Xl5fzxj3/EbK7Ic08++STDhw+ntLQULy+vuiz1PNe/vJJTZ0pY9HAfOjQJcmgtIiIi7sqlWnzmzJlDw4YN6dq1K3PnzqWsrMx2rHv37pjNZhYuXEh5eTk5OTm8++67DB482OGhByAm1A+Ao5kFDq5ERETEfblMi88jjzxCt27dCA0NZe3atUyfPp20tDTmzZsHQLNmzVi6dCmjR4/moYceory8nN69e/Ptt99e8rrFxcUUFxfbnufm5tZK/XGhfmxNzuZoVn6tXF9EREQuz6EtPtOmTTtvwPK5j8TERAAef/xx+vfvT6dOnZg0aRKvvvoqf/nLX2yhJT09nQcffJDx48ezceNGVq1ahbe3N7fddhuGYVy0htmzZxMUFGR7xMTE1Mp3jW3oD0CyWnxEREQcxmRcKhXUsoyMDDIzMy95Tnx8PN7e3ue9vnv3bjp06EBiYiKtW7fm2WefZfHixWzcuNF2zvHjx4mJiWHdunVcc801F7z+hVp8YmJiyMnJITAwsIbf7HyfbTnO4x9tp3d8Qz6YeOFaREREpGZyc3MJCgq67O+3Q7u6wsLCCAsLq9F7t23bhtlspnHjxgAUFBTYBjVX8vDwAMBqtV70OhaLBYvFUqMarkRcw4oxPslZavERERFxFJcY3Lxu3Tpef/11tm/fzqFDh/jPf/7DY489xtixYwkJCQHg5ptvZuPGjbzwwgvs37+fLVu2cO+99xIXF0fXrl0d/A0gNrSiqys1p1BT2kVERBzEJYKPxWLhww8/pF+/frRv354XX3yRxx57jLfeest2zsCBA3n//ff54osv6Nq1KzfeeCMWi4XFixfj6+vrwOorNGrgjZ+3B4YBx7IKHV2OiIiIW3LoGB9nVN0+wpq48fUfSUzP418TejCwTbhdry0iIuLOqvv77RItPvVF07Mzu7SWj4iIiGMo+NShygHOCj4iIiKOoeBTh2JtwUeLGIqIiDiCgk8dijs7s+uoprSLiIg4hIJPHars6jqeVUi5VWPKRURE6pqCTx2KCvbFy8NESbmV9NwiR5cjIiLidhR86pCH2UR0yNlxPqc0zkdERKSuKfjUsdjQs8FH43xERETqnIJPHWuqKe0iIiIOo+BTx2LPLmKYnKWuLhERkbqm4FPH4s52dR05pRYfERGRuqbgU8cqp7QnZxWgbdJERETqloJPHYs52+JzpriMrPwSB1cjIiLiXhR86piPlweRQT6AZnaJiIjUNQUfB6ic0p6smV0iIiJ1SsHHASrH+RzRZqUiIiJ1SsHHAeIqp7SrxUdERKROKfg4QGWLj8b4iIiI1C0FHweIC61o8dHqzSIiInVLwccBYs+2+Jw6U8yZ4jIHVyMiIuI+FHwcIMjXi2A/L0DjfEREROqSgo+DVG5doT27RERE6o6Cj4NUzuzSOB8REZG6o+DjIJrZJSIiUvcUfBykcvXmo1rEUEREpM4o+DiIurpERETqnoKPgzQ929WVml1ISZnVwdWIiIi4BwUfBwkLsODr5YHVgJTsQkeXIyIi4hYUfBzEZDLZxvlos1IREZG6oeDjQJUrOGsRQxERkbqh4ONAleN8NMBZRESkbij4OFDs2ZldWr1ZRESkbij4OFCcbYyPWnxERETqgoKPA1Wu3pycVYDVaji4GhERkfpPwceBooJ98TCbKCmzciKvyNHliIiI1HsKPg7k5WEmOsQX0ABnERGRuqDg42CVa/loSruIiEjtU/BxsMpxPlrEUEREpPYp+DhYXOjZzUqz1OIjIiJS2xR8HCxOqzeLiIjUGQUfB4s7u4jhkcx8DENT2kVERGqTgo+DVQ5uzisqI7ug1MHViIiI1G8KPg7m6+1B4wALoHE+IiIitU3BxwnE2TYr1cwuERGR2uRSweebb74hISEBX19fQkJCGDFiRJXjycnJ3Hzzzfj5+dG4cWOeeuopysrKHFPsFagc56MBziIiIrXL09EFVNenn37Kgw8+yEsvvcTAgQMpKytj165dtuPl5eXcfPPNREREsHbtWtLS0rjnnnvw8vLipZdecmDll1e5Wam6ukRERGqXSwSfsrIyHn30UebOncv9999ve71du3a2Py9dupQ9e/bw/fffEx4eTpcuXfjDH/7A7373O5577jm8vb0dUXq1xKqrS0REpE64RFfXli1bSElJwWw207VrVyIjIxk2bFiVFp9169bRsWNHwsPDba8NHTqU3Nxcdu/efdFrFxcXk5ubW+VR1yq7urRfl4iISO1yieBz6NAhAJ577jlmzpzJokWLCAkJoX///mRlZQGQnp5eJfQAtufp6ekXvfbs2bMJCgqyPWJiYmrpW1xc07MtPifziiksKa/zzxcREXEXDg0+06ZNw2QyXfKRmJiI1WoFYMaMGYwaNYru3buzcOFCTCYTH3/88VXVMH36dHJycmyPY8eO2eOrXZFgP28CfSp6HZM1zkdERKTWOHSMzxNPPMGECRMueU58fDxpaWlA1TE9FouF+Ph4kpOTAYiIiGDDhg1V3nvixAnbsYuxWCxYLJaalG9XcQ392ZmSw5HMfFpHBDi6HBERkXrJocEnLCyMsLCwy57XvXt3LBYLSUlJ9OnTB4DS0lKOHDlCXFwcAL179+bFF1/k5MmTNG7cGIBly5YRGBhYJTA5q9iGfuxMydGUdhERkVrkErO6AgMDmTRpErNmzSImJoa4uDjmzp0LwO233w7AkCFDaNeuHePGjePll18mPT2dmTNnMmXKFKdo0bmcynE+R7M0s0tERKS2uETwAZg7dy6enp6MGzeOwsJCEhISWLFiBSEhIQB4eHiwaNEiJk+eTO/evfH392f8+PG88MILDq68euJCNbNLRESktpkMbQleRW5uLkFBQeTk5BAYGFhnn/vzoUzufOtn4hr6seqpAXX2uSIiIvVBdX+/XWI6uzuo3K/r+OlCSsutDq5GRESkflLwcRLhAT54e5optxqkZhc6uhwREZF6ScHHSZjNJpqdXcF5b1rdrx4tIiLiDhR8nEjv5g0BWLUvw8GViIiI1E8KPk6kf+uKNY1WJmagMeciIiL2p+DjRK6Jb4iPl5n03CIS0/McXY6IiEi9o+DjRHy8PLiueSMAViaddHA1IiIi9Y+Cj5Pp36Ziu40fEjXOR0RExN4UfJxM/1YV43w2J58mp6DUwdWIiIjULwo+TiYm1I+WjRtQbjVYfUCtPiIiIvak4OOEBpzt7lqRqHE+IiIi9qTg44Qqp7WvSsrAatW0dhEREXtR8HFCPeJCaWDxJDO/hJ0pOY4uR0REpN5Q8HFC3p5m+rTQtHYRERF7U/BxUgPPjvNZmaQBziIiIvai4OOk+p0d57PjeDanzhQ7uBoREZH6QcHHSYUH+tA+KhDDgB+1aamIiIhdKPg4sQGt1d0lIiJiTwo+TmxAm4rurh/3ZVBWbnVwNSIiIq5PwceJdYkJIdjPi5zCUrYdy3Z0OSIiIi5PwceJeZhNXN+yotVH09pFRESunoKPk6vs7lqp3dpFRESumoKPk+vXqjEmE+xJyyU9p8jR5YiIiLg0BR8nF+rvTZeYYAB+UHeXiIjIVVHwcQH/m9au4CMiInI1FHxcQGXw+Wn/KUrKNK1dRESkphR8XED7qEAaNbCQX1LOpiNZji5HRETEZSn4uACz2UT/1prWLiIicrUUfFyEtq8QERG5ego+LqJPy0Z4mE0cOHmGY1kFji5HRETEJSn4uIggXy+6x4UAmtYuIiJSUwo+LmRgG3V3iYiIXA0FHxdSOc5n7cFTFJWWO7gaERER16Pg40JahTcgKsiHolIr6w9rWruIiMiVUvBxISaTietaNAIqWn1ERETkyij4uJhrWzQE4OeDmQ6uRERExPUo+LiY3vEVLT47U3LILSp1cDUiIiKuRcHHxUQE+dCskT9WAzYc0jgfERGRK6Hg44J6N6/o7lp3SN1dIiIiV0LBxwX1jq8IPms1zkdEROSKKPi4oGvOBp+9abmczi9xcDUiIiKuQ8HHBYUFWGgV3gCAn9XdJSIiUm0KPi6qsrtL43xERESqz6WCzzfffENCQgK+vr6EhIQwYsQI27Ht27dz1113ERMTg6+vL23btuXPf/6z44qtZb2bVy5kqOAjIiJSXZ6OLqC6Pv30Ux588EFeeuklBg4cSFlZGbt27bId37x5M40bN+a9994jJiaGtWvXMnHiRDw8PJg6daoDK68d18SHYjLBgZNnOJlXROMAH0eXJCIi4vRMhmEYji7icsrKymjatCnPP/88999/f7XfN2XKFPbu3cuKFSuq/Z7c3FyCgoLIyckhMDCwJuXWmZvfWM3u1FzeuKsrt3aOcnQ5IiIiDlPd32+X6OrasmULKSkpmM1munbtSmRkJMOGDavS4nMhOTk5hIaGXvKc4uJicnNzqzxchW2cj/btEhERqRaXCD6HDh0C4LnnnmPmzJksWrSIkJAQ+vfvT1bWhVcvXrt2Lf/973+ZOHHiJa89e/ZsgoKCbI+YmBi7119bbAsZapyPiIhItTg0+EybNg2TyXTJR2JiIlarFYAZM2YwatQounfvzsKFCzGZTHz88cfnXXfXrl0MHz6cWbNmMWTIkEvWMH36dHJycmyPY8eO1cp3rQ29moXiYTZxJLOA1OxCR5cjIiLi9Bw6uPmJJ55gwoQJlzwnPj6etLQ0ANq1a2d73WKxEB8fT3JycpXz9+zZw6BBg5g4cSIzZ868bA0WiwWLxXLlxTuBAB8vOjQJYvuxbNYdzGRU92hHlyQiIuLUHBp8wsLCCAsLu+x53bt3x2KxkJSURJ8+fQAoLS3lyJEjxMXF2c7bvXs3AwcOZPz48bz44ou1VrczubZ5Q7Yfy2atgo+IiMhlucQYn8DAQCZNmsSsWbNYunQpSUlJTJ48GYDbb78dqOjeGjBgAEOGDOHxxx8nPT2d9PR0MjIyHFl6rasc4PzzoUxcYIKeiIiIQ7nMOj5z587F09OTcePGUVhYSEJCAitWrCAkJASATz75hIyMDN577z3ee+892/vi4uI4cuSIg6qufT2ahuDlYSIlu5BjWYXENvRzdEkiIiJOyyXW8alLrrSOT6Xb/76WjUdOM2dkR+7sFevockREROpcvVrHRy5N+3aJiIhUj4JPPfDLfbvUgCciInJxCj71QNfYYLw9zWTkFXMwI9/R5YiIiDgtBZ96wMfLg+6xFYO8tX2FiIjIxSn41BPXNtc4HxERkctR8KknKvft+vlQFlarxvmIiIhciIJPPdEpOhg/bw+y8ktIOpHn6HJERESckoJPPeHtaaZH01BAu7WLiIhcjIJPPVK5ns9aBR8REZELUvCpRyoHOK8/nEm5xvmIiIicR8GnHmkfFUiAxZO8ojJ2p+Y4uhwRERGno+BTj3h6mEmI1zgfERGRi1HwqWeu0b5dIiIiF6XgU89Uruez4XAWpeVWB1cjIiLiXBR86pm2EYEE+3lRUFLO5qOnHV2OiIiIU1HwqWfMZhM3tA0HYOGaw1d1reTMAj7dfFwzxEREpN5Q8KmHJl4fD8DSPSc4mHGmRtcoK7dy79sbeOLj7fxz9SF7liciIuIwCj71UMvwAAa3Dccw4K1VNQstX21P5WBGPgB/XXGAU2eK7VmiiIiIQyj41FOT+1e0+ny+NYUTuUVX9N6ycit/Xr4fqNgKI6+4jHnL9tm9RhERkbqm4FNPdY8LpWfTEErKrfzrpysb6/PZ1hSOZhbQ0N+bN8d1B+DDDckkpufWRqkiIiJ1RsGnHpvUrzkA/1mfTE5habXeU1Jm5Y2zrT2T+jVnQOvG3NQxAqsBf1i0B8PQQGcREXFdCj712IDWjWkV3oAzxWX8Z/3Rar3nk83HOX66kEYNLIy9Jg6AaTe2xdvDzJoDmSzfe7I2SxYREalVCj71mNls4qHrK1p9/vXTEYpKyy95fnFZOX9dUdHa85v+zfH19gAgtqEf9/VpBsBL3+6lpEwLI4qIiGtS8Knnbu0SRVSQD6fOFPPpluOXPPejjcdIzSkiPNDC3QmxVY5NGdCcRg28OXQqn3d/rl7rkYiIiLNR8KnnvDzMPNC3YobXP348dNHFCItKy/nrygMATB3QAh8vjyrHA3y8eHJIawD+/P0+TueX1GLVIiIitUPBxw3c2SuGYD8vjmQWsHhX+gXPeX99Midyi4kK8mF0z5gLnnN7jxjaRgaSW1TG699reruIiLgeBR834OftyT29mwLw91UHz5uZVVhSzv/9cBCAqQNbYvH0OPcSAHiYTTz7q7YAvLc+mf0n8mqvaBERkVqg4OMmJlzbFB8vMztTclh7MLPKsfd+PsqpM8XEhPpye4/oS17n2uaNGNIunHKrwR+/2VubJYuIiNidgo+bCPX35s6eFQOW559t3QHILy7j76sqnj88sCVeHpf/K/HMTW3x8jCxal8GK5M0vV1ERFxHjYJPfn6+veuQOnB/n2Z4mE38dOAUO4/nAPDvdUfIzC8hrqEfI7s2qdZ1mjbyZ8K1TQF48Zu9lJZreruIiLiGGgWf8PBw7rvvPn766Sd71yO1KCbUj1s6RQLw9x8PkldUyls/Vmxi+uiglnhWo7Wn0tSBLQn19+bAyTO8vz65VuoVERGxtxoFn/fee4+srCwGDhxIq1atmDNnDqmpqfauTWrBQ2e3sfhuZxovfL2H7IJS4sP8Gd6leq09lYJ8vXj8hlYAvPb9PnIKqrclhoiIiCPVKPiMGDGCL774gpSUFCZNmsT7779PXFwcv/rVr/jss88oKyuzd51iJ20jA+nfOgyrAR9vrljQ8LeDW+FhNl3xte7sGUOr8AZkF5Ty301q9REREed3VYObw8LCePzxx9mxYwfz5s3j+++/57bbbiMqKorf//73FBQU2KtOsaPKzUsBWoU34OaOkTW6jqeHmfFnx/p8tiXFHqWJiIjUqqsKPidOnODll1+mXbt2TJs2jdtuu43ly5fz6quv8tlnnzFixAg7lSn2lNAslF5NQwF4YkjrGrX2VPpVxyi8PcwkpuexJzXXXiWKiIjUCs+avOmzzz5j4cKFLFmyhHbt2vGb3/yGsWPHEhwcbDvn2muvpW3btvaqU+zIZDLxj/E9SM4soGN00FVdK8jPi0FtG/PdrnQ+33qcdlHt7FSliIiI/dWoxefee+8lKiqKNWvWsG3bNqZOnVol9ABERUUxY8YMe9QotSDI1+uqQ0+lkd0qFj38YlsqZZraLiIiTqxGLT5paWn4+fld8hxfX19mzZpVo6LEtfRrFUaInxcZecX8dOAU/Vs3dnRJIiIiF1SjFp+ysjJyc3PPe+Tl5VFSol273Y23p5lbO0cB8PlWDXIWERHnVaPgExwcTEhIyHmP4OBgfH19iYuLY9asWVit6vZwF5XdXUt2p5NXpDV9RETEOdWoq+vtt99mxowZTJgwgV69egGwYcMG/v3vfzNz5kwyMjJ45ZVXsFgsPPPMM3YtWJxTp+gg4sP8OZSRz3e70hndI8bRJYmIiJynRsHn3//+N6+++iqjR4+2vXbLLbfQsWNH3nzzTZYvX05sbCwvvviigo+bMJlMjOoWzdwlSXy+JUXBR0REnFKNurrWrl1L165dz3u9a9eurFu3DoA+ffqQnKzVfN3JiLObnK47lMnx01e2eOWmI1k88/lOdZOJiEitqlHwiYmJYcGCBee9vmDBAmJiKv6ln5mZSUhIyNVVJy6lSbAv18RXLIz45bbq792WkVfMA+9s4v31yVf0PhERkStVo+Dzyiuv8Nprr9G5c2ceeOABHnjgAbp06cLrr7/Oq6++CsDGjRu544477FrsN998Q0JCAr6+voSEhFx0ZejMzEyio6MxmUxkZ2fbtQa5tMpBzp9tOY5hGNV6z3Nf7Sb77CanRzPza602ERGRGgWfW2+9laSkJG666SaysrLIyspi2LBhJCYm8qtf/QqAyZMnM2/ePLsV+umnnzJu3Djuvfdetm/fzpo1a7j77rsveO79999Pp06d7PbZUn3DOkTg42XmYEY+O47nXPb873am8c3ONNvz5Czt7yYiIrXnigc3l5aWcuONN/L3v/+d2bNn10ZN5ykrK+PRRx9l7ty53H///bbX27U7f3uE+fPnk52dze9//3u+++67OqlP/ifAx4sh7SL4ansqn205TueY4Iueezq/hGe/3AVA19hgtiZnczRTwUdERGrPFbf4eHl5sWPHjtqo5aK2bNlCSkoKZrOZrl27EhkZybBhw9i1a1eV8/bs2cMLL7zAO++8g9lcva9WXFx83kKMcnVGdqsY5Pz1jjRKyi6+ltMLi/Zw6kwJrcIb8OKIjgAcyyqodheZiIjIlapRV9fYsWMvOLi5thw6dAiA5557jpkzZ7Jo0SJCQkLo378/WVlZQEWAueuuu5g7dy6xsbHVvvbs2bMJCgqyPSoHZ0vN9WnRiLAAC1n5Jazal3HBc5bvPcHnW1Mwm+Dl2zoTH+aPyQT5JeVk5Wv1bxERqR013rJi/vz59OjRg4ceeojHH3+8yqO6pk2bhslkuuQjMTHRtgL0jBkzGDVqFN27d2fhwoWYTCY+/vhjAKZPn07btm0ZO3bsFX2X6dOnk5OTY3scO3bsit4v5/P0MDP87BYWn205ft7xnMJSnvl8JwAP9o2nS0wwPl4eRAT6ABrnIyIitadGCxju2rWLbt26AbBv374qx0wmU7Wv88QTTzBhwoRLnhMfH09aWsXg11+O6bFYLMTHx9vWClqxYgU7d+7kk08+AbB1lzRq1IgZM2bw/PPPX/D6FosFi8VS7ZqlekZ2i+afPx1m+d6T5BSUEuTnZTv24jd7OJFbTHwjfx67oZXt9ZhQP9JyikjOKqBrrJZCEBER+6tR8Fm5cqVdPjwsLIywsLDLnte9e3csFgtJSUn06dMHqBhkfeTIEeLi4oCKWV+FhYW292zcuJH77ruP1atX07x5c7vUK9XXLiqQNhEBJKbnsWhnKmMSKv53+nFfBh9tOo7JBC/f1gkfLw/be2JD/dhwOItkDXAWEZFaUqPgU+nAgQMcPHiQ66+/Hl9fXwzDuKIWn+oKDAxk0qRJzJo1i5iYGOLi4pg7dy4At99+O8B54ebUqVMAtG3bluDgYLvXJJc3slsTXvo2kc+2pDAmIY4zxWVM/6yii2t876b0aBpa5fy4UD9AXV0iIlJ7ahR8MjMzGT16NCtXrsRkMrF//37i4+O5//77CQkJsS1iaE9z587F09OTcePGUVhYSEJCAitWrNDq0E5seJcmzPkukc1HT3M0M59/rD5ESnYhMaG+PH1j6/POj22o4CMiIrWrRoObH3vsMby8vEhOTsbPz8/2+h133MHixYvtVtwveXl58corr3DixAlyc3NZtmwZ7du3v+j5/fv3xzAMtfY4UHigD9e1aATAs1/u5r2fK8Zj/WlkJ/y8z8/cMWrxERGRWlajFp+lS5eyZMkSoqOjq7zesmVLjh49apfCpH4Y1S2a1ftP8ePZae13J8Ry7dkwdK7Ys8EnPbeIotLyKuN/RERE7KFGLT75+flVWnoqZWVlaYaUVDGkfTj+3hUBJirIh+nD2lz03Ib+3vh7e2AYkJJdeNHzREREaqpGwadv37688847tucmkwmr1crLL7/MgAED7FacuD4/b0/G9o7D4mnm5ds6E+DjddFzTSbT/7q7NLNLRERqQY26ul5++WUGDRrEpk2bKCkp4emnn2b37t1kZWWxZs0ae9coLm76sLY8cUNrvD0vn7NjQ/1ITM/TOB8REakVNWrx6dChA/v27aNPnz4MHz6c/Px8Ro4cydatW7VmjlxQdUIPQJxmdomISC2q8To+QUFBzJgxw561iNgGOCv4iIhIbahx8MnOzmbDhg2cPHnStpdWpXvuueeqCxP3pDE+IiJSm2oUfL7++mvGjBnDmTNnCAwMrLJas8lkUvCRGvtli09trQQuIiLuq0ZjfJ544gnuu+8+zpw5Q3Z2NqdPn7Y9srKy7F2juJHoED9MJigsLefUmRJHlyMiIvVMjYJPSkoKjzzyyAXX8hG5Gt6eZqKCfAGN8xEREfurUfAZOnQomzZtsnctIgDEhFYGn3wHVyIiIvVNjcb43HzzzTz11FPs2bOHjh074uVVdVG6W2+91S7FiXuKDfXj50NZJGdq9WYREbGvGgWfBx98EIAXXnjhvGMmk4ny8vKrq0rcWlxDf0BdXSIiYn81Cj7nTl8Xsaf/7dKuri4REbGvKxrjc9NNN5GTk2N7PmfOHLKzs23PMzMzadeund2KE/ekRQxFRKS2XFHwWbJkCcXFxbbnL730UpXp62VlZSQlJdmvOnFLlcHnRG4xRaXqNhUREfu5ouBjGMYln4vYQ4ifFwGWil7Y46fV6iMiIvZTo+nsIrXJZDLZxvkc1dYVIiJiR1cUfEwm03lbCGhLAakNGucjIiK14YpmdRmGwYQJE7BYLAAUFRUxadIk/P0rph//cvyPyNWIbajgIyIi9ndFwWf8+PFVno8dO/a8c7RBqdhDrHZpFxGRWnBFwWfhwoW1VYdIFerqEhGR2qDBzeKUfhl8NHtQRETsRcFHnFKTEF/MJigus5KRp7FjIiJiHwo+4pS8PMxEBVfs0n5U3V0iImInCj7itDTAWURE7E3BR5yWBjiLiIi9KfiI06pcy+eYgo+IiNiJgo84rcoWH43xERERe1HwEaelri4REbE3BR9xWpXBJyOvmMKScgdXIyIi9YGCjzitYD9vAn0qFhdXq4+IiNiDgo84NW1WKiIi9qTgI07Nmcf5pGYXMv+Hg2QXlDi6FBERqSYFH3FqMaHOO6X91aX7+NPiRH7zny2UW7WfmIiIK1DwEacWF+oPwNHMfAdXcr4tyacBWHswk7+tPODgakREpDoUfMSpOWtXV3ZBCYdP/S+Mvf79PtYfynRgRSIiUh0KPuLUKoPPsdOFWJ2oO2n78RwAmjb0Y2S3JlgNePTDbWTla7yPiIgzU/ARpxYV7IOH2URJmZWTecWOLsdmW3I2AJ1jgvnD8A7Eh/mTnlvEkx9vxzCcJ6CJiEhVCj7i1Dw9zDQJ9gWca5zP9uPZAHSJCcbf4slf7+qGt6eZFYknWfDTYccWJyIiF6XgI07P2cb5GIbBtmPZQEWLD0C7qECe/VU7AP60ONF2vDqy8kv4fs8Jysqtdq5URETOpeAjTs/ZprQfP11IVn4JXh4m2kUG2l4fmxDLTR0jKC03ePiDLeQWlV7yOkWl5fx91UH6vbySB97ZxN9WHqzt0kVE3J6Cjzi9uIbOtUt7ZWtO28hAfLw8bK+bTCZmj+xEdIgvx7IKmfbpjguO9zEMg6+2pzLo1VXM+S6RvOIyAL7anlIn9YuIuDOXCj7ffPMNCQkJ+Pr6EhISwogRI8475+2336ZTp074+PjQuHFjpkyZUveFil05W1dXZfDpcrab65eCfL34693d8DSb+HZnOv9Zn1zl+KYjWYz4v7U88sFWUrILiQj04cVfd8Dbw8zBjHz2n8irg28gIuK+PB1dQHV9+umnPPjgg7z00ksMHDiQsrIydu3aVeWcefPm8eqrrzJ37lwSEhLIz8/nyJEjjilY7CbWybq6tleO74kOvuDxLjHB/O7GNrz47V5eWLSH7nEh+Hl7MOe7RL7blQ6An7cHk/s154G+8fh6e/D9nhOsTMpg8a50WoYH1NE3ERFxPy4RfMrKynj00UeZO3cu999/v+31du3a2f58+vRpZs6cyddff82gQYNsr3fq1KlOaxX7qxzjc+pMCfnFZfhbHPfXtrTcys6UijV8usQGX/S8+/s0Y+3BU6xMymDcgg3kFJZQWm5gNsEdPWN47IZWNA7wsZ0/rEMkK5My+G5XOg8PalnbX0NExG25RFfXli1bSElJwWw207VrVyIjIxk2bFiVFp9ly5ZhtVpJSUmhbdu2REdHM3r0aI4dO3bJaxcXF5Obm1vlIc4lyNeLYD8vwPHdXUnpeRSXWQnw8aRZQ/+Lnmc2m3h1dBfCAy2cOlNMabnB9a3C+O7R65k9slOV0AMwuF04HmYTe9JySc50jpYtEZH6yCWCz6FDhwB47rnnmDlzJosWLSIkJIT+/fuTlZVlO8dqtfLSSy/x+uuv88knn5CVlcUNN9xAScnFV9OdPXs2QUFBtkdMTEydfCe5Ms4yzueX43vMZtMlzw3192bB+J6M7NaEf9/Xi3fu60XriAt3Y4X6e5PQLBSAxbvT7FqziIj8j0ODz7Rp0zCZTJd8JCYmYrVWrG8yY8YMRo0aRffu3Vm4cCEmk4mPP/4YAKvVSmlpKW+88QZDhw7lmmuu4YMPPmD//v2sXLnyojVMnz6dnJwc2+NyLUTiGM4ypf1y43vO1aFJEPNGd6Ffq7DLnntjhwgAFp8dByQiIvbn0DE+TzzxBBMmTLjkOfHx8aSlVfwL+JdjeiwWC/Hx8SQnV8yaiYyMPO+csLAwGjVqZDvnQiwWCxaLpaZfQepIZYvPUQd3A11qRtfVGto+gt9/uZstydmk5xQREeRz+TeJiMgVcWjwCQsLIyzs8v8S7t69OxaLhaSkJPr06QNAaWkpR44cIS4uDoDrrrsOgKSkJKKjowHIysri1KlTtnPEdcU5QVdXXlEpBzLOAP9bsdmewgN96BYbzJbkbJbuSeee3k3t/hkiIu7OJcb4BAYGMmnSJGbNmsXSpUtJSkpi8uTJANx+++0AtGrViuHDh/Poo4+ydu1adu3axfjx42nTpg0DBgxwZPliB84wpX3n8RwMA5oE+xIWUDuthMM6VLRcqrtLRKR2uETwAZg7dy533nkn48aNo2fPnhw9epQVK1YQEhJiO+edd94hISGBm2++mX79+uHl5cXixYvx8vJyYOViD5VjfI6fLqTceundz62XOV5T236xMWltGdq+YpzP+sNZZOVffFC+iIjUjEus4wPg5eXFK6+8wiuvvHLRcwIDA1mwYAELFiyow8qkLkQF++JpNlFSbiU9t4gmwb4YhsGxrEJ2peawOzWH3am57E7NJaeglHl3dOZXnaLsWsO25GygdoNPbEM/2kUGsictl+/3nGB0T80yFBGxJ5cJPuLePMwmokN8OZJZwB8X7SErv4Q9abnkFZVd8PzH/ruNIF8v+ra8/Biy6tp+tsWnNsb3/NKwDhHsScvlu11pCj4iInbmMl1dIrFnFwz8blc66w9nkVdUhreHmY5NgrizZwx/GN6eTyf35pbOUZSWGzz07mbb9POrlZZTyIncYjzMJjo0Cbz8G65C5bT2NQcyL7vDu4iIXBm1+IjL+E3/5nh7mIkJ9aV9VBDtowJp0bgBXh5V83vHJsFkF5Swev8p7n17Ix9P6k3zsAZX9dmVAapVeAB+3rX7f5uW4QE0D/PnYEY+KxNPMrxLk1r9vHMZhsGetFxahQecd29FRFyd/qsmLuOa+Ib8c3wPZt3Sntu6R9M2MvCCP8zenmbmj+1Op+ggsvJLuGfBBk7kFl3VZ287dnZ/rlru5qrkyMUM//jNXm5+4yf+suJAnX+2iEhtU/CReqmBxZOFE3oS38iflOxC7lmwgZyCmncbbTt2GoAuMUH2KvGSKqe1/5CUQWFJeZ18JlQErQU/HQbgi60pGEbtzJATEXEUBR+ptxo2sPDv+3oRHmgh6UQe9/97Y41CRLnVYOfxyhafkMucbR/towKJDvGlsLScVfsy6uQzkzMLeOqT7f97nlVAYnpenXy2iEhdUfCRei0m1I9/39eLQB9PNh09zdT3t1BWbr2iaxw4eYb8knL8vT1o0fjqxgpVl8lk4sb2ld1dtb9paXFZOVM/2EJeURndYoNte4st3X2i1j9bRKQuKfhIvdcmIpAFE3pi8TSzPPEk0z/beUVdOJUDmztGB+FxmR3Z7alynM/yvScpKbuysHalZn+byI7jOQT7efGXu7txc8eKrrale7SCtIjULwo+4hZ6Ng3lr3d3w8Ns4uPNx/nT4qRqv3dr5Y7sdTSwuVK32BDCAizkFZex9uCpWvuc73am8fbaIwC8entnmgT7MqhtY8wm2J2ay/HTjt0YVkTEnhR8xG3c0C6c2b/uCMDfVx3ko43HqvW+yhafrnUcfMxmE0PbhwO1N7srObOApz/ZAcBD18czqG3F5zVsYKFH01AAlu1Rd5eI1B8KPuJWRveM4fEbWgHw3Ne7OZqZf8nzC0vKSTpRMcC3rlt8AG5sX9nldOKye5RdqeKycqa8v4W84jK6x4Xw5NDWVY4PaVcRgjTOR0TqEwUfcTtTBrQgoVkoBSXlPPbfbZcc7LwrNYdyq0F4oIXIIN86rLJCQnwowX5eZOWXsOFwll2v/dI3e9mZcnZcz11dz1sTaUi7ijFGG45kcVobpopIPaHgI27Hw2zi1dGdaWDxZEtyNm/+eOii51ZuTNo5OrhuijuHl4eZwWe7n5bstl9317c70/j3uqMAvDa6C1HB54e62IZ+tIkIoNxqsDzxpN0+W0TEkRR8xC1Fh/jx/K3tAXht2T52peRc8LxtZzcm7RIbXEeVnW/YL1Zxttqhu+toZj6/OzuuZ1K/5gxo0/ii5w45O6V+qR1Dl4iIIyn4iNsa2a0JwzpEUGY1+O1/t1FUev7ihpUtPl0c1OIDcF2LRvh7e5CeW2TbIb6mikr/N66nR1wITwxpdcnzKwdX/7i/bleQFhGpLQo+4rZMJhMv/bojjQMsHDh5hjnfJVY5npFXTEp2ISZTxRo+juLj5cHAs91dn29NuaprvbnqELtScgnx8+KNC4zrOVe7yECaBPtSVGpl9f66WUFaRKQ2KfiIWwvx9+bl2zoB8PbaI1V+3CunsbcIa0CAj5cjyrO5s2cMAB9tOkbmmeIaXaOgpIyFayv24Xru1vYXHNdzLpPJxJCzrT5LNa1dROoBBR9xe/1bN2bcNXEAPPnxdrILKmYwVXYr1dWO7JdybfOGdGwSRFGplXfODkq+Uh9uOEZ2QSlNG/rxq05R1X5f5eyu5XtPXPF2HyIizkbBRwSYflMb4hv5cyK3mGe/3A3ANget2HwhJpOJh/rFA/DvdUcoKCm7oveXllv55+qK2WsTr29+RVtv9GwaQrCfF6cLStl45PQVfa6IiLNR8BEB/Lw9mXdHFzzMJr7ensoXW1NsXV3O0OIDMKxDJHEN/cguKOW/1Vx1utJX21JJzSkiLMDCyG5Nrui9nh5mBrWp7O7S7C4RcW0KPiJndYkJ5pGBLQH43ac7yC0qw+JppnVEgIMrq+BhNvFg34pWn3+uPkxpNbudrFaDv686CMD9fZrh4+VxxZ9tG+ez+8QVbfAqIuJsFHxEfmHKgOZ0jgmm+Oxu6B2bBF125lNduq17NI0aeJOSXciiHanVes/yxJPsP3mGAIsndyfE1uhzr28Zho+XmZTsQvak5dboGiIizsB5/osu4gQ8Pcy8NrozvmdbRZxhfM8v+Xh5cO91zYCKqemXa30xDIP5PxwAYGzvOAJrODvN19uD61uGAdq7S0Rcm4KPyDniwxrwyu2daR8VyOgeMY4u5zxjE+Lw9/YgMT2PH5IuvbbOxiOn2ZKcjbenmXuva3pVn2tbxVnT2kXEhSn4iFzAzZ0i+eaRvk4zvueXgvy8uKtXRZdV5didi6ls7bmtezSNA3yu6nMHtWmM2QR703I5llVwVdcSEXEUBR8RF3R/32Z4eZhYfziLrckXnmK+Ny2XlUkZmE0w8eyg6KsR4u9Nr2ahgH03TBURqUsKPiIuKDLIl+FdKqalX6zV582zr9/UMZKmjfzt8rmVixmqu0tEXJWCj4iLmnR2QcOle05wMONMlWPHsgr4ekfa2fOa2+0zK6e1bzqSVeOtM0REHEnBR8RFtWgcwOC24RgGvLXqUJVj/1x9iHKrQd+WjejQxH4brEaH+NE+KhCrUTFNXkTE1Sj4iLiwyf0rWn0+35rCidwiAE6dKebDsys7T+5vv9aeSrbuLk1rFxEXpOAj4sK6x4XSIy6EknIr//qpYuf1f689QnGZlc7RQfSOb2j3z6zs7lq9P+OK9ww7V1FpOd/sSCOnoNQepYmIXJaCj4iLqxzD85/1yaTnFNl2b5/cvzkmU/U3I62uNhEBxIT6Ulxm5cd9l15H6FJKy61Mem8zU97fwpgFP1NSpp3fRaT2KfiIuLiBbRrTsnEDzhSXMXbBenIKS4kP87d1SdmbyWRi6Nlrz/4ukbScwiu+hmEYTPt0p20Bxl0pucxdkmjXOkVELkTBR8TFmc0mHjrb6nPgZMXsroeuj8dstn9rT6WJ18cTE+rL0cwC7nzr5ysOPy8vSeLTLcfxMJt4oE/FFhz/WH2YlUkaMC0itUvBR6QeuLVzFJFBFSszhwdaGNG1Sa1+XuNAHz6c2LtG4WfhmsPM/6FijaHZIzsy81ftGN87DoAnP9rOybyiWqtbRETBR6Qe8PY088SQ1gA8cUNrLJ4etf6ZTYJ9q4Sfu6oRfr7ensoLi/YA8NTQ1ra90Kbf1JY2EQFk5pfwxEfbsVovvfmqiEhNKfiI1BO3dY8m6Y83Mrpn3W2s2iTYlw8evIaYUF+OnA0/6TkXbrFZe+AUT3y0HcOAe3rH8ZtfTLX38fLgr3d3xcfLzOr9p3hr9aELXkNE5Gop+IjUI3XR0nOu6BA/PnjwGqJDKsLPnW+tOy/87E7NYeK7mykpt3JTxwhm3dL+vBlnLRoH8Nwt7QF4ZUkS245l19VXEBE3ouAjIlctOsSPDydeOPwcyypgwsKNnCkuI6FZKPNGd8HjIgOv7+gZw80dIymzGjzywVbyirS+j4jYl4KPiNjFueHnrn/8zJ7UXO751wYy8oppExHAP8b3wMfr4q1SJpOJl0Z2pEmwL8lZBcz8YheGofE+ImI/JkP/VakiNzeXoKAgcnJyCAwMdHQ5Ii7n+OmKWV7HT/9voHOTYF8++821hAf6VOsam49mMfrNnym3Grxye2du6x590XNLy61sP5bNluTT+Hl70iTYl8hgHyKDfAn08ayVRRzt5URuESF+3nh76t+gIlerur/fCj7nUPARuXrHsipafI6fLiTEz4tPJl9L87AGV3SNv67YzytL9+Hn7cGih/sQf/b9hmFw6FQ+P+0/xer9p/j5UCZnii+8dUYDiyeRQT5EBvvSJNiHqCBfBrUNp12U4/+//c66I8z6ajexoX68PKoTCbWwvYiIO1HwqSEFHxH7SMku5J11RxjRpQltI6/8/0vlVoOx/1zPukOZtI8K5KF+zflpfwY/7T9F6jmDp4P9vEhoFkq51SA1u4jUnEKyL7H/13UtGvJA33j6tQyr1YUeL+bddUd49svdVV6bcG1Tnr6xNX7ennVej0h9oOBTQwo+Is4jPaeIYX/+kdPnhBhvDzPd40Lo07IRfVs2on1U0HkDpgtKykjNLiItp5DU7EJSs4tITM/l+70nKT+7TlCLxg14oE8zRnRtcsmxR/b0y9Bzf59mFJSU8cGGYwAVrT+3deKaq2z9MQyDvOIyTueXkJVfwumCEnIKS+nVrCFNgn2v+juIOKN6GXy++eYbXnjhBXbs2IGPjw/9+vXjiy++sB3fuHEj06ZNY/PmzZhMJnr16sXLL79M586dq/0ZCj4izuWHpJNMfX8r0SG+9GnRiD4tG5HQrCG+3jULKseyCnh77RH+u/GYrYusob8343rHMe6aOBo2sNiz/Cre/fkoz36xC6jYVmTasDaYTCZ+3JfBtE932FqyxveO4+kb2+BvuXTrT3pOEav3Z7DuYCapOYWczi8lq6CE0/kllF1gEcjoEF+WPdavxvdOxJnVu+Dz6aef8uCDD/LSSy8xcOBAysrK2LVrF6NHjwbgzJkzxMXFceuttzJt2jTKysqYNWsWP/30E8eOHcPLy6tan6PgI+IecotK+e+GYyxcc9gWOCyeZkZ2a8IN7cJpExFIZJCP3QZH/zL0TLw+nulnQ0+lvKJSXvo2kQ82JAMQE+rLy6M607v5/1p/CkvK2XAkix/3ZbB6fwb7Tpy55Gf6eXsQ4udNqL83x04XkF1QyqR+zZk2rI1dvpOIM6lXwaesrIymTZvy/PPPc//991/wnE2bNtGzZ0+Sk5OJialYuXbnzp106tSJ/fv306JFi2p9loKPiHspLbfy3a50/rn6EDuO51Q5FujjSZuIQNpEBtA6IoA2EYG0jgigwWVaYs713s9HmXmJ0PNLq/dn8LtP/tf6M+6aOGJCfflx3yk2HMmipMxqO9dkgk7RwfRt0YhWEQGE+nkT4u9FqL83IX7eVbrvlu05wYPvbMLTbOKbR/rSOiLgir6DiLOrV8Fnw4YNJCQk8K9//Ys33niD9PR0unTpwty5c+nQoQMAeXl5NGvWjKlTp/LMM89QXl7O9OnTWbp0KTt27MDT88L/oSouLqa4uNj2PDc3l5iYGAUfETdjGAYbj5zmww3J7ErN4WBGvm0s0LliQ/3o1SyUgW0a07dlIwJ8Lt6i/J/1R5nxeUXoebBvM565qe1lW5HObf35pcggH65vGUbfVo24rnkjQvy9q/0dJ76ziaV7TtAjLoSPHupdo4HdB07m0dDfckWfK1IX6lXw+fDDD7nrrruIjY1l3rx5NG3alFdffZWlS5eyb98+QkNDAdi1axcjRozg8OHDALRs2ZIlS5YQFxd30Ws/99xzPP/88+e9ruAj4t6Ky8o5eDKfxPRcEtPzKh5puZzMK65ynqfZZAtBA9o0Jr6Rvy3Y/DL0PNCnGTNuvnzo+aWf9p/ijRX78fP24PqWYVzfqhHNwxrUuPstNbuQG+atIr+knDkjO3Jnr9grev+X21J49MNttApvwLeP9MXTQ+sPifNwieAzbdo0/vSnP13ynL1797JlyxbGjBnDm2++ycSJE4GKlpro6Gj++Mc/8tBDD1FYWEj//v1p06YNU6dOpby8nFdeeYXExEQ2btyIr++FZzKoxUdErkRWfgm7UnJYtS+DlYknOXQqv8rxuIZ+DGjdmCBfL/68fD9Qs9BTW/65+hB//GYvQb5eLH+iH42qOZh7w+Esxv5zPSXlFV1ts0d25K4rDE4itcklgk9GRgaZmZmXPCc+Pp41a9YwcOBAVq9eTZ8+fWzHEhISGDx4MC+++CILFizgmWeeIS0tDbO54l8hJSUlhISEsGDBAu68885q1aQxPiJyJY6cymdF4klWJp1k/aEsWzCodH+fZsx0ktADUFZu5da/rmFPWi4juzZh3h1dLvueQxlnGDl/LdkFpTQJ9iUlu5DwQAs/PDlAM8TEaVT399uhK2WFhYURFhZ22fO6d++OxWIhKSnJFnxKS0s5cuSIrRuroKAAs9lc5T8ulc+tVusFrysicrWaNvLnvj7NuK9PM/KLy/jpwClWJp7k50OZ/KpTFE8MaeU0oQfA08PMSyM78uv/W8NnW1O4rXs017ZodNHzM88Uc+/bG8kuKKVzTDDv3NeLm/68mpTsQt5ee4TJ/ZvXYfUiV88lOmgDAwOZNGkSs2bNYunSpSQlJTF58mQAbr/9dgBuuOEGTp8+zZQpU9i7dy+7d+/m3nvvxdPTkwEDBjiyfBFxE/4WT4a2j2DOqE788NQAnhza2qlCT6UuMcGMu6biH40zv9hFcVn5Bc8rKi3nwXc2cTSzgJhQX/55Tw+CfL14/IZWAMz/4QA5l1ghW8QZuUTwAZg7dy533nkn48aNo2fPnhw9epQVK1YQEhICQJs2bfj666/ZsWMHvXv3pm/fvqSmprJ48WIiIyMdXL2IiHN5cmhrwgIsHDqVz/wfDp533Go1eOKj7WxJzibQx5OFE3oSFlAxHmhE1ya0Dg8gt6iM+avOf6+IM3OJWV11SWN8RMRdfL09lYc/2Iq3h5nFv+1r2wgWYM53ifx91UG8PEy8c19ClYUUAZbvPcH9/96ExdPMqqcGEBHkU9fli1RR3d9vl2nxERER+/pVp0iubxVGSbmVmV/sovLfwR9sSObvZ1ty/jSq03mhB2Bgm8b0iAuhuMxqm70m4goUfERE3JTJZOKPwztg8TSz9mAmX2xLYdW+DNsq078d3JKR3aIv+t7KrS8+2nSMgxmX3j5DxFko+IiIuLHYhn48MqglAH9YtJcp/9lCudVgZLcmPHr29Yvp0TSUwW0bU241eHVpUl2UK3LVFHxERNzcg33jadm4AVn5JZwpLqN3fEPmjOxUrRlpFTPX4Nud6Ww/ll37xYpcJQUfERE35+1ZsbaPp9lEy8YN+PvY7nh7Vu/noU1EIL/u2gSAl5ck1maZInbh0AUMRUTEOfRsGsqqpwcQ4ueFn/eV/TQ8NrgVi7anseZAJqv3Z9C35eUXprWHcqvBwjWHKSgpp0OTQDpEBdE4ULPL5NIUfEREBIAmwRfe0/ByYkL9GHNNLAvXHOFPixO5rnmjS+78XrnCdXwjf1qGB9ToM4vLyvnth9v4bld6ldfDAix0iAqkQ5Mg2kcF0T4qkOgQX6dcSFIcQ8FHRESu2tQBLfho4zF2peTyzc40bukcVeV4fnEZyxNP8s2OVH5IyqC4zIq3h5k/jGjPHT2vbLPTM8VlPPTuJtYcyMTbw8wN7cPZfyKPAyfPkJFXzMqkDFYmZdjOD/L1Yuw1sTw5xDlX0pa6peAjIiJXrWEDCw9eH8/r3+/n1aVJ3NghguIyK8v3nuDbnWm2sFMpxM+L0wWl/O7TnWw/nsOsW9ph8bz8hqdZ+SVMWLiBHcdz8Pf24B/39LDtNVZYUs7e9Fx2p+SwKyWXXak57DuRR05hKX9beZAOUUEM66iV/N2dVm4+h1ZuFhGpmTPFZfR7eSWZ+SV0jg4iMT2vSthp2tCPmzpGcnOnSNpGBPJ/Pxzg1WX7MIyK/cPmj+1GZNDFu9tSswsZu2A9hzLyCfX35u17e9IpOviSNRWXlTNv6T7e/PEQYQEWvn+sH0F+Xvb6yuJEqvv7reBzDgUfEZGae3vNYZ77eo/tedOGftzcKZKbOkbSLjLwvK6mH5JO8uiH28gpLKVRA2/+enc3rok/f6XoAyfPcM+C9aTmFBEV5MM79yfQonGD8867kKLScm56YzWHMvK5s2cMc0Z1urovKU5JwaeGFHxERGqutNzKa8v2YTaZGNYx4oJh51zJmQU89N5m9qbl4mE2MeOmttx7XVPb+7Yfy2bCwg2cLiileZg/796fQNQVDsTecDiL0W+uA+CDB6+54DYc4toUfGpIwUdEpO4VlpQz7bMdfLktFYDhXaKYM7ITW5JPM/GdTeSXlNM5OoiF9/Yi1N+7Rp/xzOc7eX99Ms0a+fPdo33x8br8mCJxHQo+NaTgIyLiGIZh8PbaI/zxm72UWw3iw/w5nlVISbmV61o05M1xPWhgqfmcnNyiUga/uoqTecX8pn9znr6xjR2rF0fT7uwiIuJSTCYT917XjPcfSKBRA28OZeRTUm7lxvYR/GtCz6sKPQCBPl78YUQHAN788RB7UnPtUba4GAUfERFxKgnxDVn0cF9u7hTJb/o3529julVrqnt1DG0fwY3tIyi3Gkz7bAflVnV6uBsFHxERcToRQT787e5uPH1jGzwusQp0TTw/vD0BPp7sOJ7DwjWH7XptcX4KPiIi4lbCA3145qa2ALy6dB/HsgocXJHUJQUfERFxO3f0iCGhWSiFpeXM+GIXmufjPhR8RETE7ZjNJmaP7Ii3p5kf92XwxbYUR5ckdUTBR0RE3FJ8WAMeHdQSgBe+3kPmmWIHVyR1QcFHRETc1sTr42kTEcDpglKe/XKXZnm5AQUfERFxW14eZuaM6oTZBN/uTGfSe5spLCl3dFlSixR8RETErXWJCeYvd3XD29PMsj0nuPOtdWTk1U631960XF78Zg//+knT6B3l6pbBFBERqQdu7hRJ40ALD76zie3Hcxg5fw1v39uL5mHV2wH+UgpLyvl6RyofbEhma3K27fU2EQFc26LRVV9froz26jqH9uoSEXFfhzLOMGHhRpKzCgjy9eIf9/SgV7PQGl0rMT2XD9Yn89nWFPKKygDwNJuIDfXj0Kl82kUG8vXDfey+QKO70ialNaTgIyLi3jLPFPPAO5vYmpyNt4eZV0Z35tbOUdV6b2FJOYvOtu5s+UXrTmyoH3f2iuH27jF4mE30m7uSvKIyXr6tE6N7xNTSN3EvCj41pOAjIiJFpeX89sNtLN6dDsDTN7Zmcr/mmExVW2cKS8rZknya9YcyWX84i63HsikpswIVrTtD2odzV69YrmveCPMvWnb+8eMhXvx2L40DLKx8sj/+V7kBqyj41JiCj4iIAJRbDV76di8Lzg5EvqtXLL+7sTVbj2Wz4XAW6w9lsuN4DmXnTIGPCfXlzp6x3N4jmsYBPhe8dnFZOTfM+5HkrAIeGdSSx29oVevfp75T8KkhBR8REfmlhWsO88KiPVzs1zIyyIeEZqH0ataQhPhQ4hv5n9cydCHf7Uxj8n+24ONlZuWT/YkM8rVz5e6lur/falsTERG5hHuva0aTYF8e+XArRaVWYkP9zgadUK6Jb0h0iG+1gs65buwQQa+moWw4ksXcJUnMG93F/sXLedTicw61+IiIyIVk5BVTbjWICLpw91VN7Dieza1/XQPAV1Ovo1N0sN2u7W6q+/utBQxFRESqISzAYtfQA9ApOpiRXZsA8Mdv9mqX+Dqg4CMiIuJATw5tjY+XmQ2Hs1iy+4Sjy6n3FHxEREQcKCrYl4l94wGY/d1eisu0V1htUvARERFxsIf6NScswMLRzALeXXfU0eXUawo+IiIiDuZv8eSpIa0B+PPy/WTllzi4ovpLwUdERMQJjOoeTdvIQPKKynhj+X5Hl1NvKfiIiIg4AQ+ziZk3twXg3Z+PcuDkGQdXVD8p+IiIiDiJ61o0YnDbxpRbDZ7+ZDtLd6eTU1jq6LLqFS1geA4tYCgiIo50MOMMN77+I6XlFT/PZhN0aBJE7/iG9G7ekJ5NQ7Wp6QVor64aUvARERFH23Ysm082H2PtwUwOZeRXOeZpNtE5JphrmzdkYJvGdIkJrtGWGfWNgk8NKfiIiIgzOZFbxLqDmaw9eIo1BzJJyS6scrxTdBATrm3KzZ0isXh6OKhKx6tXweeHH35gwIABFzy2YcMGevbsCcCOHTuYMmUKGzduJCwsjIcffpinn376ij5LwUdERJzZsawC1h48xer9p1i65wQlZVYAGjWwMCYhljHXxNI4wL5ba7iCehV8SkpKyMrKqvLas88+y/Llyzl48CAmk4nc3FxatWrF4MGDmT59Ojt37uS+++7j9ddfZ+LEidX+LAUfERFxFVn5JXywIZl31x0lPbcIAC8PEzd3jOTe65rROSb4qj/DMAw+3ZJCbmEpE65titnsnN1q9Sr4nKu0tJQmTZrw8MMP8+yzzwIwf/58ZsyYQXp6Ot7e3gBMmzaNL774gsTExGpfW8FHRERcTWm5lSW701m45gibj562vd41NphJ/ZoztH1Eja878/Nd/HfTMQDGXhPLH4Z3cMoxRfV6d/avvvqKzMxM7r33Xttr69at4/rrr7eFHoChQ4eSlJTE6dOnL3QZAIqLi8nNza3yEBERcSVeHmZ+1SmKTydfy1dTr2NktyZ4e5jZmpzNQ+9u5vGPtpFfXHZF18wtKuW+tzfy303HMJvAZIL3fk5mzuJEl95F3iWDz4IFCxg6dCjR0dG219LT0wkPD69yXuXz9PT0i15r9uzZBAUF2R4xMTG1U7SIiEgd6BQdzLzRXVgzbSCT+jXHbILPtqRwy19+YldKTrWukZJdyO3z17F6/yl8vTx4a1wPXvp1RwDeXHWIv608UJtfoVY5NPhMmzYNk8l0yce53VTHjx9nyZIl3H///XapYfr06eTk5Ngex44ds8t1RUREHCkswMK0YW34cGJvIoN8OHQqn5H/t5aFaw5fssVm5/EcRvxtDUkn8mgcYOHjSb0Z3C6cu3rF2laWfmXpPv710+G6+ip25dAVkJ544gkmTJhwyXPi4+OrPF+4cCENGzbk1ltvrfJ6REQEJ06cqPJa5fOIiIv3bVosFiwWyxVULSIi4jp6NQvl20f68vSnO1i25wTPf72HNQdO8fJtnQn1965y7rI9J3jkg60UlpbTJiKAf03oSVSwr+34A33jOVNcxuvf7+eFRXtoYPFkdE/X6ilxaPAJCwsjLCys2ucbhsHChQu555578PLyqnKsd+/ezJgxg9LSUtuxZcuW0bp1a0JCQuxat4iIiCsJ8ffmrXHdeffno/xx0V6+33uSm/68mtfv7MI18Q0BWLjmMC8s2oNhQN+Wjfi/Md0I8PE671qPDmpJfnEZ/1h9mN99tgNfbw9u6RxV11+pxlxqVtfy5csZPHgwe/fupU2bNlWO5eTk0Lp1a4YMGcLvfvc7du3axX333cdrr72m6ewiIiJn7U7N4eEPtnIoIx+zCR4e2JKcwlLeXnsEgLt6xfDC8A54eVx8NIxhGMz4Yhfvr0/G02zi72O7M7hd+EXPP/e9tTErrF5OZ7/77rs5evQoa9asueDxXy5g2KhRIx5++GF+97vfXdFnKPiIiEh9l19cxqyvdvPJ5uNVXp82rA0PXR9frWBitRo8/tE2vtiWirenmYUTenJdi0ZVziktt3Iw4wy7U3LZnZrL7tQckk7ksXbaQPy87dvpVC+DT11Q8BEREXfxxdYUZny+k1KrwbzRnflVpyvrsiott/Kb/2xh2Z4T+Hl7MHtkR3KLytiTmsPu1FwS0/NsK0v/0qeTr6V7nH2HoSj41JCCj4iIuJPsghJKyq013uaiuKycB/69idX7T13weAOLJ+0iA2kXFUj7qEDaRwXRMrzBJbvSaqK6v9/a115ERMSNBft5X/6kS7B4evDmuO5M+c8W9qTl0ibifwGnfVQgsaF+TrXNhYKPiIiIXBU/b08W3tvL0WVUi0uu3CwiIiJSEwo+IiIi4jYUfERERMRtKPiIiIiI21DwEREREbeh4CMiIiJuQ8FHRERE3IaCj4iIiLgNBR8RERFxGwo+IiIi4jYUfERERMRtKPiIiIiI21DwEREREbeh4CMiIiJuw9PRBTgbwzAAyM3NdXAlIiIiUl2Vv9uVv+MXo+Bzjry8PABiYmIcXImIiIhcqby8PIKCgi563GRcLhq5GavVSmpqKgEBAZhMJrtdNzc3l5iYGI4dO0ZgYKDdritV6T7XHd3ruqH7XDd0n+tGbd5nwzDIy8sjKioKs/niI3nU4nMOs9lMdHR0rV0/MDBQ/6eqA7rPdUf3um7oPtcN3ee6UVv3+VItPZU0uFlERETchoKPiIiIuA0FnzpisViYNWsWFovF0aXUa7rPdUf3um7oPtcN3ee64Qz3WYObRURExG2oxUdERETchoKPiIiIuA0FHxEREXEbCj4iIiLiNhR86sjf/vY3mjZtio+PDwkJCWzYsMHRJbm0H3/8kVtuuYWoqChMJhNffPFFleOGYfD73/+eyMhIfH19GTx4MPv373dMsS5s9uzZ9OzZk4CAABo3bsyIESNISkqqck5RURFTpkyhYcOGNGjQgFGjRnHixAkHVeya5s+fT6dOnWyLuvXu3ZvvvvvOdlz3uHbMmTMHk8nEb3/7W9trutdX77nnnsNkMlV5tGnTxnbc0fdYwacO/Pe//+Xxxx9n1qxZbNmyhc6dOzN06FBOnjzp6NJcVn5+Pp07d+Zvf/vbBY+//PLLvPHGG/z9739n/fr1+Pv7M3ToUIqKiuq4Ute2atUqpkyZws8//8yyZcsoLS1lyJAh5Ofn28557LHH+Prrr/n4449ZtWoVqampjBw50oFVu57o6GjmzJnD5s2b2bRpEwMHDmT48OHs3r0b0D2uDRs3buTNN9+kU6dOVV7XvbaP9u3bk5aWZnv89NNPtmMOv8eG1LpevXoZU6ZMsT0vLy83oqKijNmzZzuwqvoDMD7//HPbc6vVakRERBhz5861vZadnW1YLBbjgw8+cECF9cfJkycNwFi1apVhGBX31cvLy/j4449t5+zdu9cAjHXr1jmqzHohJCTE+Oc//6l7XAvy8vKMli1bGsuWLTP69etnPProo4Zh6O+zvcyaNcvo3LnzBY85wz1Wi08tKykpYfPmzQwePNj2mtlsZvDgwaxbt86BldVfhw8fJj09vco9DwoKIiEhQff8KuXk5AAQGhoKwObNmyktLa1yr9u0aUNsbKzudQ2Vl5fz4Ycfkp+fT+/evXWPa8GUKVO4+eabq9xT0N9ne9q/fz9RUVHEx8czZswYkpOTAee4x9qktJadOnWK8vJywsPDq7weHh5OYmKig6qq39LT0wEueM8rj8mVs1qt/Pa3v+W6666jQ4cOQMW99vb2Jjg4uMq5utdXbufOnfTu3ZuioiIaNGjA559/Trt27di2bZvusR19+OGHbNmyhY0bN553TH+f7SMhIYG3336b1q1bk5aWxvPPP0/fvn3ZtWuXU9xjBR8RqZYpU6awa9euKn31Yj+tW7dm27Zt5OTk8MknnzB+/HhWrVrl6LLqlWPHjvHoo4+ybNkyfHx8HF1OvTVs2DDbnzt16kRCQgJxcXF89NFH+Pr6OrCyCurqqmWNGjXCw8PjvBHrJ06cICIiwkFV1W+V91X33H6mTp3KokWLWLlyJdHR0bbXIyIiKCkpITs7u8r5utdXztvbmxYtWtC9e3dmz55N586d+fOf/6x7bEebN2/m5MmTdOvWDU9PTzw9PVm1ahVvvPEGnp6ehIeH617XguDgYFq1asWBAwec4u+zgk8t8/b2pnv37ixfvtz2mtVqZfny5fTu3duBldVfzZo1IyIioso9z83NZf369brnV8gwDKZOncrnn3/OihUraNasWZXj3bt3x8vLq8q9TkpKIjk5Wff6KlmtVoqLi3WP7WjQoEHs3LmTbdu22R49evRgzJgxtj/rXtvfmTNnOHjwIJGRkc7x97lOhlC7uQ8//NCwWCzG22+/bezZs8eYOHGiERwcbKSnpzu6NJeVl5dnbN261di6dasBGPPmzTO2bt1qHD161DAMw5gzZ44RHBxsfPnll8aOHTuM4cOHG82aNTMKCwsdXLlrmTx5shEUFGT88MMPRlpamu1RUFBgO2fSpElGbGyssWLFCmPTpk1G7969jd69ezuwatczbdo0Y9WqVcbhw4eNHTt2GNOmTTNMJpOxdOlSwzB0j2vTL2d1GYbutT088cQTxg8//GAcPnzYWLNmjTF48GCjUaNGxsmTJw3DcPw9VvCpI3/5y1+M2NhYw9vb2+jVq5fx888/O7okl7Zy5UoDOO8xfvx4wzAqprQ/++yzRnh4uGGxWIxBgwYZSUlJji3aBV3oHgPGwoULbecUFhYav/nNb4yQkBDDz8/P+PWvf22kpaU5rmgXdN999xlxcXGGt7e3ERYWZgwaNMgWegxD97g2nRt8dK+v3h133GFERkYa3t7eRpMmTYw77rjDOHDggO24o++xyTAMo27alkREREQcS2N8RERExG0o+IiIiIjbUPARERERt6HgIyIiIm5DwUdERETchoKPiIiIuA0FHxEREXEbCj4iIucwmUx88cUXji5DRGqBgo+IOJUJEyZgMpnOe9x4442OLk1E6gFPRxcgInKuG2+8kYULF1Z5zWKxOKgaEalP1OIjIk7HYrEQERFR5RESEgJUdEPNnz+fYcOG4evrS3x8PJ988kmV9+/cuZOBAwfi6+tLw4YNmThxImfOnKlyzr/+9S/at2+PxWIhMjKSqVOnVjl+6tQpfv3rX+Pn50fLli356quvbMdOnz7NmDFjCAsLw9fXl5YtW54X1ETEOSn4iIjLefbZZxk1ahTbt29nzJgx3HnnnezduxeA/Px8hg4dSkhICBs3buTjjz/m+++/rxJs5s+fz5QpU5g4cSI7d+7kq6++okWLFlU+4/nnn2f06NHs2LGDm266iTFjxpCVlWX7/D179vDdd9+xd+9e5s+fT6NGjeruBohIzdXZdqgiItUwfvx4w8PDw/D396/yePHFFw3DqNgxftKkSVXek5CQYEyePNkwDMN46623jJCQEOPMmTO24998841hNpuN9PR0wzAMIyoqypgxY8ZFawCMmTNn2p6fOXPGAIzvvvvOMAzDuOWWW4x7773XPl9YROqUxviIiNMZMGAA8+fPr/JaaGio7c+9e/eucqx3795s27YNgL1799K5c2f8/f1tx6+77jqsVitJSUmYTCZSU1MZNGjQJWvo1KmT7c/+/v4EBgZy8uRJACZPnsyoUaPYsmULQ4YMYcSIEVx77bU1+q4iUrcUfETE6fj7+5/X9WQvvr6+1TrPy8urynOTyYTVagVg2LBhHD16lG+//ZZly5YxaNAgpkyZwiuvvGL3ekXEvjTGR0Rczs8//3ze87Zt2wLQtm1btm/fTn5+vu34mjVrMJvNtG7dmoCAAJo2bcry5cuvqoawsDDGjx/Pe++9x+uvv85bb711VdcTkbqhFh8RcTrFxcWkp6dXec3T09M2gPjjjz+mR48e9OnTh//85z9s2LCBBQsWADBmzBhmzZrF+PHjee6558jIyODhhx9m3LhxhIeHA/Dcc88xadIkGjduzLBhw8jLy2PNmjU8/PDD1arv97//Pd27d6d9+/YUFxezaNEiW/ASEeem4CMiTmfx4sVERkZWea1169YkJiYCFTOuPvzwQ37zm98QGRnJBx98QLt27QDw8/NjyZIlPProo/Ts2RM/Pz9GjRrFvHnzbNcaP348RUVFvPbaazz55JM0atSI2267rdr1eXt7M336dI4cOYKvry99+/blww8/tMM3F5HaZjIMw3B0ESIi1WUymfj8888ZMWKEo0sRERekMT4iIiLiNhR8RERExG1ojI+IuBT1zovI1VCLj4iIiLgNBR8RERFxGwo+IiIi4jYUfERERMRtKPiIiIiI21DwEREREbeh4CMiIiJuQ8FHRERE3IaCj4iIiLiN/wdq9ngHJF+zaQAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "cudaq.set_target('nvidia')\n",
    "start_time = timeit.default_timer()\n",
    "result = minimize(cost,\n",
    "                  x0,\n",
    "                  method='COBYLA',\n",
    "                  callback=callback,\n",
    "                  options={'maxiter': 50})\n",
    "end_time = timeit.default_timer()\n",
    "\n",
    "print('UCCSD-VQE energy =  ', result.fun)\n",
    "print('Total number of qubits = ', qubit_count)\n",
    "print('Total number of parameters = ', parameter_count)\n",
    "print('Total number of terms in the spin hamiltonian = ',\n",
    "      spin_ham.term_count)\n",
    "print('Total elapsed time (s) = ', end_time - start_time)\n",
    "\n",
    "plt.plot(exp_vals)\n",
    "plt.xlabel('Epochs')\n",
    "plt.ylabel('Energy')\n",
    "plt.title('VQE')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The result of this procedure is an estimate of the ground state energy of water. However, the convergence behavior is not perfect, more iterations would greatly improve the result, but would take a few minutes to run."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Parallel Parameter Shift Gradients"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "One way to accelerate VQE is to use an optimizer that accepts a gradient. This can drastically lower the number of VQE iterations required at the cost of computing the gradient on the quantum side of the algorithm.\n",
    "\n",
    "The parameter shift rule is a common technique to compute the gradient for parameterized circuits. It is obtained by computing two expectation values for each parameter corresponding to a small forward and backward shift in the ith parameter. These results are used to estimate  finite difference contribution to the gradient.\n",
    "\n",
    "![parametershift.png](./images/parametershift.png)\n",
    "\n",
    "This procedure can become cost prohibitive as the number of parameters becomes large."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Each of the expectation values needed to evaluate a parameter shift gradient can be computed independently.  The CUDA-Q `nvidia-mqpu` backend is designed for parallel computations across multiple simulated QPUs. The function below uses `cudaq.observe_asynch` to distribute all of the expectation values evaluations across as many GPUs that are available. First, try it with `num_qpus` set to 1."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.random.seed(42)\n",
    "x0 = np.random.normal(0, 1, parameter_count)\n",
    "\n",
    "cudaq.set_target(\"nvidia\", option=\"mqpu\")\n",
    "\n",
    "num_qpus = cudaq.get_target().num_qpus()\n",
    "\n",
    "epsilon = np.pi / 4\n",
    "\n",
    "\n",
    "def batched_gradient_function(kernel, parameters, hamiltonian, epsilon):\n",
    "\n",
    "    x = np.tile(parameters, (len(parameters), 1))\n",
    "\n",
    "    xplus = x + (np.eye(x.shape[0]) * epsilon)\n",
    "\n",
    "    xminus = x - (np.eye(x.shape[0]) * epsilon)\n",
    "\n",
    "    g_plus = []\n",
    "    g_minus = []\n",
    "    gradients = []\n",
    "\n",
    "    qpu_counter = 0  # Iterate over the number of GPU resources available\n",
    "    for i in range(x.shape[0]):\n",
    "\n",
    "        g_plus.append(\n",
    "            cudaq.observe_async(kernel,\n",
    "                                hamiltonian,\n",
    "                                qubit_count,\n",
    "                                electron_count,\n",
    "                                xplus[i],\n",
    "                                qpu_id=qpu_counter%num_qpus))\n",
    "        qpu_counter += 1\n",
    "\n",
    "        g_minus.append(\n",
    "            cudaq.observe_async(kernel,\n",
    "                                hamiltonian,\n",
    "                                qubit_count,\n",
    "                                electron_count,\n",
    "                                xminus[i],\n",
    "                                qpu_id=qpu_counter%num_qpus))\n",
    "        qpu_counter += 1\n",
    "\n",
    "    gradients = [\n",
    "        (g_plus[i].get().expectation() - g_minus[i].get().expectation()) /\n",
    "        (2 * epsilon) for i in range(len(g_minus))\n",
    "    ]\n",
    "\n",
    "    assert len(gradients) == len(\n",
    "        parameters) == x.shape[0] == xplus.shape[0] == xminus.shape[0]\n",
    "\n",
    "    return gradients"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The cost function needs to be slightly updated to make use of the gradient in the optimization procedure and allow for a gradient based optimizer like L-BFGS-B to be used."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "gradient = batched_gradient_function(kernel, x0, spin_ham, epsilon)\n",
    "\n",
    "exp_vals = []\n",
    "\n",
    "\n",
    "def objective_function(parameter_vector: list[float], \\\n",
    "                       gradient=gradient, hamiltonian=spin_ham, kernel=kernel):\n",
    "\n",
    "    get_result = lambda parameter_vector: cudaq.observe\\\n",
    "        (kernel, hamiltonian, qubit_count, electron_count, parameter_vector).expectation()\n",
    "\n",
    "    cost = get_result(parameter_vector)\n",
    "    exp_vals.append(cost)\n",
    "    gradient_vector = batched_gradient_function(kernel, parameter_vector,\n",
    "                                                spin_ham, epsilon)\n",
    "\n",
    "    return cost, gradient_vector"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Run the code below. Notice how the result is converged to a lower energy using only 10% of the steps as optimization above without a gradient. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "VQE-UCCSD energy=  -73.19436446992752\n",
      "Total elapsed time (s) =  72.83047100389376\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAj4AAAHHCAYAAAC/R1LgAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjMsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvZiW1igAAAAlwSFlzAAAPYQAAD2EBqD+naQAAR5pJREFUeJzt3Xd4VGXC/vF7JmXSC2kEEhISEqogREQQlCIgVlYsSEBAXAQBUdD9obiC7lpBbKuo74W4ithFFFeKIKgU6R1CDS1AGqSSOvP7IxANNYQkJzPz/VzXudaZc2ZyT959zb3PeeZ5TDabzSYAAAAnYDY6AAAAQG2h+AAAAKdB8QEAAE6D4gMAAJwGxQcAADgNig8AAHAaFB8AAOA0KD4AAMBpUHwAAIDToPgAAACnQfEBYBfuuOMOeXl5KScn54LXJCYmyt3dXRkZGZKkvLw8/etf/1Lr1q3l5eUlf39/denSRZ988onOt1uPyWS64DFixIga+2wAao+r0QEAoDISExP1ww8/aM6cOXrggQfOOZ+fn6+5c+fq5ptvVlBQkI4fP64ePXpox44d6t+/v0aPHq2CggJ98803euCBBzR//nx98sknMpsr/u+/nj17nvf94+Pja+yzAag9FB8AduGOO+6Qr6+vZs+efd5iMnfuXOXl5SkxMVGSNHjwYO3YsUNz5szRHXfcUX7do48+qieffFJTp07V1VdfrSeffLLC+8THx2vgwIE1+2EAGIZbXQDsgqenp+666y4tXrxYqamp55yfPXu2fH19dccdd2jVqlVasGCBhgwZUqH0nPHSSy8pLi5OL7/8sk6dOlUb8QHUERQfAHYjMTFRJSUl+vLLLys8n5mZqQULFuhvf/ubPD099cMPP0jSeUeGJMnV1VUDBgxQZmamVqxYUeFcQUGB0tPTzzmKiopq5kMBqFUUHwB2o3v37goPD9fs2bMrPP/VV1+puLi4/DbX9u3bJUlt2rS54HudOXfm2jNmzJihkJCQc45vv/22Oj8KAIMwxweA3XBxcVH//v31+uuvKzk5WdHR0ZLKbnOFhYWpR48eklT+zS9fX98LvteZc2d/S+zOO+/U6NGjz7n+qquuqo6PAMBgFB8AdiUxMVGvv/66Zs+eraefflqHDx/Wb7/9pkcffVQuLi6SKpaagICA877PmcITGhpa4fmIiAjddNNNNfcBABiKW10A7EpCQoKaNWumzz77TJL02WefyWazld/mkqQWLVpIkjZv3nzB9zlzLiYmpgbTAqhrKD4A7E5iYqK2bt2qzZs3a/bs2YqLi1P79u3Lz99+++2SpI8//vi8ry8tLS2/PXbDDTfUSmYAdQPFB4DdOTO68+yzz2rjxo0VRnsk6brrrlOvXr00c+ZMzZs375zXT5w4Ubt27dI//vEPubpyxx9wJibb+dZtB4A67vrrry//Kvru3bvVpEmTCuePHz+u7t27a+fOnRowYIC6dOmiwsJCffvtt1q6dKkGDhyojz/+WCaTqfw1JpPpgis3h4WFqWfPnjX7oQDUOIoPALv07rvvatSoUbr22mv1xx9/nPea3NxcTZs2TV9++aX27t2rgoICSdI///lPPf/88+dc/9cSdLYbb7xRS5curZbsAIxD8QHgNI4cOaJOnTqppKREK1euVKNGjYyOBKCWMccHgNNo2LCh5s+fr4KCAvXp00cnTpwwOhKAWsaIDwAAcBqM+AAAAKdB8QEAAE6D4gMAAJwGxQcAADgNliw9i9VqVUpKinx9fS+6pgcAAKg7bDabcnJy1KBBA5nNFx7XoficJSUlRZGRkUbHAAAAVXDo0CFFRERc8DzF5yy+vr6Syn5xfn5+BqcBAACVkZ2drcjIyPK/4xdC8TnLmdtbfn5+FB8AAOzMpaapMLkZAAA4DYoPAABwGhQfAADgNCg+AADAaVB8AACA06D4AAAAp0HxAQAAToPiAwAAnAbFBwAAOA2KDwAAcBoUHwAA4DQoPgAAwGlQfGpJUYlVv+9ONzoGAABOjeJTC3IKitVj2lIN+vAP7T6eY3QcAACcFsWnFvh6uKlFuJ9sNmnaol1GxwEAwGlRfGrJuJ5NZTJJP209pq1HsoyOAwCAU6L41JKm9X11R5sGkqTXFiYZnAYAAOdE8alFj98ULxezSb8kpWndgUyj4wAA4HQoPrUoOthb9yRESJKmLmCuDwAAtY3iU8vG9IiTu4tZK/dlaPkevt4OAEBtovjUsoYBnhrQoZEkacqCJNlsNoMTAQDgPCg+BnikW6w83MzaeOikluxMNToOAABOg+JjgFBfDw3uFC1Jmrpwl6xWRn0AAKgNFB+DjLghVr4WV+04mq3/bT1qdBwAAJwCxccggd7uGtalsaSy1ZxLSq0GJwIAwPFRfAw0rHNjBXi5aV9anr7bmGJ0HAAAHB7Fx0C+Hm4acWOsJOmNn3epqIRRHwAAahLFx2CDO0Yr2MeiwydO6cu1h4yOAwCAQ6P4GMzT3UWju5WN+ry9ZLcKiksNTgQAgOOi+NQB93dopIYBnjqeXahZqw4YHQcAAIdF8akDLK4uerRHE0nSu0v3Kq+wxOBEAAA4JrsqPj/++KM6dOggT09PBQYGqm/fvuXnNm3apPvvv1+RkZHy9PRU8+bN9eabbxoX9jLd1S5C0UFeyswr0szl+42OAwCAQ3I1OkBlffPNN/r73/+uF198Ud27d1dJSYm2bt1afn7dunUKDQ3VrFmzFBkZqRUrVmj48OFycXHR6NGjDUxeOW4uZj3eM15jP9+o93/dp0HXRcvfy83oWAAAOBSTzQ52ySwpKVF0dLSee+45DRs2rNKvGzVqlHbs2KElS5ZU+jXZ2dny9/dXVlaW/Pz8qhK3yqxWm/q8+ZuSjudodLcmeqJ301r9+QAA2KvK/v22i1td69ev15EjR2Q2m9W2bVuFh4erT58+FUZ8zicrK0v16tWrpZRXzmw26fGe8ZKkD5fvV3puocGJAABwLHZRfPbt2ydJmjx5sp555hnNmzdPgYGB6tq1qzIzM8/7mhUrVuiLL77Q8OHDL/rehYWFys7OrnAYqXfLMLWO8Fd+UammL91raBYAAByNocVnwoQJMplMFz127twpq7VsReOJEyeqX79+SkhI0MyZM2UymfTVV1+d875bt27VnXfeqUmTJqlXr14XzfDSSy/J39+//IiMjKyRz1pZJpNJ43uV3eL6ZNUBHcsqMDQPAACOxNDJzePHj9eQIUMuek1MTIyOHi3bvbxFixblz1ssFsXExOjgwYMVrt++fbt69Oih4cOH65lnnrlkhqeeekrjxo0rf5ydnW14+bkhLljtowO1JvmE3l6yWy/87SpD8wAA4CgMLT4hISEKCQm55HUJCQmyWCxKSkpS586dJUnFxcVKTk5WVFRU+XXbtm1T9+7dNXjwYL3wwguVymCxWGSxWKr2AWqIyWTSE72a6r4PVumLNYf08A2xahTkZXQsAADsnl3M8fHz89OIESM0adIkLVy4UElJSRo5cqQk6Z577pFUdnurW7du6tWrl8aNG6djx47p2LFjSktLMzJ6lXWICVKXuGCVWG16c/Fuo+MAAOAQ7KL4SNKUKVPUv39/DRo0SO3bt9eBAwe0ZMkSBQYGSpK+/vprpaWladasWQoPDy8/2rdvb3Dyqjsz12fOhsPak5pjcBoAAOyfXazjU5uMXMfnfP7+8Vot2n5ct14VrncS2xkdBwCAOsmh1vFxZuN7xctkkn7cclTbUrKMjgMAgF2j+NRxzer76bbWDSRJ0xbuMjgNAAD2jeJjBx6/KU4uZpMW70zV+oMnjI4DAIDdovjYgZgQH/Vr11CS9NrCJIPTAABgvyg+dmJM9zi5uZi0fE+GVuxNNzoOAAB2ieJjJyLreen+axtJkqYuSBJfxgMA4PJRfOzI6G5NZHE1a/3Bk1qaZJ8LMwIAYCSKjx0J9fPQ4E7RkqSpC5NktTLqAwDA5aD42JkRN8bKx+KqbSnZmr/tmNFxAACwKxQfO1PP210Pdm4sSZq2aJdKGfUBAKDSKD526KEujeXv6aY9qbmau/GI0XEAALAbFB875OfhpodvjJEkvfHzbhWXWg1OBACAfaD42KkhnaIV7OOug5n5+mrtYaPjAABgFyg+dsrL3VWPdG0iSXp7yW4VFJcanAgAgLqP4mPHBnRopHB/Dx3NKtCnfxw0Og4AAHUexceOebi56NEecZKk6Uv3KK+wxOBEAADUbRQfO3d3QoSigryUnlukj1YkGx0HAIA6jeJj59xczHrsprJRn/eX7VXWqWKDEwEAUHdRfBzAHW0aKi7UR9kFJZrx2z6j4wAAUGdRfByAi9mkcT3jJUkzft+vjNxCgxMBAFA3UXwcxM2t6qtVQz/lFZXqvWV7jY4DAECdRPFxECaTSeN7NZUkfbzygI5nFxicCACAuofi40C6xocoISpQhSVW/WfJHqPjAABQ51B8HIjJZNITp0d9Pl9zUIcy8w1OBABA3ULxcTAdY4PUuUmwikttemvxbqPjAABQp1B8HND4XmXf8Ppm/WHtTcs1OA0AAHUHxccBtW0UqJuah8pqk15ftMvoOAAA1BkUHwc1rmfZXJ95m49qx9Fsg9MAAFA3UHwcVIsGfrq1dbgk6bWFjPoAACBRfBza4zfFy2ySft5xXBsOnjA6DgAAhqP4OLAmoT66q12EJGkac30AAKD4OLqxPeLk5mLSb7vTtWpfhtFxAAAwFMXHwUXW89J97SMlSVMXJMlmsxmcCAAA41B8nMCY7nGyuJq19sAJLduVZnQcAAAMQ/FxAmF+Hhp0XZSksm94MeoDAHBWFB8nMbJrrLzdXbTlSJYWbDtmdBwAAAxB8XESQT4WPdi5saSyb3iVWhn1AQA4H4qPE3moS4z8PFy163iuftiUYnQcAABqHcXHifh7uunhG2MlSa//vEvFpVaDEwEAULsoPk5mSKdoBfu460BGvr5Zd9joOAAA1CqKj5PxtrhqZNcmkqS3Fu9WYUmpwYkAAKg9FB8nlNihker7eSglq0Cz/zhodBwAAGoNxccJebi5aEyPslGfd37Zq/yiEoMTAQBQOyg+TureayLVqJ6X0nML9d8VB4yOAwBAraD4OCk3F7PG9oiTJL23bK+yC4oNTgQAQM2j+Dixvm0bqkmoj7JOFWvGb/uNjgMAQI2zq+Lz448/qkOHDvL09FRgYKD69u173usyMjIUEREhk8mkkydP1mpGe+JiNmlcz3hJ0ozf9+tEXpHBiQAAqFl2U3y++eYbDRo0SEOHDtWmTZu0fPlyDRgw4LzXDhs2TK1bt67lhPbp5pb11SLcT7mFJXpv2V6j4wAAUKPsoviUlJRo7NixmjJlikaMGKH4+Hi1aNFC99577znXTp8+XSdPntQTTzxhQFL7Yzab9ETvslGf/65MVmp2gcGJAACoOXZRfNavX68jR47IbDarbdu2Cg8PV58+fbR169YK123fvl3PP/+8Pv74Y5nNlftohYWFys7OrnA4m25NQ9WuUYAKiq1655c9RscBAKDG2EXx2bdvnyRp8uTJeuaZZzRv3jwFBgaqa9euyszMlFRWYO6//35NmTJFjRo1qvR7v/TSS/L39y8/IiMja+Qz1GUmk0lP9GoqSZq9+qAOn8g3OBEAADXD0OIzYcIEmUymix47d+6U1Vq2mebEiRPVr18/JSQkaObMmTKZTPrqq68kSU899ZSaN2+ugQMHXlaGp556SllZWeXHoUOHqv1z2oNOTYLVKTZIxaU2vb2YUR8AgGNyNfKHjx8/XkOGDLnoNTExMTp69KgkqUWLFuXPWywWxcTE6ODBsi0XlixZoi1btujrr7+WJNlsNklScHCwJk6cqOeee+6872+xWGSxWK70oziE8b2aasX0Ffp6/WGN6BqrxsHeRkcCAKBaGVp8QkJCFBIScsnrEhISZLFYlJSUpM6dO0uSiouLlZycrKioKEll3/o6depU+WvWrFmjBx98UL/99ptiY2Nr5gM4mISoQHVvFqolO1P1+qJdeuv+tkZHAgCgWhlafCrLz89PI0aM0KRJkxQZGamoqChNmTJFknTPPfdI0jnlJj09XZLUvHlzBQQE1Gpeeza+V7yW7EzVD5tT9Ei3WDWr72d0JAAAqo1dTG6WpClTpqh///4aNGiQ2rdvrwMHDmjJkiUKDAw0OppDadnAX7deFS6bTZq2cJfRcQAAqFYm25nJMJAkZWdny9/fX1lZWfLzc87Rjj2pOer1+q+y2qS5o65Xm8gAoyMBAHBRlf37bTcjPqg9TUJ91bdtQ0nSa4sY9QEAOA6KD87rsR7xcjWb9OuuNK3en2l0HAAAqgXFB+fVKMhL97YvW8xx6oIkcUcUAOAIKD64oDHdm8jd1azVyZn6dXe60XEAALhiFB9cULi/pwZdV7ZO0msLGfUBANg/ig8uamTXWHm5u2jz4Swt3H7c6DgAAFwRig8uKtjHoqHXR0sqW9en1MqoDwDAflF8cEnDu8TK18NVScdzNG9zitFxAACoMooPLsnfy03Du8RIkt74ebdKSq0GJwIAoGooPqiUoZ0bq563u/an5+mb9YeNjgMAQJVQfFApPhZXPdK1bCPYtxbvUWFJqcGJAAC4fBQfVNrA66IU5mfRkZOn9PnqQ0bHAQDgslF8UGkebi4a3T1OkvSfX/boVBGjPgAA+0LxwWW575pIRQR6Ki2nUB+vTDY6DgAAl4Xig8vi7mrW2B5loz7Tl+1VTkGxwYkAAKg8ig8u29/aNlRMiLdO5hdrxu/7jY4DAEClUXxw2VxdzBrXM16SNOO3/TqRV2RwIgAAKofigyq5pVW4mof7KaewRO//us/oOAAAVArFB1ViNps0/vSoz0cr9is1p8DgRAAAXBrFB1XWo3moro4MUEGxVe/+stfoOAAAXBLFB1VmMpn0RK+mkqTZfxzUkZOnDE4EAMDFUXxwRa5vEqTrYuqpqNSqtxfvNjoOAAAXRfHBFfnrqM9X6w4rOT3P4EQAAFwYxQdX7JroeuraNESlVpve+HmX0XEAALggig+qxZlRn7mbUpR0LMfgNAAAnB/FB9WiVUN/9WlVXzab9PoiRn0AAHUTxQfV5vGe8TKZpPnbjmnL4Syj4wAAcA6KD6pNfJiv+l7dUJI0dWGSwWkAADgXxQfV6rGb4uRiNmnZrjStSc40Og4AABVQfFCtooK8de81EZKkqQuSZLPZDE4EAMCfKD6odmO6x8ndxaw/9mfq9z3pRscBAKAcxQfVrkGApxKvayRJmrpwF6M+AIA6g+KDGvFI1ybydHPRpkMn9fOOVKPjAAAgieKDGhLia9GQ66MlSa8tTJLVyqgPAMB4FB/UmIdviJGvxVU7j+Xoxy1HjY4DAADFBzUnwMtdD3WJkVS2mnNJqdXgRAAAZ0fxQY16sHO0Ar3ctC89T99uOGJ0HACAk6P4oEb5erhpZNdYSdKbP+9WUQmjPgAA41B8UOMGXRetEF+Ljpw8pS/WHDQ6DgDAiVF8UOM83V00pnsTSdLbS/boVFGpwYkAAM6K4oNacV/7SDUM8FRqTqFmrTpgdBwAgJOi+KBWWFxdNLZHnCRp+rK9yi0sMTgRAMAZUXxQa+5q11Axwd7KzCvSh7/vNzoOAMAJUXxQa1xdzHqsZ7wk6f9+3aeT+UUGJwIAOBuKD2rVbVeFq1l9X+UUluiDX/cZHQcA4GTsqvj8+OOP6tChgzw9PRUYGKi+ffuec81HH32k1q1by8PDQ6GhoRo1alTtB8UFmc0mjTs96jNzebLScgoNTgQAcCauRgeorG+++UZ///vf9eKLL6p79+4qKSnR1q1bK1wzbdo0vfbaa5oyZYo6dOigvLw8JScnGxMYF9SzRZjaRPhr0+EsTV+6V8/e3sLoSAAAJ2Gy2Wx1ftvskpISRUdH67nnntOwYcPOe82JEyfUsGFD/fDDD+rRo0eVf1Z2drb8/f2VlZUlPz+/Kr8PLu7XXWl64MPVcnc1a9mTXRXu72l0JACAHavs32+7uNW1fv16HTlyRGazWW3btlV4eLj69OlTYcRn0aJFslqtOnLkiJo3b66IiAjde++9OnTo0EXfu7CwUNnZ2RUO1LwuccG6tnE9FZVY9dbiPUbHAQA4CbsoPvv2lU2CnTx5sp555hnNmzdPgYGB6tq1qzIzM8uvsVqtevHFF/XGG2/o66+/VmZmpnr27Kmiogt/e+ill16Sv79/+REZGVkrn8nZmUwmPdGrqSTpq7WHdCAjz+BEAABnYGjxmTBhgkwm00WPnTt3ymot29hy4sSJ6tevnxISEjRz5kyZTCZ99dVXkiSr1ari4mK99dZb6t27t6677jp99tln2r17t3755ZcLZnjqqaeUlZVVflxqhAjV59rG9XRDfIhKrDa9+fNuo+MAAJyAoZObx48fryFDhlz0mpiYGB09elSS1KLFn5NgLRaLYmJidPBg2aaX4eHh51wTEhKi4ODg8mvOx2KxyGKxVPUj4Ao90Stev+5K05yNRzSya6ziwnyNjgQAcGCGFp+QkBCFhIRc8rqEhARZLBYlJSWpc+fOkqTi4mIlJycrKipKknT99ddLkpKSkhQRESFJyszMVHp6evk1qHtaRwSod8swLdh2XK//vEvvJiYYHQkA4MDsYo6Pn5+fRowYoUmTJmnhwoVKSkrSyJEjJUn33HOPJCk+Pl533nmnxo4dqxUrVmjr1q0aPHiwmjVrpm7duhkZH5cwrmdTmUzS/7Yc09YjWUbHAQA4MLsoPpI0ZcoU9e/fX4MGDVL79u114MABLVmyRIGBgeXXfPzxx+rQoYNuvfVW3XjjjXJzc9P8+fPl5uZmYHJcStP6vrqjTQNJ0msLkwxOAwBwZHaxjk9tYh0fY+xPz9NN05ap1GrTNyM7KiGqntGRAAB2xKHW8YHjaxzsrbvblc3Nmrpgl8FpAACOiuKDOuPRm+Lk7mLWyn0ZWr4n3eg4AAAHRPFBndEwwFMDOjSSJE1dmCTuwgIAqhvFB3XKI91i5eFm1oaDJ7VkZ6rRcQAADobigzol1NdDgztFS5KmLtwlq5VRHwBA9alS8cnLY18l1JwRN8TKx+KqHUez9dPWY0bHAQA4kCoVn7CwMD344IP6/fffqzsPoEBvdw3r3FiSNG1RkkoZ9QEAVJMqFZ9Zs2YpMzNT3bt3V3x8vF5++WWlpKRUdzY4sYe6NFaAl5v2puVpzoYjRscBADiIKhWfvn376rvvvtORI0c0YsQIzZ49W1FRUbrtttv07bffqqSkpLpzwsn4erhpxI2xkqQ3F+9SUYnV4EQAAEdwRZObQ0JCNG7cOG3evFnTpk3Tzz//rLvvvlsNGjTQs88+q/z8/OrKCSf0QMcoBftYdCjzlL5ce8joOAAAB3BFxef48eN69dVX1aJFC02YMEF33323Fi9erNdee03ffvut+vbtW00x4Yy83F01ulvZqM/bS3aroLjU4EQAAHvnWpUXffvtt5o5c6YWLFigFi1a6JFHHtHAgQMVEBBQfk2nTp3UvHnz6soJJ3V/h0b64Nd9Sskq0KxVB/RQlxijIwEA7FiVRnyGDh2qBg0aaPny5dq4caNGjx5dofRIUoMGDTRx4sTqyAgnZnF10aM94iRJ05fuVV4h88cAAFVXpd3Z8/Pz5eXlVRN5DMfu7HVPcalVPactU3JGvp7oFa/R3eOMjgQAqGNqdHf2kpISZWdnn3Pk5OSoqKioyqGB83FzMevxnvGSpPd/3aes/GKDEwEA7FWVik9AQIACAwPPOQICAuTp6amoqChNmjRJVitfQUb1uK11A8WH+SinoET/99s+o+MAAOxUlYrPRx99pAYNGujpp5/Wd999p++++05PP/20GjZsqOnTp2v48OF666239PLLL1d3XjgpF7NJ43o2lSR9uHy/0nMLDU4EALBHVfpW13//+1+99tpruvfee8ufu/3223XVVVfp/fff1+LFi9WoUSO98MILevrpp6stLJxb75Zhuqqhv7YcydJ7S/fqmdtaGB0JAGBnqjTis2LFCrVt2/ac59u2bauVK1dKkjp37qyDBw9eWTrgL0wmk8b3Kpvr8/GqAzqWVWBwIgCAvalS8YmMjNSMGTPOeX7GjBmKjIyUJGVkZCgwMPDK0gFnuTE+RO2jA1VUYtXbS3YbHQcAYGeqdKtr6tSpuueee/TTTz+pffv2kqS1a9dq586d+vrrryVJa9as0X333Vd9SQGVjfo80aup7vtglb5Yc0gjboxVZD3HXFoBAFD9qrSOjyQlJyfr/fffV1JSkiSpadOmevjhhxUdHV2d+Wod6/jYh0Ez/tBvu9PVr12EXru3jdFxAAAGq+zf78suPsXFxbr55pv13nvvKS7O8RaSo/jYh42HTqrvO8tlNkkLH79BTUJ9jY4EADBQjS1g6Obmps2bN19ROOBKXR0ZoJ4twmS1Sa//zFwfAEDlVGly88CBA887uRmoTeN6xstkkn7cfFTbUrKMjgMAsANVmtxcUlKiDz/8UD///LMSEhLk7e1d4fy0adOqJRxwMc3D/XRb6wb6YVOKpi3cpRlD2hsdCQBQx1Wp+GzdulXt2rWTJO3atavCOZPJdOWpgEp6/KY4/W/LUS3emar1B0+oXSOWUAAAXFiVis8vv/xS3TmAKokJ8VG/dg315drDem1hkj596DqjIwEA6rAqzfE5Y8+ePVqwYIFOnTolSariN+OBKzKme5zcXExavidDK/amGx0HAFCHVan4ZGRkqEePHoqPj9ctt9yio0ePSpKGDRum8ePHV2tA4FIi63np/msbSZJeW7iLAg4AuKAqFZ/HH39cbm5uOnjwoLy8/lw197777tP8+fOrLRxQWaO7NZHF1ax1B05oaVKa0XEAAHVUlYrPwoUL9corrygiIqLC83FxcTpw4EC1BAMuR6ifhwZ3ipYkTV2YJKuVUR8AwLmqVHzy8vIqjPSckZmZKYvFcsWhgKoYcWOsfCyu2paSra/XHzY6DgCgDqpS8enSpYs+/vjj8scmk0lWq1WvvvqqunXrVm3hgMtRz9tdI7vGSpKenbtV21OyDU4EAKhrqvR19ldffVU9evTQ2rVrVVRUpH/84x/atm2bMjMztXz58urOCFTaiBtjtXp/ppbtStOIWev0w+jO8vdyMzoWAKCOqNKIT6tWrbRr1y517txZd955p/Ly8nTXXXdpw4YNio2Nre6MQKW5mE16s//Vigj01MHMfD32xQbm+wAAyl327uyOjt3ZHcPWI1nqN32FCkusGtsjTo/3jDc6EgCgBlX273eVbnVJ0smTJ7V69WqlpqbKarVWOPfAAw9U9W2BatGqob9e+NtVeuKrTXpz8W61ifRX92ZhRscCABisSiM+P/zwgxITE5Wbmys/P78K+3OZTCZlZmZWa8jaxIiPY3nmuy2ateqg/Dxc9cOYzooK8r70iwAAdqeyf7+rNMdn/PjxevDBB5Wbm6uTJ0/qxIkT5Yc9lx44nmdva6m2jQKUXVCiEbPW61RRqdGRAAAGqlLxOXLkiB599NHzruUD1CXurmZNT0xQsI+7dhzN1tNztrClBQA4sSoVn969e2vt2rXVnQWoEfX9PfSfAe3kYjZpzoYj+nglq4sDgLOq0uTmW2+9VU8++aS2b9+uq666Sm5uFddJueOOO6olHFBdrosJ0lN9munfP+7Qv+ZtV8sGfromup7RsQAAtaxKk5vN5gsPFJlMJpWW2u88CiY3Oy6bzaYxn23QvM1HFepr0bxHOyvU18PoWACAalCjk5utVusFD3suPXBsJpNJr/RrrbhQH6XmFGr0pxtUXGq99AsBAA7jsorPLbfcoqysrPLHL7/8sk6ePFn+OCMjQy1atKi2cGf78ccf1aFDB3l6eiowMFB9+/atcH7NmjXq0aOHAgICFBgYqN69e2vTpk01lgf2x9viqvcGJcjH4qrVyZl68X87jI4EAKhFl1V8FixYoMLCwvLHL774YoWvr5eUlCgpKan60v3FN998o0GDBmno0KHatGmTli9frgEDBpSfz83N1c0336xGjRrpjz/+0O+//y5fX1/17t1bxcXFNZIJ9ik2xEev3dtGkjRzebLmbjxicCIAQG25rMnNZ08Hqq2vBZeUlGjs2LGaMmWKhg0bVv78X0eXdu7cqczMTD3//POKjIyUJE2aNEmtW7fWgQMH1KRJk1rJCvvQu2V9PdI1Vu8u3asJ32xR0/q+alafOV0A4OiqNMentq1fv15HjhyR2WxW27ZtFR4erj59+mjr1q3l1zRt2lRBQUGaMWOGioqKdOrUKc2YMUPNmzdXdHT0Bd+7sLBQ2dnZFQ44h/G9mqpLXLBOFZdqxCfrlHWKkUEAcHSXVXxMJlOF7SnOPFfT9u3bJ0maPHmynnnmGc2bN0+BgYHq2rVr+a02X19fLV26VLNmzZKnp6d8fHw0f/58/fTTT3J1vfDA1ksvvSR/f//y48xoERxf2U7ubdUwwFPJGfka/+VGdnIHAAd3WV9nN5vN6tOnjywWi6SyPbu6d+8ub++y/Y8KCws1f/78Sn+za8KECXrllVcues2OHTu0fv16JSYm6v3339fw4cPLf1ZERIT+/e9/6+GHH9apU6fUtWtXNWvWTKNHj1ZpaammTp2qnTt3as2aNfL09Dzv+xcWFlaYt5Sdna3IyEi+zu5EthzOUr/3VqioxKrxPeM1pkec0ZEAAJepRnZnHzx4cIXHAwcOPOeay9mZffz48RoyZMhFr4mJidHRo0clVZzTY7FYFBMTo4MHD0qSZs+ereTkZK1cubJ8naHZs2crMDBQc+fOVf/+/c/7/haLpbzIwTldFeGvf9/ZSv/4ZrOm/bxLV0X4q2vTUKNjAQBqwGUVn5kzZ1brDw8JCVFISMglr0tISJDFYlFSUpI6d+4sSSouLlZycrKioqIkSfn5+TKbzRVuvZ15bLWyVgsu7t72kdpw6KQ+W31QYz/fqHljOiuyHnvRAYCjsYvJzX5+fhoxYoQmTZqkhQsXKikpSSNHjpQk3XPPPZKknj176sSJExo1apR27Nihbdu2aejQoXJ1dVW3bt2MjA87MfmOFmoTGaCsU8V6+JN1KihmMU4AcDR2UXwkacqUKerfv78GDRqk9u3b68CBA1qyZIkCAwMlSc2aNdMPP/ygzZs3q2PHjurSpYtSUlI0f/58hYeHG5we9sDi6qLpie1Uz9td249ma+KcrezkDgAOpkp7dTky9urCij3pGjjjD1lt0r/6ttKg66KMjgQAuIQa3asLcGSdmgTrHzc3kyQ9/8M2rTtwwuBEAIDqQvEBzuPhG2LUp1V9FZfa9Min65SWU3jpFwEA6jyKD3AeJpNJU+5po9gQbx3PLtSYz9arhJ3cAcDuUXyAC/CxuOr9QdfI291Fq/Zl6pX5O42OBAC4QhQf4CKahPpo6j1lO7n/32/7NW9zisGJAABXguIDXEKfq8L18I0xkqR/fL1Zu47nGJwIAFBVFB+gEp7s1VSdYoOUX1S2k3t2ATu5A4A9ovgAleDqYtbb97dVuL+H9qXn6YkvN7GTOwDYIYoPUElBPhZNH5ggdxezFm4/runL9hodCQBwmSg+wGW4OjJAk+9oKUl6bWGSftudZnAiAMDloPgAl+n+ayN17zURstqkRz/boMMn8o2OBACoJIoPcJlMJpOev7OVrmrorxP5xRo5az07uQOAnaD4AFXg4eai6QPbKdDLTVuOZGnS3G1GRwIAVALFB6iiiEAvvXV/W5lN0hdrD+mz1QeNjgQAuASKD3AFusSFaHyvppKkSXO3aeOhk8YGAgBcFMUHuEKPdI1VrxZhKiq16pFZ65SRy07uAFBXUXyAK2QymTT13jaKCfZWSlaBxny2gZ3cAaCOovgA1cDPw03vDUqQl7uLVuzN0JSFSUZHAgCcB8UHqCbxYb569e7WkqT3l+3TT1uOGpwIAHA2ig9QjW5r3UAPdW4sSXriq03ak8pO7gBQl1B8gGo2oU8zdWhcT3lFpXr4k3XKLSwxOhIA4DSKD1DNXF3M+s+Adqrv56G9aXl68qtNstnYyR0A6gKKD1ADQnwtendgO7m5mPTT1mP64Nd9RkcCAIjiA9SYdo0C9eztZTu5vzJ/p1bsSTc4EQCA4gPUoIEdGqlfu7Kd3Ed/tkEpJ08ZHQkAnBrFB6hBJpNJL/ytlVqE+ykzr0gjP12vwhJ2cgcAo1B8gBrm4eai9wclyN/TTZsOndTk77cbHQkAnBbFB6gFkfW89Gb/q2UySZ+tPqgv1xwyOhIAOCWKD1BLujYN1eM3xUuSnpm7VZsPnzQ2EAA4IYoPUItGd2uim5qHqqjEqpGz1iszr8joSADgVCg+QC0ym0167d6rFR3kpSMnT2ns5xtUamVxQwCoLRQfoJb5e5bt5O7p5qLfdqdr2iJ2cgeA2kLxAQzQrL6fXu53lSTpnV/2auG2YwYnAgDnQPEBDHLn1Q019PpoSdL4LzdpX1qusYEAwAlQfAADPX1Lc7WPDlROYYke/mSd8tjJHQBqFMUHMJCbi1nvDGinEF+Ldqfm6h/fbGYndwCoQRQfwGChfh6anthOrmaTftx8VDN+3290JABwWBQfoA64Jrqenrm1uSTppZ92auXeDIMTAYBjovgAdcTgTtHqe3UDlVptGvPZeh3NYid3AKhuFB+gjjCZTHrprtZqVt9X6blFeuTT9SoqsRodCwAcCsUHqEM83ct2cvfzcNWGgyf1r3ns5A4A1YniA9QxUUHeeqP/1ZKkT1Yd0DfrDhsbCAAcCMUHqIO6NwvT2B5xkqSn52zR1iNZBicCAMdA8QHqqLE94tStaYgKS6wa+ek6ncxnJ3cAuFIUH6COMptNev2+qxVZz1OHMk9p7Ocb2ckdAK6QXRSfpUuXymQynfdYs2ZN+XWbN29Wly5d5OHhocjISL366qsGpgauXICXu94bmCCLq1nLdqXpzZ93GR0JAOyaXRSfTp066ejRoxWOhx56SI0bN9Y111wjScrOzlavXr0UFRWldevWacqUKZo8ebI++OADg9MDV6ZlA3+9dFfZTu5vLdmjn7cfNzgRANgvuyg+7u7uql+/fvkRFBSkuXPnaujQoTKZTJKkTz/9VEVFRfrwww/VsmVL9e/fX48++qimTZtmcHrgyt3VLkIPdIySJD3+5UYlp+cZnAgA7JNdFJ+zff/998rIyNDQoUPLn1u5cqVuuOEGubu7lz/Xu3dvJSUl6cSJExd8r8LCQmVnZ1c4gLromVtbKCEqUDkFJRoxa53yi9jJHQAul10WnxkzZqh3796KiIgof+7YsWMKCwurcN2Zx8eOHbvge7300kvy9/cvPyIjI2smNHCF3F3NejexnYJ9LNp5LEdPfbuFndwB4DIZWnwmTJhwwUnLZ46dO3dWeM3hw4e1YMECDRs2rFoyPPXUU8rKyio/Dh06VC3vC9SEMD8PvTOgrVzMJs3dmKKPViQbHQkA7IqrkT98/PjxGjJkyEWviYmJqfB45syZCgoK0h133FHh+fr16+v48YqTPs88rl+//gXf32KxyGKxXEZqwFgdYoL09C3N9a952/XCjzvUqqG/2kfXMzoWANgFQ4tPSEiIQkJCKn29zWbTzJkz9cADD8jNza3CuY4dO2rixIkqLi4uP7do0SI1bdpUgYGB1ZobMNqD10dr46GT+mFTih75dL1+HNNZoX4eRscCgDrPrub4LFmyRPv379dDDz10zrkBAwbI3d1dw4YN07Zt2/TFF1/ozTff1Lhx4wxICtQsk8mkV/pdpaZhvkrLKWQndwCoJLsqPjNmzFCnTp3UrFmzc875+/tr4cKF2r9/vxISEjR+/Hg9++yzGj58uAFJgZrn5e6q9wYlyNfiqrUHTujF/+0wOhIA1HkmG18LqSA7O1v+/v7KysqSn5+f0XGAS1q0/bj+/vFaSdLr97XR39pGXOIVAOB4Kvv3265GfACcq2eLMI3u1kSS9NS3W7Q9hbWoAOBCKD6AA3i8Z7xuiA9RQbFVI2atU1Z+sdGRAKBOovgADsDFbNKb912tiEBPHczM1+NfbpSVndwB4BwUH8BBBHr/uZP7kp2penvJHqMjAUCdQ/EBHEirhv76d99WkqQ3Fu/SL0mpBicCgLqF4gM4mHuuiVRih0ay2aSxn23QwYx8oyMBQJ1B8QEc0LO3t9DVkQHKLijRw7PW6VRRqdGRAKBOoPgADsji6qLpA9spyNtdO45ma+IcdnIHAIniAziscH9PvT2grcwm6dsNR/TJqgNGRwIAw1F8AAfWKTZYE/qUbfHy/A/bte5ApsGJAMBYFB/Awf29S4xuvSpcJVabRs5ar9ScAqMjAYBhKD6AgzOZTHrl7tZqEuqj1JxCjZ69QcWl7OQOwDlRfAAn4GNx1fuDEuRjcdXq/Zl6+aedRkcCAENQfAAnERvio6n3tJEkzfh9v77flGJwIgCofRQfwInc3Kq+RnaNlST9v683K+lYjsGJAKB2UXwAJ/NEr6bq3CRYp4pLNWLWOmUXsJM7AOdB8QGcjIvZpDf7X60G/h7an56ncV9sYid3AE6D4gM4oSAfi6YPTJC7i1k/7ziud5eykzsA50DxAZxUm8gAPX9nS0nSa4t2admuNIMTAUDNo/gATqz/tY3Uv31k2U7un2/QoUx2cgfg2Cg+gJObfEdLtY7w18n8Yo2YtU4FxezkDsBxUXwAJ+fh5qLpAxNUz9td21Ky9cx3W9nJHYDDovgAUMMAT719f9lO7l+vO6zZqw8aHQkAagTFB4Ak6fomwXqyd9lO7pO/36YNB08YnAgAqh/FB0C5ETfGqHfLMBWXlu3knp5baHQkAKhWFB8A5Uwmk6be00YxId46ll2g0bPXq4Sd3AE4EIoPgAp8Pdz0/sAEebm7aNW+TL26IMnoSABQbSg+AM4RF+arKXeX7eT+wa/79OPmowYnAoDqQfEBcF63tg7X8BtiJElPfr1Ju4+zkzsA+0fxAXBB/+jdVB1jgpRfVKqHP1mnHHZyB2DnKD4ALsjVxay3B7RVuL+H9qXn6YmvNrG4IQC7RvEBcFHBPha9m9hO7i5mLdh2XO8t22d0JACoMooPgEtq2yhQk+5oIUmasmCnlu9JNzgRAFQNxQdApQy4tpHuSYiQ1SaN+WyDjpw8ZXQkALhsFB8AlWIymfSvvq3UqqGfMvOKNJKd3AHYIZONmYoVZGdny9/fX1lZWfLz8zM6DlDnHMrM1+3/+V0n84vl6eaia6ID1Sk2WJ1ig9Sqob9czCajIwJwQpX9+03xOQvFB7i0FXvT9djnG5WaU3EvL18PV3VoXE8dTxehpmG+MlOEANQCik8VUXyAyrFabdqVmqOVezO0Ym+GVu3LUE5BSYVr6nm7q2NMkK6LDVKn2CDFBHvLZKIIAah+FJ8qovgAVVNqtWlbSlZ5EVqTnKn8oopzgML8LOoUG6yOMUHqGBukyHpeBqUF4GgoPlVE8QGqR1GJVZsPn9SKvRlauTdD6w6eUFFJxZ3eI+t5qlNMsDrGlhWhMD8Pg9ICsHcUnyqi+AA1o6C4VOsPnCgrQvsytOnQSZVYK/7rJzbEu2xEKDZI18UEqZ63u0FpAdgbik8VUXyA2pFbWKI1yZlaeXpEaGtKls7+t1HzcD91ig1Sx5ggXRtTT34ebsaEBVDnUXyqiOIDGCMrv1ir9meUF6Gks3aDN5ukqyIC1On0ROlrourJ093FoLQA6hqKTxVRfIC6IS2nUKv2ZZyeI5Su5Iz8CufdXExqGxlYPj+obaMAWVwpQoCzovhUEcUHqJtSTp4q/8bYyr3pSskqqHDew82sa6LqqePpEaGrGvrL1YXF6QFn4VDFZ+nSperWrdt5z61evVrt27fX0qVL9frrr2v16tXKzs5WXFycnnzySSUmJl7Wz6L4AHWfzWbTwcx8rfhLEUrPLapwjY/FVdc2rlc2Ryg2SM3r+7GYIuDAHKr4FBUVKTMzs8Jz//znP7V48WLt3btXJpNJL774ok6dOqU+ffooLCxM8+bN07hx4zR37lzddtttlf5ZFB/A/thsNu1OzT09IpSuVfsylXWquMI1AV5uuq5xkDo1KRsRig3xYTFFwIE4VPE5W3FxsRo2bKgxY8bon//85wWvu/XWWxUWFqYPP/yw0u9N8QHsX6nVph1Hs8uL0Or9mco7azHFEF+LOsYEnZ4sHazIep4UIcCOVfbvt2stZqo233//vTIyMjR06NCLXpeVlaXmzZtf9JrCwkIVFv6531B2dna1ZARgHBezSa0a+qtVQ3/9/YYYFZdatflw1unJ0ulam3xCaTmF+n5Tir7flCJJahjgWT4/qGNskML9PQ3+FABqgl2O+Nxyyy2SpP/9738XvObLL7/UoEGDtH79erVs2fKC102ePFnPPffcOc8z4gM4roLiUm04eFIr95XND9pw8NzFFGOCvcv3GLsuJkjBPhaD0gKoDLu41TVhwgS98sorF71mx44datasWfnjw4cPKyoqSl9++aX69et33tf88ssvuu222zR9+nQ98MADF33/8434REZGUnwAJ5JXWKK1B06cXkMoXVuOZOmsHqRm9X113elbYx1iguTvyWKKQF1iF8UnLS1NGRkZF70mJiZG7u5/Llv/r3/9S2+//baOHDkiN7dz/8WzbNky3XrrrZo2bZqGDx9+2ZmY4wMg61SxVu/P1Iq96Vq5N0M7j527mGKrhv7lm622j64nb4tdzhwAHIZdFJ/LZbPZFBsbq7vuuktTp0495/zSpUt122236ZVXXtGoUaOq9DMoPgDOlpFbqFX7ThehfRnal5ZX4byr2aSrIwPKF1Ns1yhQHm4spgjUJocsPosXL9ZNN910zu0v6c/bW2PHjtWjjz5a/ry7u7vq1atX6Z9B8QFwKceyCrRyX7pW7ClbR+jIyVMVzru7mnVNVGD5ROnWEQFyYzFFoEY5ZPEZMGCADhw4oOXLl59zbsiQIfrvf/97zvM33nijli5dWumfQfEBcLkOZeZrxd7004spZig1p7DCeS93lz8XU4wJVosGfnJhMUWgWjlk8akNFB8AV8Jms2lvWp5WnilC+zJ0Mr/iYop+Hq7lE6U7xgYrPozFFIErRfGpIooPgOpktdq081hO+UTpP/ZnKrewpMI1wT7up4tQsDrGBik6yIsiBFwmik8VUXwA1KSSUqu2pmSXF6E1yZkqKLZWuCbc3+P0YoplRahhAIspApdC8akiig+A2lRYUqqNpxdTXLE3QxsOnlBxacV/LXu5uyjU16JQXw+F+FkU6mtRyOnHob4WhfpZFOJjUaCXOxuxwmlRfKqI4gPASKeKSrXuwInyydKbD588ZzHFC3FzMSnY50wx8jhdjsqK0ZmSFHL64FtmcDQUnyqi+ACoS04VlepYdoFSswuUmlOotJxCpeYUKjWnoOyfswuVlluozLyiy3rfet7uFYrQ2aNHoX5lj1mYEfbCoTcpBQBn4enuosbB3moc7H3R64pKrErPPV2KsguUlltWisrKUkF5YUrLKVSJ1abMvCJl5hWdsyr12bzdXcqL0ZnbbKHnGU0K9HJjQjbsAsUHAByAu6tZDQI81eASE6GtVptO5BedM3p0ZuQoLfv045xC5ReVKq+oVHkZ+UrOyL/o+559m+3PkaOKo0nBPtxmg7EoPgDgRMxmk4J8LArysah5+MWvzS0sOX07reB0QfrzFtuZ22ypOQU6kV+s4lKbjmYV6GhWgaSsi77vX2+znW/06MxjbrOhJvDfKgDAeflYXOVjcb3s22x/HU1KOz16lJpdqPTcy7/NFurn8Zd5SGfNRfLlNhsuH8UHAHBFqnKb7c9ydPo2218nbP/lNtv+9DztT8+76Pu6uZgU4nNmorbH6ZGjc0eTuM0GieIDAKgll3ubLTW74C/zkE4Xo9Nzkc6+zZaSVaCUKtxm++tcpMhALzUJ9eEWm4Pj/7oAgDrHx+IqnxAfxYT4XPS6ohJr2aTss+YipV3BbbaIQE/Fh/kqLsxH8aG+ig/zVZNQH3m6u1T3x4QBWMfnLKzjAwCO5+zbbH/9yn9aTqGOZxcoOSNP6bnnXw/JZJIiA70UH+ajuDBfNT1djGJDfOThRiGqC1jAsIooPgDgvDLzirTreI52H8/RruO5Zf+cmnvBBSLNJikqyFtxoT5/jhKF+SomxFsWVwpRbaL4VBHFBwBwtvTcwtOFKPfP/0zN0cn84vNe72I2KSrI6/StsrJRovgwXzUO9pa7KxOsawLFp4ooPgCAyrDZbErLLSwvQ7uO52r38RwlHc9RTkHJeV/jajYpOti7rAydnj8UH+aj6GBvvnF2hSg+VUTxAQBcCZvNpuPZhafL0J+jQ7uP5yq38PyFyM3FpMbB3qeL0J+jRFH1vORKIaoUik8VUXwAADXBZitb3fqvt8x2peZqz/Ec5RWVnvc17i5mxYR4VyhD8WG+alTPSy5mFm38K4pPFVF8AAC1yWq1KSXrVMVbZqdHiE4Vn78QWVzNig3xqVCG4sN8FBnoJbOTFiKKTxVRfAAAdYHVatPhE6dOjwz9OUq0JzVXhSXW877Gw82sJqFl6w/Fhfmqaf2yuUQNAzwdvhBRfKqI4gMAqMtKrTYdyswv/6r9mVGivWm5KrpAIfJyd1Fc6JnRoT9HiRr4ezjMPmcUnyqi+AAA7FFJqVUHM/PLv122K7XsP/em5aq49Px/6n0srmUjRGFn1iEqK0b1/eyvEFF8qojiAwBwJCWlViVn5Ff8ltnxHO1Pz1OJ9fwVwNfD9S+LMvqWF6NQX0udLUQUnyqi+AAAnEFRiVXJGXkV1iDadTxHyRn5Kr1AIfL3dPvzVtlfilGwj7vhhYjiU0UUHwCAMyssKdX+9LwKZWj38VwlZ+TpAn1IgV5uFUaG4k6vWB3kY6m13BSfKqL4AABwroLiUu1Ly9Pu1JzyUaJdx3N0MDNfF2oSQd7u5fuX/XWUKNDbvdrzUXyqiOIDAEDlnSoq1d60itt27ErN0aHMUxd8zTcjOykhKrBac1T277drtf5UAADgVDzdXdSqob9aNfSv8Hx+UYn2pOZWuGW263iujpw8pSYhPgalpfgAAIAa4OXuqtYRAWodEVDh+bzCEnlbjKsf7HwGAABqjZGlR6L4AAAAJ0LxAQAAToPiAwAAnAbFBwAAOA2KDwAAcBoUHwAA4DQoPgAAwGlQfAAAgNOg+AAAAKdB8QEAAE6D4gMAAJwGxQcAADgNig8AAHAaxm6RWgfZbDZJUnZ2tsFJAABAZZ35u33m7/iFUHzOkpOTI0mKjIw0OAkAALhcOTk58vf3v+B5k+1S1cjJWK1WpaSkyNfXVyaTqdreNzs7W5GRkTp06JD8/Pyq7X1REb/n2sPvunbwe64d/J5rR03+nm02m3JyctSgQQOZzReeycOIz1nMZrMiIiJq7P39/Pz4f6pawO+59vC7rh38nmsHv+faUVO/54uN9JzB5GYAAOA0KD4AAMBpUHxqicVi0aRJk2SxWIyO4tD4Pdcefte1g99z7eD3XDvqwu+Zyc0AAMBpMOIDAACcBsUHAAA4DYoPAABwGhQfAADgNCg+teSdd95RdHS0PDw81KFDB61evdroSA7n119/1e23364GDRrIZDLpu+++MzqSw3nppZfUvn17+fr6KjQ0VH379lVSUpLRsRzO9OnT1bp16/JF3jp27KiffvrJ6FgO7+WXX5bJZNJjjz1mdBSHM3nyZJlMpgpHs2bNDMlC8akFX3zxhcaNG6dJkyZp/fr1atOmjXr37q3U1FSjozmUvLw8tWnTRu+8847RURzWsmXLNGrUKK1atUqLFi1ScXGxevXqpby8PKOjOZSIiAi9/PLLWrdundauXavu3bvrzjvv1LZt24yO5rDWrFmj999/X61btzY6isNq2bKljh49Wn78/vvvhuTg6+y1oEOHDmrfvr3+85//SCrbDywyMlJjxozRhAkTDE7nmEwmk+bMmaO+ffsaHcWhpaWlKTQ0VMuWLdMNN9xgdByHVq9ePU2ZMkXDhg0zOorDyc3NVbt27fTuu+/q3//+t66++mq98cYbRsdyKJMnT9Z3332njRs3Gh2FEZ+aVlRUpHXr1ummm24qf85sNuumm27SypUrDUwGXLmsrCxJZX+UUTNKS0v1+eefKy8vTx07djQ6jkMaNWqUbr311gr/nkb12717txo0aKCYmBglJibq4MGDhuRgk9Ialp6ertLSUoWFhVV4PiwsTDt37jQoFXDlrFarHnvsMV1//fVq1aqV0XEczpYtW9SxY0cVFBTIx8dHc+bMUYsWLYyO5XA+//xzrV+/XmvWrDE6ikPr0KGDPvroIzVt2lRHjx7Vc889py5dumjr1q3y9fWt1SwUHwBVMmrUKG3dutWw+/SOrmnTptq4caOysrL09ddfa/DgwVq2bBnlpxodOnRIY8eO1aJFi+Th4WF0HIfWp0+f8n9u3bq1OnTooKioKH355Ze1fvuW4lPDgoOD5eLiouPHj1d4/vjx46pfv75BqYArM3r0aM2bN0+//vqrIiIijI7jkNzd3dWkSRNJUkJCgtasWaM333xT77//vsHJHMe6deuUmpqqdu3alT9XWlqqX3/9Vf/5z39UWFgoFxcXAxM6roCAAMXHx2vPnj21/rOZ41PD3N3dlZCQoMWLF5c/Z7VatXjxYu7Xw+7YbDaNHj1ac+bM0ZIlS9S4cWOjIzkNq9WqwsJCo2M4lB49emjLli3auHFj+XHNNdcoMTFRGzdupPTUoNzcXO3du1fh4eG1/rMZ8akF48aN0+DBg3XNNdfo2muv1RtvvKG8vDwNHTrU6GgOJTc3t8L/eti/f782btyoevXqqVGjRgYmcxyjRo3S7NmzNXfuXPn6+urYsWOSJH9/f3l6ehqcznE89dRT6tOnjxo1aqScnBzNnj1bS5cu1YIFC4yO5lB8fX3PmZ/m7e2toKAg5q1VsyeeeEK33367oqKilJKSokmTJsnFxUX3339/rWeh+NSC++67T2lpaXr22Wd17NgxXX311Zo/f/45E55xZdauXatu3bqVPx43bpwkafDgwfroo48MSuVYpk+fLknq2rVrhednzpypIUOG1H4gB5WamqoHHnhAR48elb+/v1q3bq0FCxaoZ8+eRkcDquTw4cO6//77lZGRoZCQEHXu3FmrVq1SSEhIrWdhHR8AAOA0mOMDAACcBsUHAAA4DYoPAABwGhQfAADgNCg+AADAaVB8AACA06D4AAAAp0HxAYCzmEwmfffdd0bHAFADKD4A6pQhQ4bIZDKdc9x8881GRwPgANiyAkCdc/PNN2vmzJkVnrNYLAalAeBIGPEBUOdYLBbVr1+/whEYGCip7DbU9OnT1adPH3l6eiomJkZff/11hddv2bJF3bt3l6enp4KCgjR8+HDl5uZWuObDDz9Uy5YtZbFYFB4ertGjR1c4n56err/97W/y8vJSXFycvv/++/JzJ06cUGJiokJCQuTp6am4uLhzihqAuoniA8Du/POf/1S/fv20adMmJSYmqn///tqxY4ckKS8vT71791ZgYKDWrFmjr776Sj///HOFYjN9+nSNGjVKw4cP15YtW/T999+rSZMmFX7Gc889p3vvvVebN2/WLbfcosTERGVmZpb//O3bt+unn37Sjh07NH36dAUHB9feLwBA1dkAoA4ZPHiwzcXFxebt7V3heOGFF2w2m80myTZixIgKr+nQoYNt5MiRNpvNZvvggw9sgYGBttzc3PLzP/74o81sNtuOHTtms9lstgYNGtgmTpx4wQySbM8880z549zcXJsk208//WSz2Wy222+/3TZ06NDq+cAAahVzfADUOd26ddP06dMrPFevXr3yf+7YsWOFcx07dtTGjRslSTt27FCbNm3k7e1dfv7666+X1WpVUlKSTCaTUlJS1KNHj4tmaN26dfk/e3t7y8/PT6mpqZKkkSNHql+/flq/fr169eqlvn37qlOnTlX6rABqF8UHQJ3j7e19zq2n6uLp6Vmp69zc3Co8NplMslqtkqQ+ffrowIED+t///qdFixapR48eGjVqlKZOnVrteQFUL+b4ALA7q1atOudx8+bNJUnNmzfXpk2blJeXV35++fLlMpvNatq0qXx9fRUdHa3FixdfUYaQkBANHjxYs2bN0htvvKEPPvjgit4PQO1gxAdAnVNYWKhjx45VeM7V1bV8AvFXX32la665Rp07d9ann36q1atXa8aMGZKkxMRETZo0SYMHD9bkyZOVlpamMWPGaNCgQQoLC5MkTZ48WSNGjFBoaKj69OmjnJwcLV++XGPGjKlUvmeffVYJCQlq2bKlCgsLNW/evPLiBaBuo/gAqHPmz5+v8PDwCs81bdpUO3fulFT2javPP/9cjzzyiMLDw/XZZ5+pRYsWkiQvLy8tWLBAY8eOVfv27eXl5aV+/fpp2rRp5e81ePBgFRQU6PXXX9cTTzyh4OBg3X333ZXO5+7urqeeekrJycny9PRUly5d9Pnnn1fDJwdQ00w2m81mdAgAqCyTyaQ5c+aob9++RkcBYIeY4wMAAJwGxQcAADgN5vgAsCvcnQdwJRjxAQAAToPiAwAAnAbFBwAAOA2KDwAAcBoUHwAA4DQoPgAAwGlQfAAAgNOg+AAAAKdB8QEAAE7j/wMFEFxzqUXY1QAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "np.random.seed(42)\n",
    "init_params = np.random.normal(0, 1, parameter_count)\n",
    "\n",
    "start_time = timeit.default_timer()\n",
    "result_vqe = minimize(objective_function,\n",
    "                      init_params,\n",
    "                      method='L-BFGS-B',\n",
    "                      jac=True,\n",
    "                      tol=1e-8,\n",
    "                      options={'maxiter': 5})\n",
    "end_time = timeit.default_timer()\n",
    "\n",
    "print('VQE-UCCSD energy= ', result_vqe.fun)\n",
    "print('Total elapsed time (s) = ', end_time - start_time)\n",
    "\n",
    "plt.plot(exp_vals)\n",
    "plt.xlabel('Epochs')\n",
    "plt.ylabel('Energy')\n",
    "plt.title('VQE')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, run the code again (the three previous cells) and specify `num_qpus` to be more than one if you have access to multiple GPUs and notice resulting speedup.  Thanks to CUDA-Q, this code could be used without modification in a setting where multiple physical QPUs were available."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Using an Active Space"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Performing electronic structure computations with all electrons and orbitals is often prohibitively expensive and unnecessary.  Most of the interesting chemistry can be modeled by restricting simulations to the highest energy occupied molecular orbitals and lowest energy unoccupied molecular orbitals. This is known as the active space approximation. \n",
    "\n",
    "Below is an example of STO-3G water modeled with a 4 electron 3 orbital active space simulated with UCCSD-VQE. Using an active space means you can run VQE for the same molecule using fewer qubits and a more shallow circuit.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![cas.png](./images/cas.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The molecule is defined the same way, expect for you now include variables `nele_cas` and `norb_cas` to define the active space. The `ncore` "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "geometry = [('O', (0.1173, 0.0, 0.0)), ('H', (-0.4691, 0.7570, 0.0)),\n",
    "            ('H', (-0.4691, -0.7570, 0.0))]\n",
    "basis = 'sto3g'\n",
    "multiplicity = 1\n",
    "charge = 0\n",
    "ncore = 3\n",
    "nele_cas, norb_cas = (4, 3)\n",
    "\n",
    "molecule = openfermionpyscf.run_pyscf(\n",
    "    openfermion.MolecularData(geometry, basis, multiplicity, charge))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The Hamiltonian is now constrcuted with the same steps, but only models the active space."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/tmp/ipykernel_225414/1900341958.py:9: ComplexWarning: Casting complex values to real discards the imaginary part\n",
      "  spin_ham = cudaq.SpinOperator(qubit_hamiltonian)\n"
     ]
    }
   ],
   "source": [
    "molecular_hamiltonian = molecule.get_molecular_hamiltonian(\n",
    "    occupied_indices=range(ncore),\n",
    "    active_indices=range(ncore, ncore + norb_cas))\n",
    "\n",
    "fermion_hamiltonian = get_fermion_operator(molecular_hamiltonian)\n",
    "\n",
    "qubit_hamiltonian = jordan_wigner(fermion_hamiltonian)\n",
    "\n",
    "spin_ham = cudaq.SpinOperator(qubit_hamiltonian)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Similarly, the kernel is defined only by the orbitals and electrons in the active space. Notice how this means you only need to optimize 8 parameters now. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "8\n"
     ]
    }
   ],
   "source": [
    "electron_count = nele_cas\n",
    "qubit_count = 2 * norb_cas\n",
    "\n",
    "\n",
    "@cudaq.kernel\n",
    "def kernel(qubit_num: int, electron_num: int, thetas: list[float]):\n",
    "    qubits = cudaq.qvector(qubit_num)\n",
    "\n",
    "    for i in range(electron_num):\n",
    "        x(qubits[i])\n",
    "\n",
    "    cudaq.kernels.uccsd(qubits, thetas, electron_num, qubit_num)\n",
    "\n",
    "\n",
    "parameter_count = cudaq.kernels.uccsd_num_parameters(electron_count,\n",
    "                                                     qubit_count)\n",
    "\n",
    "print(parameter_count)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "def cost(theta):\n",
    "\n",
    "    exp_val = cudaq.observe(kernel, spin_ham, qubit_count, electron_count,\n",
    "                            theta).expectation()\n",
    "    thetas = theta\n",
    "    return exp_val\n",
    "\n",
    "\n",
    "exp_vals = []\n",
    "\n",
    "\n",
    "def callback(xk):\n",
    "    exp_vals.append(cost(xk))\n",
    "\n",
    "\n",
    "# Initial variational parameters.\n",
    "np.random.seed(42)\n",
    "x0 = np.random.normal(0, 1, parameter_count)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The VQE procedure below is much faster using an active space compared to inclusion of all orbitals and electrons."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "UCCSD-VQE energy =   -74.96447648834257\n",
      "Total number of qubits =  6\n",
      "Total number of parameters =  8\n",
      "Total number of terms in the spin hamiltonian =  62\n",
      "Total elapsed time (s) =  3.2149466909468174\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAksAAAHHCAYAAACvJxw8AAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjMsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvZiW1igAAAAlwSFlzAAAPYQAAD2EBqD+naQAAWFVJREFUeJzt3Xd4VGX+Pv77TM2kTsqEVEooBqSIIkgwuIIS0F3ELuBCaH5RxF11XWHFykcRC7u/dVl32Y1BlxUVBYMNBCKuktCbtNBCSQgkTJKZ1MmU8/tjMicMaZOQ5Mwk9+u6zmVOzfswkNw+z3OeI4iiKIKIiIiIGqSQuwAiIiIib8awRERERNQEhiUiIiKiJjAsERERETWBYYmIiIioCQxLRERERE1gWCIiIiJqAsMSERERURMYloiIiIiawLBERERE1ASGJSLqtCZOnAh/f3+UlZU1eszUqVOh0WhgNBoBABUVFVi8eDEGDx4Mf39/hISEIDk5Gf/5z3/Q0NuhBEFodJk7d2673RsRdRyV3AUQEbWXqVOn4quvvsK6deswbdq0evsrKyuRkZGB8ePHIzw8HJcuXcLYsWNx9OhRPPLII3jyySdRXV2NL774AtOmTcOGDRvwn//8BwqF+/9n3nnnnQ1ev1+/fu12b0TUcRiWiKjTmjhxIoKCgvDxxx83GGYyMjJQUVGBqVOnAgCmT5+Oo0ePYt26dZg4caJ03FNPPYXnnnsO77zzDm644QY899xzbtfp168fHn300fa9GSKSDbvhiKjT0ul0uO+++7BlyxYUFhbW2//xxx8jKCgIEydOxPbt27Fx40akpqa6BSWXJUuWoG/fvnjzzTdRVVXVEeUTkZdgWCKiTm3q1Kmw2Wz47LPP3LYXFxdj48aNuPfee6HT6fDVV18BQIMtUACgUqkwZcoUFBcXIysry21fdXU1Ll++XG+pqalpn5siog7FsEREndqYMWMQHR2Njz/+2G37mjVrYLVapS64I0eOAACGDBnS6LVc+1zHuqSlpcFgMNRb1q5d25a3QkQy4ZglIurUlEolHnnkEfz5z3/GmTNn0LNnTwDOLrhu3bph7NixACA9MRcUFNTotVz7rn667p577sGTTz5Z7/hBgwa1xS0QkcwYloio05s6dSr+/Oc/4+OPP8af/vQn5OXl4aeffsJTTz0FpVIJwD0I6fX6Bq/jCkmRkZFu2+Pi4nDHHXe03w0QkazYDUdEnd5NN92ExMRErF69GgCwevVqiKIodcEBwIABAwAABw8ebPQ6rn0JCQntWC0ReRuGJSLqEqZOnYpDhw7h4MGD+Pjjj9G3b1/cfPPN0v7f/OY3AICPPvqowfPtdrvUdTd69OgOqZmIvAPDEhF1Ca5WpJdeegn79+93a1UCgFtuuQXjxo1Deno6vv7663rnv/DCCzh+/Dj++Mc/QqXiCAairkQQG5q/n4ioExo1apT02P+JEyfQp08ft/2XLl3CmDFjcOzYMUyZMgXJycmwWCxYu3Yttm7dikcffRQfffQRBEGQzhEEodEZvLt164Y777yzfW+KiNodwxIRdRl///vfMW/ePAwfPhw7duxo8Jjy8nIsW7YMn332GU6dOoXq6moAwIsvvojXXnut3vFXBqer3Xbbbdi6dWub1E5E8mFYIiJqQn5+PpKSkmCz2ZCdnY3u3bvLXRIRdTCOWSIiakJsbCw2bNiA6upqTJgwASUlJXKXREQdjC1LRERERE1gyxIRERFRExiWiIiIiJrAsERERETUBIYlIiIioiZwGto24HA4cOHCBQQFBTU55woRERF5D1EUUVZWhpiYGCgUjbcfMSy1gQsXLiA+Pl7uMoiIiKgVzp8/j7i4uEb3+0RY2rp1K26//fYG9+3cudPtZZgAcPLkSQwdOhRKpRKlpaWNXvfAgQN488038fPPP+Py5cvo2bMn5s6di9/97nctqi8oKAiA8w87ODi4RecSERGRPMxmM+Lj46Xf443xibCUlJSEgoICt20vvvgitmzZgmHDhrltt1qtmDx5MpKTk6V3QDVmz549iIyMxKpVqxAfH4+srCw89thjUCqVePLJJz2uz9X1FhwczLBERETkY5obQuMTYUmj0SAqKkpat1qtyMjIwPz58+vd4KJFi5CYmIixY8c2G5Zmzpzptp6QkIDs7GysXbu2RWGJiIiIOi+ffBpu/fr1MBqNmDFjhtv2zMxMrFmzBsuXL2/1tU0mE8LCwq61RCIiIuokfKJl6WppaWlISUlxG4xlNBqRmpqKVatWtborLCsrC59++im++eabJo+zWCywWCzSutlsbtX3IyIiIu8na8vSggULIAhCk8uxY8fczsnLy8PGjRsxa9Yst+1z5szBlClTMHr06FbVcujQIdxzzz14+eWXMW7cuCaPXbJkCUJCQqSFT8IRERF1XrK+SLeoqAhGo7HJYxISEqDRaKT1xYsX47333kN+fj7UarW0Xa/Xo7y8XFoXRREOhwNKpRIrVqyoNz7pSkeOHMHtt9+O2bNn4/XXX2+27oZaluLj42EymTjAm4iIyEeYzWaEhIQ0+/tb1m44g8EAg8Hg8fGiKCI9PR3Tpk1zC0oAkJ2dDbvdLq1nZGRg6dKlyMrKQmxsbKPXPHz4MMaMGYPp06d7FJQAQKvVQqvVelw3ERER+S6fGrOUmZmJ3NxczJ49u96+/v37u63v3r0bCoUCAwcOlLatW7cOCxculLr2Dh06hDFjxiAlJQXPPPMMLl68CABQKpUtCnFERETUefnU03BpaWlISkpCYmJiq843mUzIycmR1j///HMUFRVh1apViI6OlparJ7kkIiKirkvWMUudhad9nkREROQ9PP397VMtS0REREQdjWGJiIiIqAkMS0RERERN8Kmn4bqaQnM1LDYHDEFa+KmVcpdDRETUJbFlyYs9vGI7kt/6Ab/km+QuhYiIqMtiWPJiSoUAALDZ+cAiERGRXBiWvJiqNizZHQxLREREcmFY8mIqpTMsWR0OmSshIiLquhiWvJhK4fx47OyGIyIikg3DkhdzdcPZ2LJEREQkG4YlLyYN8OaYJSIiItkwLHkxtbK2G45hiYiISDYMS17M1bJk5ZglIiIi2TAseTG10jV1AMcsERERyYVhyYuxZYmIiEh+DEteTMUxS0RERLJjWPJiKj4NR0REJDuGJS9W9244jlkiIiKSC8OSF1PXzuDNliUiIiL5MCx5MaXS1bLEsERERCQXhiUvplZw6gAiIiK5MSx5MWVtN5yV3XBERESyYVjyYippUkqGJSIiIrkwLHkxaeoAjlkiIiKSDcOSF6ubZ4ljloiIiOTCsOTFXDN4c+oAIiIi+TAseTFOSklERCQ/hiUvplbydSdERERyY1jyYq6pA/g0HBERkXwYlrwYn4YjIiKSH8OSF1Mp+TQcERGR3BiWvBhbloiIiOTHsOTFVApOHUBERCQ3hiUvxm44IiIi+TEseTElu+GIiIhk5xNhaevWrRAEocFl165d9Y4/efIkgoKCoNfrm7yu0WjE+PHjERMTA61Wi/j4eDz55JMwm83tdCcto+LUAURERLLzibCUlJSEgoICt2X27Nno1asXhg0b5nas1WrF5MmTkZyc3Ox1FQoF7rnnHqxfvx7Hjx/HypUrsXnzZsydO7e9bqVFXAO8rQxLREREslHJXYAnNBoNoqKipHWr1YqMjAzMnz8fgiC4Hbto0SIkJiZi7NixyMrKavK6oaGhePzxx6X1Hj164IknnsDbb7/dtjfQSq4xS3aOWSIiIpKNT4Slq61fvx5GoxEzZsxw256ZmYk1a9Zg//79WLt2bYuve+HCBaxduxa33XZbk8dZLBZYLBZpvb267aSn4ThmiYiISDY+0Q13tbS0NKSkpCAuLk7aZjQakZqaipUrVyI4OLhF15s8eTL8/f0RGxuL4OBg/Pvf/27y+CVLliAkJERa4uPjW3UfzVHx3XBERESykzUsLViwoNGB267l2LFjbufk5eVh48aNmDVrltv2OXPmYMqUKRg9enSL6/jzn/+MvXv3IiMjA6dOncIzzzzT5PELFy6EyWSSlvPnz7f4e3rCNWaJA7yJiIjkI4iiKNtv4qKiIhiNxiaPSUhIgEajkdYXL16M9957D/n5+VCr1dJ2vV6P8vJyaV0URTgcDiiVSqxYsQIzZ870qKaff/4ZycnJuHDhAqKjoz06x2w2IyQkBCaTqcWtWk3Zd64E9/49C3GhOvz8/Jg2uy4RERF5/vtb1jFLBoMBBoPB4+NFUUR6ejqmTZvmFpQAIDs7G3a7XVrPyMjA0qVLkZWVhdjYWI+/h6N2MPWVY5LkolZy6gAiIiK5+dQA78zMTOTm5mL27Nn19vXv399tfffu3VAoFBg4cKC0bd26dVi4cKHUtfftt9/i0qVLuPnmmxEYGIjDhw/jueeew6hRo9CzZ892vRdPuCaltHKANxERkWx8KiylpaUhKSkJiYmJrTrfZDIhJydHWtfpdPjXv/6Fp59+GhaLBfHx8bjvvvuwYMGCtir5mqg5dQAREZHsZB2z1Fm015il3MsVuP2drQjSqvDLqyltdl0iIiLy/Pe3T04d0FW4nobj1AFERETyYVjyYnUzeDMsERERyYVhyYtJA7w5ZomIiEg2DEteTF37uhNRBBxsXSIiIpIFw5IXUyrrXhLM1iUiIiJ5MCx5MVfLEsBxS0RERHJhWPJirjFLAJ+IIyIikgvDkhdTXRmWOIs3ERGRLBiWvJhCIcCVl2wcs0RERCQLhiUvp6odt8SWJSIiInkwLHk5TkxJREQkL4YlL6fkK0+IiIhkxbDk5dRKVzccxywRERHJgWHJy7FliYiISF4MS17ONX0AB3gTERHJg2HJy7kGeHPqACIiInkwLHk5aeoAdsMRERHJgmHJy7EbjoiISF4MS17ONcCb8ywRERHJg2HJy7nGLFk5ZomIiEgWDEtezjVmyc5uOCIiIlkwLHk5acwSW5aIiIhkwbDk5eqmDmDLEhERkRwYlryc1A3HsERERCQLhiUvJw3w5pglIiIiWTAseTmVNHUAxywRERHJgWHJy7nmWWLLEhERkTwYlrycSskxS0RERHJiWPJyKqllid1wREREcmBY8nJ8Go6IiEheDEterm5SSoYlIiIiOTAseTmla1JKDvAmIiKSBcOSl1Nz6gAiIiJZMSx5OWXtmCUru+GIiIhk4RNhaevWrRAEocFl165d9Y4/efIkgoKCoNfrPf4eRqMRcXFxEAQBpaWlbVf8NVIrXS1LDEtERERy8ImwlJSUhIKCArdl9uzZ6NWrF4YNG+Z2rNVqxeTJk5GcnNyi7zFr1iwMHjy4LctuE65JKTlmiYiISB4+EZY0Gg2ioqKkJTw8HBkZGZgxYwYEQXA7dtGiRUhMTMRDDz3k8fXff/99lJaW4g9/+ENbl37NXJNS2jhmiYiISBY+EZautn79ehiNRsyYMcNte2ZmJtasWYPly5d7fK0jR47gtddew0cffQSFwvv+ODh1ABERkbxUchfQGmlpaUhJSUFcXJy0zWg0IjU1FatWrUJwcLBH17FYLJg8eTLefvttdO/eHadPn/b4PIvFIq2bzeaW3UAL1HXDsWWJiIhIDrI2pSxYsKDRgduu5dixY27n5OXlYePGjZg1a5bb9jlz5mDKlCkYPXq0x99/4cKF6N+/Px599NEW1b1kyRKEhIRIS3x8fIvObwnXAG+2LBEREclDEEVRtt/CRUVFMBqNTR6TkJAAjUYjrS9evBjvvfce8vPzoVarpe16vR7l5eXSuiiKcDgcUCqVWLFiBWbOnFnv2jfccAN++eUXadzTlee88MILePXVVxusqaGWpfj4eJhMJo9btTyV9nMuFn99BBOHxOCvk4c2eWxxRQ3mrtqDB2+Kw4PD2i/AERERdQZmsxkhISHN/v6WtRvOYDDAYDB4fLwoikhPT8e0adPcghIAZGdnw263S+sZGRlYunQpsrKyEBsb2+D1vvjiC1RVVUnru3btwsyZM/HTTz+hd+/ejdah1Wqh1Wo9rvtatGTqgB2njdiZWwyL1c6wRERE1EZ8asxSZmYmcnNzMXv27Hr7+vfv77a+e/duKBQKDBw4UNq2bt06LFy4UOrauzoQXb58WbpWS+Zoak/SmCUPnoartjnDYmmVtV1rIiIi6kq87/GvJqSlpSEpKQmJiYmtOt9kMiEnJ6eNq2pfqhbMs1Rjcwaq0kqGJSIiorbiU2Hp448/xrZt2zw6NjU1td5M3KmpqWhqiNavfvUriKLoNa1KAKBSuOZZ8jwsmautnPGbiIiojfhUWOqKVErPu+EstWFJFAEzu+KIiIjaBMOSl5NaljzohnOFJYDjloiIiNoKw5KXcw3w9qRbrebKsFRZ0241ERERdSUMS17ONXWA1ZOwZL8yLLFliYiIqC0wLHm5upYlD8YsWa/shmPLEhERUVtgWPJyLRmzVHPFpJxsWSIiImobDEteTtWCd8O5j1liWCIiImoLDEteTtWCAd4WDvAmIiJqcwxLXk6ldH5EVnvzY5ZqOHUAERFRm2NY8nItaVliNxwREVHbY1jycq4xS1aPBnizZYmIiKitMSx5OVUrpw4wccwSERFRm2BY8nLKlrzu5IqWpRJ2wxEREbUJhiUv52pZaunUAeZqq0fjnIiIiKhpDEtezjVmybMB3nWTUooiUFbN1iUiIqJrxbDk5VwzeFuvGLNUY3Ng//nSegHqynmWAD4RR0RE1BYYlrycqxtOFAFHbTj654+nMGn5Nny667zbsTVXhaUSDvImIiK6ZgxLXk5Z2w0H1LUuHS8sBwCcK650O9Y1dYCf2vmxcvoAIiKia8ew5OXUirqPyNXtVlLhbDGqqrG5HeuaOiAyyA8AYGI3HBER0TVjWPJySkVdy5LriThjbViqrLG7HetqWeoWrAXA98MRERG1BYYlL6e6MizVzrVUXGEBAFRZ68KS3SFKLU+uliXOtURERHTtGJa8nEIhwJWXbA4HRFFESYUzBFVd0bJ05eBuQ5CzZcnEMUtERETXjGHJB6iUdbN4l1tsUnfbld1wlivmWIpkNxwREVGbYVjyAXXvhxNRXFEXgK7shnO1LCkEICKgNiyxZYmIiOiaqeQugJrnGuRttTvcAlCVW8uSMyxpVAro/dUAOGaJiIioLbBlyQeoa7vh7A5RmjYAACqtdVMHSGFJqUBYgAYAYCy3dGCVREREnRPDkg9QXvEyXeOV3XA17q9AAQCtWon4MH8AwIXSqnqzehMREVHLMCz5ALUrLNmvGrN0xaSUrkHfGqUCkUFa+KkVcIhAfmlVxxZLRETUyTAs+QDXK09sDodbWKq02iGKzrmVpJYllQKCIKBneAAA4IyxooOrJSIi6lwYlnyAqvaVJ7arnoYTxbqxSq6pAzQq57E9wp1dcWcvu4elaqsd5moO/CYiIvIUw5IPUDXSDQfUPRF3ZcsSAKll6exVL9udtHwbbnvrB7cn6YiIiKhxDEs+QHnFPEvGq8JSpdU9LNW1LNWGJWNdWKq22nHsYhlKKq0cy0REROQhhiUf4Jo6wOpwSO+Fc5FaluxXhyVnN9yVY5aKyurOLWNXHBERkUcYlnyA1LJkr3svnIsrLFmsrm44JYC6sHS+uFJ6we7l8ivDkg1ERETUPIYlH6CufRqu3GJDucUZciICnRNPVtZOH2C5YuoAAIgO0UGjVMBqF3GhtsvtcnldFx7DEhERkWd8Iixt3boVgiA0uOzatave8SdPnkRQUBD0en2z127omp988kk73EXrdQ9zjj/amlMIwNnSFBnkB6Du/XBXj1lSKgTEh+kA1I1bcm9ZYjccERGRJ3wiLCUlJaGgoMBtmT17Nnr16oVhw4a5HWu1WjF58mQkJyd7fP309HS3a0+aNKmN7+DajO4XAQD4/sglAECovwYBWmd3m9QNd9XUAQDqzbV0+YoxS5w+gIiIyDM+8SJdjUaDqKgoad1qtSIjIwPz58+HIAhuxy5atAiJiYkYO3YssrKyPLq+Xq93u763GdXHGZYqa4NReIAGOo3KbdvVUwcAVz4RVxuWOGaJiIioxXyiZelq69evh9FoxIwZM9y2Z2ZmYs2aNVi+fHmLrjdv3jxERERg+PDh+OCDD6RZsRtjsVhgNpvdlvYUEajFgOhgaT0sQAN/dW3LUiPdcADQM8L1RJyrG45jloiIiFrKJ1qWrpaWloaUlBTExcVJ24xGI1JTU7Fq1SoEBwc3cba71157DWPGjIG/vz++//57PPHEEygvL8dTTz3V6DlLlizBq6++ek330FLJ/SJwpMAZysICNFIounpSSk0TLUtF5eyGIyIiailZW5YWLFjQ6MBt13Ls2DG3c/Ly8rBx40bMmjXLbfucOXMwZcoUjB49ukU1vPjiixg1ahSGDh2K559/Hn/84x/x9ttvN3nOwoULYTKZpOX8+fMt+p6tMbqvQfo6LEADncbZslQpjVlynzoAAHqG17UsORwiu+GIiIhaQdaWpWeffRapqalNHpOQkOC2np6ejvDwcEycONFte2ZmJtavX4933nkHACCKIhwOB1QqFVasWIGZM2d6VNOIESOwePFiWCwWaLXaBo/RarWN7msvN/UIhZ9agWqrA2EBGlTUTiFwdTfclWOWYvU6qJUCamwO5JdWuQ/wrmLLEhERkSdkDUsGgwEGg6H5A2uJooj09HRMmzYNarXabV92djbs9rr3nWVkZGDp0qXIyspCbGysx99j//79CA0N7fAw1Bw/tRJJvSOQeawQsaE6nK9951tV7TxLNVfNswQAKqUCPcIDcLKwHDkXy2C+ojWJLUtERESe8akxS5mZmcjNzcXs2bPr7evfv7/b+u7du6FQKDBw4EBp27p167Bw4UKpa++rr77CpUuXcMstt8DPzw+bNm3CG2+8gT/84Q/teyOt9H+TBuK7Qxcx6YZY/Pvn0wCu7IarP3UAACREOMPSrrPFbtvLLGxZIiIi8oRPhaW0tDQkJSUhMTGxVeebTCbk5ORI62q1GsuXL8fTTz8NURTRp08fLFu2DHPmzGmrkttUjF6HWbf2AoBGn4bTXh2WDIEALmFn7lVhiS1LREREHvGpsPTxxx97fGxqamq98VBXbxs/fjzGjx/fRtV1LNcA76qrBnjXa1kyOJ+I+yXPBMD5mpTL5TUoq7ZBFMV681QRERGRO5+cZ4nQ6KSUV4el3rVhyVb7Mt2EiEAAgN0hSucSERFR4xiWfJTuqm64hqYOAOrCkUtcmA5KhbM1iV1xREREzWNY8lH+V3XDNdayFBqggd6/7slBQ6AWQX7OVim+TJeIiKh5DEs+SpqU0tr41AEuCREB0tcRV4QlM1uWiIiImsWw5KOkbrgaZ0iSnoZTNxCWDHVdcRFBGgT7OVua2LJERETUPIYlH1XXDedsHZLmWWqoZcnAliUiIqLWYljyUXXdcHaIotjoPEuA+yBvZ1hiyxIREZGnGJZ8lKsbThSdT8I1NsAbqJs+AHBvWeLTcERERM3zqUkpqY6/pu6jq6qxNzp1AAD0CA9ATIgf1CoFwgI4ZomIiKglGJZ8lFIhQKNSoMbmQLnFJk062VDLkkalwOZnb4MAAUqFwJYlIiKiFmBY8mE6tRI1NgdMVXUtRA2FJcC9JaquZYlhiYiIqDkcs+TDXE/EXRmWGhrgfTXpabgqdsMRERE1h2HJh7meiCutdIYeQQBUiuZfjBvEliUiIiKPMSz5MNcTcaVVNQCccywJgidhyTXPEluWiIiImtOqsFRRUdHWdVAr+F/VstTYeKWrcYA3ERGR51oVlrp164aZM2fi559/but6qAV0tYO2XWOWGpo2oCHBOk4dQERE5KlWhaVVq1ahuLgYY8aMQb9+/fDmm2/iwoULbV0bNUNX+x640kpnN5wng7uBK1qWLDY4aqccICIiooa1KixNmjQJX375JfLz8zF37lx8/PHH6NGjB379619j7dq1sNnYvdMR/K9qWfK0G841dYAoAhU17p/VD8cK8ev3fsKhfFMbVkpEROS7rmmAt8FgwDPPPIODBw9i2bJl2Lx5Mx544AHExMTgpZdeQmVlZVvVSQ1wPQ1XUunqhvPs49SqFFArnQPBrx639FH2GRzKN+ONb4+2YaVERES+65rC0qVLl/DWW29hwIABWLBgAR544AFs2bIF7777LtauXYtJkya1UZnUENfTcKYWDvAWBEGaPuDqJ+KOFJgBAFmnjDiYV9pGlRIREfmuVs3gvXbtWqSnp2Pjxo0YMGAAnnjiCTz66KPQ6/XSMUlJSejfv39b1UkNcD0NV2CqAuCcOsBTMXo/FFfU4FRhBRKjggEAl8stuGS2SMf848dT+PvUm9qwYiIiIt/TqpalGTNmICYmBtu2bcP+/fvx5JNPugUlAIiJicELL7zQFjVSI5L7GgAA5tquNE9blgBgaHwoAGDvuRJp29HaVqXg2gHg3x26iNzLnCaCiIi6tlaFpYKCAvzzn//EzTff3OgxOp0OL7/8cqsLo+YN7xWGJ37VW1r3dMwSANzYQw8A2HdFWDpywRmWkvsaMCYxEqIIfLEnr22KJSIi8lGt6oaz2Wwwm831tguCAK1WC41Gc82FkWeeubMfdp8twc7cYgRoPf84XS1Lh/LNsNjs0KqU0nilATHBEAQg81ghLpqr26VuIiIiX9GqliW9Xo/Q0NB6i16vh06nQ48ePfDyyy/D4XC0db10FZVSgb9PvRGzb+2Fubf1bv6EWj3C/REWoEGN3YHDtS1KrpalAdHBCPN3Bt6Sipq2L5qIiMiHtKplaeXKlXjhhReQmpqK4cOHAwB27tyJDz/8EIsWLUJRURHeeecdaLVa/OlPf2rTgqm+iEAtFv16QIvOEQQBQ+P12HKsEPvOlWJAdDBOFZUDAK6PCYbF5gy6xZUMS0RE1LW1Kix9+OGHePfdd/HQQw9J237zm99g0KBB+Oc//4ktW7age/fueP311xmWvNiNPUKx5Vgh9p4rwbAeoXCIQESgBoYgLcICnC1LrvfOERERdVWt6obLysrC0KFD620fOnQosrOzAQC33norzp07d23VUbsa2l0PANh/rlQar9Q/OhiCICAswDkPUzG74YiIqItrVViKj49HWlpave1paWmIj48HABiNRoSGhl5bddSuhsTpoRCA/NIq/OPHUwCcg7sBQF87ZslUZYXNzrFnRETUdbWqG+6dd97Bgw8+iO+++06aPmD37t04duwYPv/8cwDArl278PDDD7ddpdTmArQqDIoNwYE8E84ana+mGd4zDACg16ml40xVVoQHamWpkYiISG6tCksTJ05ETk4O/vnPfyInJwcAMGHCBHz55Zfo2bMnAODxxx9vsyKp/bz70BBsPlqIUH81EgyBGNbD2RqoUioQolPDVGVFSWUNwxIREXVZLQ5LVqsV48ePxz/+8Q8sWbKkPWqiDtQnMgh9IoMa3Bfq7wxLxRUc5E1ERF1Xi8csqdVqHDx4sD1qIS8TWvtEXAmnDyAioi6sVQO8H3300QYHeFPnwokpiYiIruF1Jx988AE2b96Mm266CQEBAW77ly1b1ibFuWzduhW33357g/t27txZ7x11J0+exNChQ6FUKlFaWtrs9VeuXIlly5bh+PHjCA4OxoMPPojly5e3Rek+zfVEHCemJCKirqxVYenQoUO48cYbAQDHjx932ycIwrVXdZWkpCQUFBS4bXvxxRexZcsWDBs2zG271WrF5MmTkZycjKysrGavvWzZMrz77rt4++23MWLECFRUVODMmTNtWb7Pcs21xIkpiYioK2tVWPrhhx/auo4maTQaREVFSetWqxUZGRmYP39+vXC2aNEiJCYmYuzYsc2GpZKSEixatAhfffUVxo4dK20fPHhw296Aj3KNWeLElERE1JW1asySy8mTJ7Fx40ZUVVUBAERRbJOimrN+/XoYjUbMmDHDbXtmZibWrFnjcRfapk2b4HA4kJ+fj/79+yMuLg4PPfQQzp8/3+R5FosFZrPZbemMQjlmiYiIqHVhyWg0YuzYsejXrx/uuusuqYts1qxZePbZZ9u0wIakpaUhJSUFcXFxbjWlpqZi5cqVCA4O9ug6p0+fhsPhwBtvvIG//OUv+Pzzz1FcXIw777wTNTWNB4QlS5YgJCREWlyzlnc2UljimCUiIurCWhWWnn76aajVapw7dw7+/v7S9ocffhgbNmzw+DoLFiyAIAhNLseOHXM7Jy8vDxs3bsSsWbPcts+ZMwdTpkzB6NGjPf7+DocDVqsVf/3rX5GSkoJbbrkFq1evxokTJ5rsaly4cCFMJpO0NNcS5atC/Z1jlko4ZomIiLqwVo1Z+v7777Fx40a3lh0A6Nu3L86ePevxdZ599lmkpqY2eUxCQoLbenp6OsLDwzFx4kS37ZmZmVi/fj3eeecdAM4uQYfDAZVKhRUrVmDmzJn1rh0dHQ0AGDBggLTNYDAgIiKiyZcAa7VaaLWdf0brMI5ZIiIial1YqqiocGtRcikuLm5RiDAYDDAYDB4fL4oi0tPTMW3aNKjVard92dnZsNvt0npGRgaWLl2KrKwsxMbGNni9UaNGAQBycnKk4FdcXIzLly+jR48eHtfVWbkGeJurnS/TVSmvaYgbERGRT2rVb7/k5GR89NFH0rogCHA4HHjrrbcanQ+pLWRmZiI3NxezZ8+ut69///4YOHCgtMTGxkKhUGDgwIEIDXW+72zdunVITEyUzunXrx/uuece/O53v0NWVhYOHTqE6dOnIzExsV3vw1e4XqYris6X6RIREXVFrWpZeuuttzB27Fjs3r0bNTU1+OMf/4jDhw+juLgY27Zta+saJWlpaUhKSnILPC1hMpmkF/+6fPTRR3j66adx9913Q6FQ4LbbbsOGDRvqtVx1RSqlAsF+KpirbXyZLhERdVmC2Mrn/U0mE/72t7/hwIEDKC8vx4033oh58+ZJ44C6ErPZjJCQEJhMJo+fxPMVv3r7B5wxVmLN3JG4uWeY3OUQERG1GU9/f7eqZQkAQkJC8MILL7T2dPIRen8NYKzkIG8iIuqyWh2WSktLsXPnThQWFsLhcLjtmzZt2jUXRt7B9UQcJ6YkIqKuqlVh6auvvsLUqVNRXl6O4OBgt1eOCILAsNSJ1E1MyQHeRETUNbXqabhnn30WM2fORHl5OUpLS1FSUiItxcXFbV0jyahuYkpny1JxRQ3uWPYj/rL5eFOnERERdRqtCkv5+fl46qmnGpxriToX11xLxnJnWMo+ZcTJwnJ8vidPzrKIiIg6TKvCUkpKCnbv3t3WtZAXig7xAwDkl1YCAM4VO/97yVwNh6PhByk3H7mEldtyO6ZAIiKidtaqMUt33303nnvuORw5cgSDBg2qNyfR1a8iId/VI9zZenjO6AxJ50uc/7XaRVwutyAy2M/t+Evmajzx372osTtwa18D+kQGdmzBREREbaxVYWnOnDkAgNdee63ePkEQ3F47Qr6te1gAAKDAXA2LzY7ztS1LAJBfWlUvLP37p9OosTuk/QxLRETk61rVDedwOBpdGJQ6l4hADfw1SogikFdShbySKmlfgana7djSyhp8vKPuBcSFZvf9REREvqhFYemuu+6CyWSS1t98802UlpZK60ajEQMGDGiz4kh+giCge5izK+7M5QrkldS1LF0orXI79qPss6ioqQvLhWWWjimSiIioHbUoLG3cuBEWS90vwDfeeMNtqgCbzVbv3Wvk+1zjlnaeKYbVXjeo+0JpXcuR3SFiZdYZAEBChLPrrohhiYiIOoEWhaWrXyPXytfKkY/pEe4MP9tOXnbbXmCqa1k6VVSO4ooa+GuUmDKiOwCgsIzdcERE5PtaNWaJupb42m64wxfMAABF7YTtV3bDHThfCgAYGBOCGL0OAFBoZssSERH5vhaFJUEQ3F5t4tpGnVuP2rDkakgcGBsCALhwxQDvg3nOsWyD40IQGaQFwDFLRETUObRo6gBRFJGamgqt1vnLsLq6GnPnzkVAgLOb5srxTNR5uMYsuQzvGYaDeSYUlVlgsdmhVSlxMK8UADA4Xo/IIOd0AoVl1RBFkYGaiIh8WovC0vTp093WH3300XrH8CW6nU+MXgelQoC9dsbuwfF6aFUKWGwOXDJZEBXih6MFZQCAIXEhiAyuDdNWB8zVNoTo1I1em4iIyNu1KCylp6e3Vx3kxdRKBWL1OulVJ93D/BGj1yH3cgUumKpgqrKixu6A3l+N7mH+EAQBQX4qlFXbUFRWzbBEREQ+jQO8ySNXdsV1D/OX3hl3obQKB2q74AbFhkhdbtK4JQ7yJiIiH8ewRB5xTUwZoFEi1F8tPfFWYKqWxisNidNLx9eNW2JYIiIi38awRB5xtSzF13azxVzRsnTlk3AurnFLnGuJiIh8XatepEtdz/Be4RAEYHivMACQWpb+d6II+bXvixsSr5eOZzccERF1FgxL5JEb4vXY/cIdCPXXAACia8PS+WJnUPrVdQZ0C/aTjmc3HBERdRYMS+Sx8ECt9HW/boFQCM4n5X5/Rz/MTu7ldiy74YiIqLNgWKJWiQ7R4dvfJUOv0yAqxK/efgNn8SYiok6CYYlaLTEquNF9rm64Io5ZIiIiH8en4ahduLrhyiw2VNXYZa6GiIio9RiWqF0EaVXwUzv/enHcEhER+TKGJWoXgiBIT8cdyjfLXA0REVHrMSxRuxk/MAoA8OaGo/W64tg1R0REvoJhidrNU2P6IjrED+eLq/C3H05I21fvPIf+L23At78UyFgdERGRZxiWqN0EaFV4ZeL1AIAV/zuNU0XlEEUR//rpNABg4+GLcpZHRETkEYYlalfjBnTDr64zwGoX8e+fTuPwBTNOF1UAAI5c4FgmIiLyfpxnidqVIAh44ld9sDWnCGv35sPuEKV9p4rKUW21w0+tlLFCIiKiprFlidrdzT1DMTA2GBabA5/tzpO2O0Qg52KZjJURERE1zyfC0tatWyEIQoPLrl276h1/8uRJBAUFQa/XN3ndlStXNnrdwsLCdrqbrkcQBMy6te7dcUFaFW7uGQoAOFLArjgiIvJuPhGWkpKSUFBQ4LbMnj0bvXr1wrBhw9yOtVqtmDx5MpKTk5u97sMPP1zvuikpKbjtttsQGRnZXrfTJd09KAaRte+LG3d9FIZ2d4alowxLRETk5XxizJJGo0FUVJS0brVakZGRgfnz50MQBLdjFy1ahMTERIwdOxZZWVlNXlen00Gn00nrRUVFyMzMRFpaWtveAEGjUmDRrwfgb5kn8NjoBCkkHblgxp6zxUj9YBcMwVqM7mvAtJE9kGAIlLliIiIiJ58IS1dbv349jEYjZsyY4bY9MzMTa9aswf79+7F27doWX/ejjz6Cv78/HnjggbYqla4wcUgMJg6JAQC4Mu7RAjP+75ujKLPYUFZkw+miChzIK8W6J0bJWCkREVEdn+iGu1paWhpSUlIQFxcnbTMajUhNTcXKlSsRHBzc6utOmTLFrbWpIRaLBWaz2W2hlkmICIBGpUBFjR37zpXCT63Am/cNAgAcOF8Kc7VV5gqJiIicZA1LCxYsaHSAtWs5duyY2zl5eXnYuHEjZs2a5bZ9zpw5mDJlCkaPHt2qWrKzs3H06NF6123IkiVLEBISIi3x8fGt+p5dmUqpQGJUkLQ+fWRPPDK8O3qG+8MhArtyi2WsjoiIqI4giqLY/GHto6ioCEajscljEhISoNFopPXFixfjvffeQ35+PtRqtbRdr9ejvLxcWhdFEQ6HA0qlEitWrMDMmTOb/D6zZs3C3r17sW/fvmbrtlgssFgs0rrZbEZ8fDxMJlOrW7W6ogVfHMQnu84jQKPET8+PQViABgvXHsTqnecx+9ZeWPTrAXKXSEREnZjZbEZISEizv79lHbNkMBhgMBg8Pl4URaSnp2PatGluQQlwtgzZ7XUvZ83IyMDSpUuRlZWF2NjYJq9bXl6Ozz77DEuWLPGoDq1WC61W63Hd1LCUgVH4ZNd5PH1nP4QFOAPxLQnhWL3zPLJPNx2iiYiIOopPDfDOzMxEbm4uZs+eXW9f//793dZ3794NhUKBgQMHStvWrVuHhQsX1uva+/TTT2Gz2fDoo4+2T+HUoNuvi0TO/42HVlU3g/fIhHAAzvmXSitroPfXNHY6ERFRh/CpAd5paWlISkpCYmJiq843mUzIyclp8Lr33Xdfs5NYUtu7MigBQGSwH3obAiCKwA6OWyIiIi8g65ilzsLTPk/yzKIvf8Gq7eeQmtQTr0y8Xu5yiIiok/L097dPtSxR1zAyIQIA8ENOIWx2h8zVEBFRV8ewRF5ndL8I6P3VOGusxJo9ec2fQERE1I4YlsjrBPmp8dSYvgCAd78/jnKLTeaKiIioK2NYIq/06C090CPcH5fLLVjx4ym5yyEioi6MYYm8kkalwILxzqceP9h2BlaOXSIiIpkwLJHXSrk+Cnp/NcotNvySb5K7HCIi6qIYlshrKRQCRvQKAwBkn+KM3kREJA+GJfJqt9TO6L299vUna/fm4d3vc1BcUSNnWURE1IX41OtOqOsZ2dsZlnafKcHJwjL8Yc0BOETgP9vP4oW7+uPBYfEyV0hERJ0dW5bIq/WLDEJYgAZVVjt+/+l+OERApRBQWmnFc58fxIHzpXKXSEREnRzDEnk1hULALQnOcUuH8s0AgI9mDcfYxEgAwLeHCmSrjYiIugaGJfJ6rnFLAHBjdz1GJoTj3htjAQCbDl+SqywiIuoiGJbI6428IizNva03BEHAbf0M0CgVOH25AicLy2WsjoiIOjuGJfJ6fSIDMXl4PB64KQ539O8GwPlKFNfg7++PXJSzPCIi6uQYlsjrCYKAJfcNxjsPDoFCIUjbx13vDE6bjrArjoiI2g/DEvmsO2tbmfadK0WhuVrmaoiIqLNiWCKfFRnshxvi9QCArTlF8hZDRESdFsMS+bTRfSMAANtOXZa5EiIi6qwYlsinjepTG5ZOXoYoijJXQ0REnRHDEvm0od1DoVMrcbm8BjmXyuQuh4iIOiGGJfJpGpUCI2pn+P75BLviiIio7TEskc+79YquOCIiorbGsEQ+zzVuaUduMWpsjmaPv1Bahaoae3uXRUREnQTDEvm867oFISJQg8oaO/adK2nwGFEUseHQRdz/fhaS3szE7z/d18FVEhGRr2JYIp+nUAgY3c8AAPh01/kGj1m26TjmrtqDPWedYernE3x6joiIPMOwRJ1CalJPAMD6AxdwobTKbd+/fzqN9zJPAgBm39oLSoWAiho7CsssHV0mERH5IIYl6hQGx+lxS0IYbA4R6dtype2f78nD/31zFADwXMp1WPTrAege5g8AOFVYLkutRETkWxiWqNN4bHQCAGD1zvMoNFdj4+GLeP6LgwCAOcm98MSvegMAehsCAACnihiWiIioeSq5CyBqK7/qF4m+kYE4UViO4W9sgUIAHCLw0LA4/Omu/hAEAQDQ2xCIzUcLcaqoQuaKiYjIF7BliToNhULA0gcGo09kIABnUBp/fRTeuHeQFJQAIIEtS0RE1AJsWaJO5cbuodj8zG0oqajBGWMFBsfpoVQIbsf0NjjD1Gm2LBERkQcYlqhTCg3QIDRA0+A+V1jKL61CZY0N/hr+MyAiosaxG466nNAADcJqg1TuZbYuERFR0xiWqEuqeyKuLiz9eLwIGw4VyFUSERF5KYYl6pISIpxdca65lqqtdvy//+zGE//di0vmajlLIyIiL+MTYWnr1q0QBKHBZdeuXfWOP3nyJIKCgqDX65u99q5duzB27Fjo9XqEhoYiJSUFBw4caIe7IG/SO9L9ibjTRRWotjrgEIEjF8xylkZERF7GJ8JSUlISCgoK3JbZs2ejV69eGDZsmNuxVqsVkydPRnJycrPXLS8vx/jx49G9e3fs2LEDP//8M4KCgpCSkgKr1dpet0Ne4Oon4o5fKpP2Hb3IsERERHV84jEgjUaDqKgoad1qtSIjIwPz5893mz8HABYtWoTExESMHTsWWVlZTV732LFjKC4uxmuvvYb4+HgAwMsvv4zBgwfj7Nmz6NOnT9vfDHmFvpFBAICTheWw2OxuYelYQVljpxERURfkEy1LV1u/fj2MRiNmzJjhtj0zMxNr1qzB8uXLPbrOddddh/DwcKSlpaGmpgZVVVVIS0tD//790bNnz0bPs1gsMJvNbgv5lvgwHcICNKixO3C0oAzHL9VNUJlzkWGJiIjq+GRYSktLQ0pKCuLi4qRtRqMRqampWLlyJYKDgz26TlBQELZu3YpVq1ZBp9MhMDAQGzZswHfffQeVqvFGtyVLliAkJERaXK1S5DsEQcDQeD0AYN+5EreWpVNFztYmIiIiQOawtGDBgkYHbruWY8eOuZ2Tl5eHjRs3YtasWW7b58yZgylTpmD06NEef/+qqirMmjULo0aNwvbt27Ft2zYMHDgQd999N6qqqho9b+HChTCZTNJy/vz5lt04eYWh3fUAgKxTRpwvqQQAaFQK2BwiThVy/iUiInISRFEU5frmRUVFMBqNTR6TkJAAjaZuJubFixfjvffeQ35+PtRqtbRdr9ejvLyuK0UURTgcDiiVSqxYsQIzZ86sd+20tDT86U9/QkFBARQKZ26sqalBaGgo0tLS8Mgjj3h0H2azGSEhITCZTB63apH8fj5xGY+m7YBKIcDmEBEWoEGfyEDszC3GsoeG4L4b45q/CBER+SxPf3/LOsDbYDDAYDB4fLwoikhPT8e0adPcghIAZGdnw26v6zrJyMjA0qVLkZWVhdjY2AavV1lZCYVC4TZI3LXucDhaeDfkawbHh0AQAJvD+f8L/boFol+3IOzMLea4JSIikvjUmKXMzEzk5uZi9uzZ9fb1798fAwcOlJbY2FgoFAoMHDgQoaGhAIB169YhMTFROufOO+9ESUkJ5s2bh6NHj+Lw4cOYMWMGVCoVbr/99g67L5JHsJ8afSMDpfV+3YKQGOX8P4ujDEtERFTLp8JSWloakpKS3AJPS5hMJuTk5EjriYmJ+Oqrr3Dw4EGMHDkSycnJuHDhAjZs2IDo6Oi2Kpu82ND4UOnrft2CkBjtnFLgWAGfcCQiIidZxyx1Fhyz5Ls+2XkOC9b+AgD47P+NxICYYAx8eSMAYM+iOxAeqJWzPCIiakee/v72qZYlorY2tPuVLUuBCNSq0Ke2ay7zWKFcZRERkRdhWKIurV+3QEwf2QPzx/SB3t/51OU9Q2IAAF/szZOzNCIi8hIMS9SlCYKAV+8ZiGfHXSdtu/dG59OT208X43xxpVylERGRl2BYIrpKXKg/knqHAwDW7cuXuRoiIpIbwxJRA+6vnZDyi715qLFxzi0ioq6MYYmoAeMHRsFfo8RZYyUGv7oR0z/YicvlFrnLIiIiGTAsETUgQKvCH1OuQ6i/GtVWB348XoT/bj8nd1lERCQDhiWiRqSO6oW9L96JxfdcDwDIOJAPTktGRNT1MCwRNUEQBNx7Yxz81AqcLqrAoXzO7E1E1NUwLBE1I1Crwh39uwEAvtzPp+OIiLoahiUiD9xzg3Pupa8OXIDdwa44IqKuhGGJyAO39TNA769GYZkF208b5S6HiIg6EMMSkQc0KgVSBkQBALbm8J1xRERdCcMSkYeS+jhn9d6ZWyxzJURE1JEYlog8dHPPMADAoQtmlFtsMldDREQdhWGJyEMxeh3iw3SwO0TsPVsidzlERNRBGJaIWmB4T2dX3I5cDvImIuoqGJaIWmBEL2dXHMctERF1HQxLRC0wvDYsHThvQrXVLnM1RETUERiWiFqgR7g/ugVrUWN3YP/50kaPO5RvwkVTdccVRkRE7YZhiagFBEHA8F7OcUs/Hi9q8JhTReW4Z/k2/DZtR70X7+7MLcZ3vxS0e51ERNR2GJaIWujuQc7JKf+7/WyDUwjsOF0Mu0PEicJynCoql7ZbbHbMWrkLT3y8F+eLKzusXiIiujYMS0QtNG5AFBIMATBX2/Df7Wfr7f8lv1T6emtOXevTnrMlKLPYIIpA7uWKjiiViIjaAMMSUQspFAIev603AOBfP+Xib5knMO7PP+Kz3ecBAAfzTNKxP1zxapT/Hb8sfZ1fWtVB1RIR0bViWCJqhUlDYxGr1+FyuQXvfH8cxy+V4//bfALVVjtyLpZJx+3MLZa66n46UdfKlFfCbjgiIl/BsETUCmqlAr+/oy8AoFdEALQqBfJLq7B2bz5sDhHhARr0CPeH1S4i6+RlXC634PAFs3R+XglbloiIfIVK7gKIfNWDw+Jxa98IGAK1mLtqLzYfvYS/ZZ4AAAyOC0GP8ACszDqDH3KKUHXVnEwMS0REvoMtS0TXIDpEB5VSgTsHRAIALtTOrTQoTo/brjMAAL4+eAH/+uk0AGBkgnPaAXbDERH5DoYlojYwJrEbBKFufUhcCJJ6h+P6mGCUVdtwKN/ZBTd5RHcAQGGZBRYbZwAnIvIFDEtEbcAQpMWN3UOl9UGxIdCqlFj7RBIWTkhEoFaFXhEBGDegG3RqJUQRKCjlDN9ERL6AY5aI2sidA7phz9kSRAX7ITLYDwCgVSnx/27rjRmjesEhivBTKxEbqsPJwnLklVShZ0SAzFUTEVFz2LJE1EbuvzEON/UIxWOjE+rt06gU8FMrAQBxoToAHLdEROQr2LJE1EYMQVp88XhSs8e5whInpiQi8g1sWSLqYHGh/gCc0wfsOG3EWxuONfiOOSIi8g4+EZa2bt0KQRAaXHbt2lXv+JMnTyIoKAh6vb7Za2/ZsgVJSUkICgpCVFQUnn/+edhs/MVF7cfVspRzsQxP/Hcv/r71FOb9dy9sdofMlRERUUN8IiwlJSWhoKDAbZk9ezZ69eqFYcOGuR1rtVoxefJkJCcnN3vdAwcO4K677sL48eOxb98+fPrpp1i/fj0WLFjQXrdChFi9MywdKTDDWFEDAPjxeBFe+eowRFGUszQiImqAT4QljUaDqKgoaQkPD0dGRgZmzJgB4crJbQAsWrQIiYmJeOihh5q97qefforBgwfjpZdeQp8+fXDbbbfhrbfewvLly1FWVtbs+USt4eqGc3nwpjgIArBq+zms3nlepqqIiKgxPhGWrrZ+/XoYjUbMmDHDbXtmZibWrFmD5cuXe3Qdi8UCPz8/t206nQ7V1dXYs2dPm9VLdKWIQA20Kuc/vQHRwVh6/2A8l3IdAOC1rw/jxCUGdSIib+KTYSktLQ0pKSmIi4uTthmNRqSmpmLlypUIDg726DopKSnIysrC6tWrYbfbkZ+fj9deew0AUFBQ0Oh5FosFZrPZbSHylCAIuD7G+Xf0+QmJUCgEzB3dG8l9I1BtdWD+6n2otnJ2byIibyFrWFqwYEGjA7ddy7Fjx9zOycvLw8aNGzFr1iy37XPmzMGUKVMwevRoj7//uHHj8Pbbb2Pu3LnQarXo168f7rrrLgCAQtH4H82SJUsQEhIiLfHx8S24ayLg71NvwpfzRuG2fs73xykUAt59aAjCAzQ4drEMf91yQuYKiYjIRRBlHFFaVFQEo9HY5DEJCQnQaDTS+uLFi/Hee+8hPz8farVa2q7X61FeXi6ti6IIh8MBpVKJFStWYObMmY1+D1EUUVBQgNDQUJw5cwYDBgzAzp07cfPNNzd4vMVigcVikdbNZjPi4+NhMpk8btUiasjGwxfx//6zB2qlgA2/H43ehkC5SyIi6rTMZjNCQkKa/f0t66SUBoMBBoPB4+NFUUR6ejqmTZvmFpQAIDs7G3Z7XddFRkYGli5diqysLMTGxjZ5XUEQEBMTAwBYvXo14uPjceONNzZ6vFarhVar9bhuIk+NG9ANYxIjkXmsEC9nHMYz4/ph+2kj7h4UjR7hfDUKEZEcfGoG78zMTOTm5mL27Nn19vXv399tfffu3VAoFBg4cKC0bd26dVi4cKFb197bb7+N8ePHQ6FQYO3atXjzzTfx2WefQalUtt+NEDVCEAS8/JsB+PnkZWkBgNU7z+Hr+ckI0ambuQIREbU1nxrgnZaWhqSkJCQmJrbqfJPJhJycHLdt3333HZKTkzFs2DB88803yMjIwKRJk9qgWqLW6REegKfG9AEABGiUCAvQ4HxxFZ5bc8BtHia7Q+S8TEREHUDWMUudhad9nkSeEkURJwrL0T3MH8cvleGB97NRY3dg6ojumD+mL7Ycu4Rl3x9HWIAGL/56AEb387w7m4iInDz9/c2w1AYYlqi9/Sf7DF7MONzo/nuHxuLdB4dAoRAaPYaIiNx5+vvbp7rhiLqq347sibTpw3BTj1AAQJCfCi/+egBm3doLSoWAdfvy8fmePJmrJCLqnHxqgDdRVza2fzeM7d8NORfL0C1YC72/c0qNqGA/vP7tUbz+7VGM6R+JiEA+qUlE1JbYskTkY66LCpKCEgDMGNUT/aODYaqy4vVvjspYGRFR58SwROTjVEoFltw3CIIArNuXj093nZO7JCKiToVhiagTuCFej9+P7QcAeGHdIWSdvAxRFHGysBxLvj2K2R/uRu7lCpmrJCLyTRyzRNRJPDW2D04VlWP9gQuY8u8d9fb/kl+KTx4biV4RnAmciKgl2LJE1EkIgoC3HhjsNueSUiHgjv6R6BsZiEtmCyav2I68kkoZqyQi8j2cZ6kNcJ4l8jallTWosTugUysR5KdGUZkFk/+1HScLyzGiVxhWz7mFczIRUZfHeZaIujC9vwaRQX4I8nO+S84QpMUH02+GTq3EjtxirNpxVuYKiYh8B8MSURfRPdwfCyY436v45nfHcNbY+gHfm45cwjsbc2CzO9qqPCIir8WwRNSF/PaWHhjeKwyVNXZM+dcOnCoqR+7lCqz43ylsq32Crjl7zpbg8VV78LcfTuKrgxc6oGoiInlxzFIb4Jgl8iUFpipM/fcOnC6qQIBGiUqrHa6fAiMTwjE7uReSekdAp1FK59TYHLDaHbDYHLj7rz+hwFQNAEjqHY6P59wix20QEV0zvki3AzEska8xlluQmr4Lv+SbAADDeoTiYJ4JNbXdalqVAv26BaFbsBaXy2tw5IIZNXYHFALgEIFYvQ4XTFUQReCnP96O+DB/OW+HiKhVGJY6EMMS+aJyiw1fH7iAm3qEom+3IOSVVOJf/zuNzUcLkV9a1eh5ARolPn88Cf/3zRFsO2nE78b2xdN39uvAyomI2gbDUgdiWKLORBRFnCoqx5nLlbhUVo1ArQo3xOsREahFcUUNArUqhAZokLE/H7/7ZD9i9TrMSe4FY0UNJgyMxoAY/hsgIt/AsNSBGJaoK6q22nHz65tRVm1z2z6stqUqPECD7mH+6BHuD6tdhKnKiohADXpFBCCvtAp7zpSgwFSNcosVPcID8NuRPRBcO9UBEVFHYFjqQAxL1FV98HMuPsw+g96GQKiVAjYfLYTd0bofKSE6NR68KQ7dw/0RGeSHbsFaRIX4ISJQC7VSAZvdgRq7A/4avqWJiNoGw1IHYlgiciowVWHz0UIYyy24XG7BmcuVOFdcCT+1AsF+ahSWWXC+pBJ6nRrDeoahtyEQOrUSXx28gJOF5Q1eUxCAQI0KZRZnC9ag2BAk942AUiGgrNoGc5UV5RYbdBol9Do1okJ06BHuj/AADfw1Kug0Cug0Kuh1agRory1onblcgY+yzyI+TIcpI7pDq1I2fxIReS2GpQ7EsETkuRqbA2qlAEGoe92K3SHim18KsOdMMS6ZLbhUVo1LpmoUlllga2VLVUOC/FQI9dfA7hDhEMXa/wIOUYRCAAxBfogO8UNUiB/CAzTIL6nCGWMFdBolNEoF/nfistRyFqvX4eGb49Ej3B/BfmrYHCL81ArE6HW4aKrG53vykF9ShXHXd8P9N8YhNEDjVovFZkd1jQPBOpX0Z2GzO6AQBL6KhqiDMCx1IIYlovbhcIgorqyBucqKEJ0adoeIH48XYc/ZEmhUCgT5qRDk52wxsljtKKmsQX5JFc4WV8JUZUVVjR1VVjsqa+yosbXNbOO39onAicIyXDJbWnSe3l8NvU6Niho7zFVWWGrriQr2w829wlBQWoUDeaVwiECovxr9o4Nx54BuiAzS4lRRBaqtdhiCtIgO0SExKghxoTpYaue/CtSq3MInEXmGYakDMSwReb9yiw0FpVUwV9ugVAhQCgIEAc6vFQJsdhGXyqpRUFqNi6YqXK6oQUyIH3pFBMJis6Ok0ophPUIxJF6Pqho7Ptl1Dr/km5BfUoXKGjuUCgGVNTbkl1RBqRDw6yEx6BsZiM/35OHwBXO73ptaKSDUX4OwAA2Ca0NltdVeuzicrVhWB4L8VIgM0sIhAqYqKxSCc6yYRqWA1S7C5nDAZhdhc4hSK5chSIvwQA2UCgUEOLtFFYJQ+7Xzz1BA7TYBtYvgvs3tWAEKoe46uGqb62sItdsg1B4LKRC6rnvlfkGq6cpanNuUCgFqpQJqpVDbogjYHQ7YHYBKIUCtcu1XQKNSQKNUwCGKsFgdta2OAhQK13Xrvo/0vQHpPty+f+1xwFV/Dh5+rgpBgE6jgFalhEN0fi4Ohwi7KEIUIU0mK6L5X+ONfdeGMnZjubtF12isjkZ2iCKkVl5H7Y0JcP65u/5ce4T7t3nXN8NSB2JYIiIX14/UK1t6iitqUFRmganKigCtEsF+agT7qaFWCdh3rhR7z5agW4gfbu4ZBn+NEkVlFmSduowtRwtRbbWjtyEQAVoVisosOFtciZOFZbDa+aObupYtz96G3obANr2mp7+/+VgJEVEbaqg7LCzA2erTkFF9IjCqT4Tbtm7BfhgYG4LHRvdu8JwamwOllTXQaZRQKRQoqaxBcYVzMVVZoVYq4Kd2tkj4qRXSmCtztQ2F5mqolAJCdGo4RKC00gqb3QGVUgGVQoBKKUClcLbC2BwiCsssKKmoqf0/fmcYFEVna4arNcD1tWtfvW1wtRw4v3Y4L+D82uG+H7X7Xd/DIbWgiHA43Ldd/T0ars/ZkmS1O1Bjc0CldLYkqhQCFIIAu0NETe0+a+0TlzU2Z6uaVqWAQiFI13WIrhrqQrFbHXC/X7ev4d4a5Am7Q0SV1Q6LzQ6lIEitoEpF3Zg/19+2pnphG/ueDW1uvL6GdzR0fGOXaKxtRoSzFU0huLe8Xf3nqJJxLB/DEhGRj9GoFIgM9pPWdRodYvQ6GSsi6twUchdARERE5M0YloiIiIiawLBERERE1ASGJSIiIqImMCwRERERNYFhiYiIiKgJDEtERERETWBYIiIiImoCwxIRERFRExiWiIiIiJrgE2Fp69attW+srr/s2rULAHDmzJkG92/fvr3Ja587dw533303/P39ERkZieeeew42m60jbouIiIh8gE+8Gy4pKQkFBQVu21588UVs2bIFw4YNc9u+efNmXH/99dJ6eHh4o9e12+24++67ERUVhaysLBQUFGDatGlQq9V444032vYmiIiIyCf5RFjSaDSIioqS1q1WKzIyMjB//vx6b/gODw93O7Yp33//PY4cOYLNmzejW7duuOGGG7B48WI8//zzeOWVV6DRNPyWcCIiIuo6fKIb7mrr16+H0WjEjBkz6u2bOHEiIiMjceutt2L9+vVNXic7OxuDBg1Ct27dpG0pKSkwm804fPhwo+dZLBaYzWa3hYiIiDonn2hZulpaWhpSUlIQFxcnbQsMDMS7776LUaNGQaFQ4IsvvsCkSZPw5ZdfYuLEiQ1e5+LFi25BCYC0fvHixUa//5IlS/Dqq6/W287QRERE5Dtcv7dFUWz6QFFGzz//vAigyeXo0aNu55w/f15UKBTi559/3uz1f/vb34q33npro/vnzJkjjhs3zm1bRUWFCED89ttvGz2vurpaNJlM0nLkyJFm74MLFy5cuHDh4p3L+fPnm8wTsrYsPfvss0hNTW3ymISEBLf19PR0hIeHN9padKURI0Zg06ZNje6PiorCzp073bZdunRJ2tcYrVYLrVYrrQcGBuL8+fMICgqqN4bqWpjNZsTHx+P8+fMIDg5us+tS2+Fn5N34+Xg3fj7er7N/RqIooqysDDExMU0eJ2tYMhgMMBgMHh8viiLS09OlJ9aas3//fkRHRze6f+TIkXj99ddRWFiIyMhIAMCmTZsQHByMAQMGeFyXQqFw6xJsa8HBwZ3yL2lnws/Iu/Hz8W78fLxfZ/6MQkJCmj3Gp8YsZWZmIjc3F7Nnz66378MPP4RGo8HQoUMBAGvXrsUHH3yAf//739Ix69atw8KFC3Hs2DEAwLhx4zBgwAD89re/xVtvvYWLFy9i0aJFmDdvnlvLEREREXVdPhWW0tLSkJSUhMTExAb3L168GGfPnoVKpUJiYiI+/fRTPPDAA9J+k8mEnJwcaV2pVOLrr7/G448/jpEjRyIgIADTp0/Ha6+91u73QkRERL5BEMXmhoCTXCwWC5YsWYKFCxeypctL8TPybvx8vBs/H+/Hz8iJYYmIiIioCT45KSURERFRR2FYIiIiImoCwxIRERFRExiWiIiIiJrAsOTFli9fjp49e8LPzw8jRoyoN9s4dYxXXnkFgiC4LVdOX1FdXY158+YhPDwcgYGBuP/++6WZ4Knt/e9//8NvfvMbxMTEQBAEfPnll277RVHESy+9hOjoaOh0Otxxxx04ceKE2zHFxcWYOnUqgoODodfrMWvWLJSXl3fgXXRuzX1Gqamp9f5NjR8/3u0YfkbtZ8mSJbj55psRFBSEyMhITJo0yW1aHcCzn2vnzp3D3XffDX9/f0RGRuK5556DzWbryFvpMAxLXurTTz/FM888g5dffhl79+7FkCFDkJKSgsLCQrlL65Kuv/56FBQUSMvPP/8s7Xv66afx1VdfYc2aNfjxxx9x4cIF3HfffTJW27lVVFRgyJAhWL58eYP733rrLfz1r3/FP/7xD+zYsQMBAQFISUlBdXW1dMzUqVNx+PBhbNq0CV9//TX+97//4bHHHuuoW+j0mvuMAGD8+PFu/6ZWr17ttp+fUfv58ccfMW/ePGzfvh2bNm2C1WrFuHHjUFFRIR3T3M81u92Ou+++GzU1NcjKysKHH36IlStX4qWXXpLjltpfs2+jJVkMHz5cnDdvnrRut9vFmJgYccmSJTJW1TW9/PLL4pAhQxrcV1paKqrVanHNmjXStqNHj4oAxOzs7A6qsOsCIK5bt05adzgcYlRUlPj2229L20pLS0WtViuuXr1aFEVRevH1rl27pGO+++47URAEMT8/v8Nq7yqu/oxEURSnT58u3nPPPY2ew8+oYxUWFooAxB9//FEURc9+rn377beiQqEQL168KB3z/vvvi8HBwaLFYunYG+gAbFnyQjU1NdizZw/uuOMOaZtCocAdd9yB7OxsGSvruk6cOIGYmBgkJCRg6tSpOHfuHABgz549sFqtbp9VYmIiunfvzs9KBrm5ubh48aLb5xESEoIRI0ZIn0d2djb0ej2GDRsmHXPHHXdAoVBgx44dHV5zV7V161ZERkbiuuuuw+OPPw6j0Sjt42fUsUwmEwAgLCwMgGc/17KzszFo0CB069ZNOiYlJQVmsxmHDx/uwOo7BsOSF7p8+TLsdrvbX0IA6NatGy5evChTVV3XiBEjsHLlSmzYsAHvv/8+cnNzkZycjLKyMly8eBEajQZ6vd7tHH5W8nD9mTf1b+fixYvSi7NdVCoVwsLC+Jl1kPHjx+Ojjz7Cli1bsHTpUvz444+YMGEC7HY7AH5GHcnhcOD3v/89Ro0ahYEDBwKARz/XLl682OC/M9e+zsan3g1HJIcJEyZIXw8ePBgjRoxAjx498Nlnn0Gn08lYGZFveuSRR6SvBw0ahMGDB6N3797YunUrxo4dK2NlXc+8efNw6NAht3GYVB9blrxQREQElEplvScPLl26hKioKJmqIhe9Xo9+/frh5MmTiIqKQk1NDUpLS92O4WclD9efeVP/dqKiouo9KGGz2VBcXMzPTCYJCQmIiIjAyZMnAfAz6ihPPvkkvv76a/zwww+Ii4uTtnvycy0qKqrBf2eufZ0Nw5IX0mg0uOmmm7BlyxZpm8PhwJYtWzBy5EgZKyMAKC8vx6lTpxAdHY2bbroJarXa7bPKycnBuXPn+FnJoFevXoiKinL7PMxmM3bs2CF9HiNHjkRpaSn27NkjHZOZmQmHw4ERI0Z0eM0E5OXlwWg0Ijo6GgA/o/YmiiKefPJJrFu3DpmZmejVq5fbfk9+ro0cORK//PKLW6jdtGkTgoODMWDAgI65kY4k9whzatgnn3wiarVaceXKleKRI0fExx57TNTr9W5PHlDHePbZZ8WtW7eKubm54rZt28Q77rhDjIiIEAsLC0VRFMW5c+eK3bt3FzMzM8Xdu3eLI0eOFEeOHClz1Z1XWVmZuG/fPnHfvn0iAHHZsmXivn37xLNnz4qiKIpvvvmmqNfrxYyMDPHgwYPiPffcI/bq1UusqqqSrjF+/Hhx6NCh4o4dO8Sff/5Z7Nu3rzh58mS5bqnTaeozKisrE//whz+I2dnZYm5urrh582bxxhtvFPv27StWV1dL1+Bn1H4ef/xxMSQkRNy6datYUFAgLZWVldIxzf1cs9ls4sCBA8Vx48aJ+/fvFzds2CAaDAZx4cKFctxSu2NY8mLvvfee2L17d1Gj0YjDhw8Xt2/fLndJXdLDDz8sRkdHixqNRoyNjRUffvhh8eTJk9L+qqoq8YknnhBDQ0NFf39/8d577xULCgpkrLhz++GHH0QA9Zbp06eLouicPuDFF18Uu3XrJmq1WnHs2LFiTk6O2zWMRqM4efJkMTAwUAwODhZnzJghlpWVyXA3nVNTn1FlZaU4btw40WAwiGq1WuzRo4c4Z86cev8jyM+o/TT02QAQ09PTpWM8+bl25swZccKECaJOpxMjIiLEZ599VrRarR18Nx1DEEVR7OjWLCIiIiJfwTFLRERERE1gWCIiIiJqAsMSERERURMYloiIiIiawLBERERE1ASGJSIiIqImMCwRERERNYFhiYioDQiCgC+//FLuMoioHTAsEZHPS01NhSAI9Zbx48fLXRoRdQIquQsgImoL48ePR3p6uts2rVYrUzVE1JmwZYmIOgWtVouoqCi3JTQ0FICzi+z999/HhAkToNPpkJCQgM8//9zt/F9++QVjxoyBTqdDeHg4HnvsMZSXl7sd88EHH+D666+HVqtFdHQ0nnzySbf9ly9fxr333gt/f3/07dsX69evl/aVlJRg6tSpMBgM0Ol06Nu3b71wR0TeiWGJiLqEF198Effffz8OHDiAqVOn4pFHHsHRo0cBABUVFUhJSUFoaCh27dqFNWvWYPPmzW5h6P3338e8efPw2GOP4ZdffsH69evRp08ft+/x6quv4qGHHsLBgwdx1113YerUqSguLpa+/5EjR/Ddd9/h6NGjeP/99xEREdFxfwBE1Hpyv8mXiOhaTZ8+XVQqlWJAQIDb8vrrr4ui6HzL+ty5c93OGTFihPj444+LoiiKK1asEENDQ8Xy8nJp/zfffCMqFArx4sWLoiiKYkxMjPjCCy80WgMAcdGiRdJ6eXm5CED87rvvRFEUxd/85jfijBkz2uaGiahDccwSEXUKt99+O95//323bWFhYdLXI0eOdNs3cuRI7N+/HwBw9OhRDBkyBAEBAdL+UaNGweFwICcnB4Ig4MKFCxg7dmyTNQwePFj6OiAgAMHBwSgsLAQAPP7447j//vuxd+9ejBs3DpMmTUJSUlKr7pWIOhbDEhF1CgEBAfW6xdqKTqfz6Di1Wu22LggCHA4HAGDChAk4e/Ysvv32W2zatAljx47FvHnz8M4777R5vUTUtjhmiYi6hO3bt9db79+/PwCgf//+OHDgACoqKqT927Ztg0KhwHXXXYegoCD07NkTW7ZsuaYaDAYDpk+fjlWrVuEvf/kLVqxYcU3XI6KOwZYlIuoULBYLLl686LZNpVJJg6jXrFmDYcOG4dZbb8V///tf7Ny5E2lpaQCAqVOn4uWXX8b06dPxyiuvoKioCPPnz8dvf/tbdOvWDQDwyiuvYO7cuYiMjMSECRNQVlaGbdu2Yf78+R7V99JLL+Gmm27C9ddfD4vFgq+//loKa0Tk3RiWiKhT2LBhA6Kjo922XXfddTh27BgA55Nqn3zyCZ544glER0dj9erVGDBgAADA398fGzduxO9+9zvcfPPN8Pf3x/33349ly5ZJ15o+fTqqq6vx5z//GX/4wx8QERGBBx54wOP6NBoNFi5ciDNnzkCn0yE5ORmffPJJG9w5EbU3QRRFUe4iiIjakyAIWLduHSZNmiR3KUTkgzhmiYiIiKgJDEtERERETeCYJSLq9DjagIiuBVuWiIiIiJrAsERERETUBIYlIiIioiYwLBERERE1gWGJiIiIqAkMS0RERERNYFgiIiIiagLDEhEREVETGJaIiIiImvD/A8QP4vEvm/rIAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "cudaq.set_target(\"nvidia\")\n",
    "\n",
    "start_time = timeit.default_timer()\n",
    "result = minimize(cost,\n",
    "                  x0,\n",
    "                  method='COBYLA',\n",
    "                  callback=callback,\n",
    "                  options={'maxiter': 500})\n",
    "end_time = timeit.default_timer()\n",
    "\n",
    "print('UCCSD-VQE energy =  ', result.fun)\n",
    "print('Total number of qubits = ', qubit_count)\n",
    "print('Total number of parameters = ', parameter_count)\n",
    "print('Total number of terms in the spin hamiltonian = ',\n",
    "      spin_ham.term_count)\n",
    "print('Total elapsed time (s) = ', end_time - start_time)\n",
    "\n",
    "plt.plot(exp_vals)\n",
    "plt.xlabel('Epochs')\n",
    "plt.ylabel('Energy')\n",
    "plt.title('VQE')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Gate Fusion for Larger Circuits"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "CUDA-Q simulations take advantage of a technique called gate fusion. Gate fusion is an optimization technique where consecutive gates are combined into a single gate operation to improve the efficiency of the simulation (See figure below).  By targeting the `nvidia` or `nvidia-mgpu` backend and setting the `CUDAQ_FUSION_MAX_QUBITS` or `CUDAQ_MGPU_FUSE` environment variable respectively, you can select the degree of fusion that takes place, see [link](https://nvidia.github.io/cuda-quantum/latest/using/backends/sims/svsims.html#single-gpu) for more details. \n",
    "\n",
    "![gate-fuse.png](./images/gate-fuse.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is particularly important for larger circuits and can have a significant impact on the performance of the simulation. Each system is different, so you should test different gate fusion levels to find out what is best for your system. You can find more information [here](https://developer.nvidia.com/blog/new-nvidia-cuda-q-features-boost-quantum-application-performance/)."
   ]
  }
 ],
 "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.10.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}