From e2659766868de4ce53a4dd966aaede63505d33c3 Mon Sep 17 00:00:00 2001 From: Julien Gacon Date: Fri, 7 May 2021 15:33:43 +0200 Subject: [PATCH] Add a tutorial for the VQE script and ``VQEProgram`` (#11) * add tutorial * update the VQE tutorital texts replace runtime for Qiskit Runtime. and add some context to differientate Qiskit Runtime programs, and Qiskit Runtime service * rerun notebook once to get final cell executed * add link to qiskit nature tutorials Co-authored-by: Ismael Faro --- tutorials/vqe.ipynb | 696 ++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 696 insertions(+) create mode 100644 tutorials/vqe.ipynb diff --git a/tutorials/vqe.ipynb b/tutorials/vqe.ipynb new file mode 100644 index 000000000..04b749444 --- /dev/null +++ b/tutorials/vqe.ipynb @@ -0,0 +1,696 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# VQE\n", + "\n", + "The Variational Quantum Eigensolver (VQE) is a central algorithm in many applications from e.g. quantum chemistry or optimization.\n", + "This tutorial shows you how to run the VQE as a Qiskit Runtime program. We'll start off by defining the algorithm settings, such as the Hamiltonian and ansatz, and then run a VQE both locally, on your machine, and remotely, using the Qiskit Runtime.\n", + "\n", + "**Note:** You can find tutorials on solving more comprehensive problems, such as finding the ground state of the lithium hydride molecule, using the VQE (and Qiskit Runtime) within [the tutorials of Qiskit Nature](https://github.com/Qiskit/qiskit-nature/tree/main/docs/tutorials)." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## System Hamiltonian" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Let's start by defining the operator of which we want to determine the ground state. Here we'll chose a simple diagonal Hamiltonian $\\hat H$ acting with Pauli-Z operators on the first two qubits\n", + "\n", + "$$\n", + "\\hat H = \\hat Z_0 \\otimes \\hat Z_1.\n", + "$$\n", + "\n", + "We can construct this Hamiltonian with Qiskit's `opflow` module:" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [], + "source": [ + "from qiskit.opflow import Z, I\n", + "\n", + "num_qubits = 4\n", + "hamiltonian = (Z ^ Z) ^ (I ^ (num_qubits - 2))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "This Hamiltonian has a ground state energy of -1." + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "target_energy = -1" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Parameterized Ansatz Circuit" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Next, we choose a parameterized quantum circuit $\\hat U(\\theta)$ to prepare the ansatz wavefunction\n", + "\n", + "$$\n", + "|\\psi(\\theta)\\rangle = \\hat U(\\theta)|0\\rangle.\n", + "$$\n", + "\n", + "We'll use the `EfficientSU2` circuit from Qiskit's circuit library, which is a hardware efficient, heuristic ansatz with alternating rotation and entanglement layers." + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAegAAADWCAYAAAAaVxFlAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAAzNUlEQVR4nO3dd0AUd/o/8PcuSxUVDfYepQiKIgJ2YklEzC+KLdH7egbN2Yk1JuqpQRK8GGNJYowenuWMJPYYKSY2YkNQiUYFiYIggiAq4lKEZff3h+cqbXaBWT6f2X1e/1yYmR3ecg/Pw85OkWk0Gg0IIYQQwhU56wCEEEIIqYgGNCGEEMIhGtCEEEIIh2hAE0IIIRyiAU0IIYRwiAY0IYQQwiEa0IQQQgiHaEATQgghHKIBTQghhHCIBjQhhBDCIRrQhBBCCIdoQBNCCCEcogFNCCGEcIgGNCGEEMIhGtCEEEIIh2hAE0IIIRyiAU0IIYRwiAY0IYQQwiEa0IQQQgiHFKwDEKKPuLg4ndts2bIFU6dOFdzG09NTrEiE1BrVNRFC76CJ0QgNDWUdgRDRUV2bLhrQhBBCCIdoQBNCCCEcogFNjMaOHTtYRyBEdFTXposGNCGEEMIhGtDEaEyaNIl1BEJER3VtuugyKx0uLNuGR9dTmHzvxq4d4B0cUKPXsspdm8yEED5Q36seQ/U9GtA6PLqegqzzN1jHqDap5iaEsCfV/iHV3FWhQ9zEaHzwwQesIxAiOqpr00UDmhgNXXdbIkSKqK5NFw1oYjT8/PxYRyBEdFTXposGNDEaOTk5rCMQIjqqa9NFJ4mJxHd/EJp4OEKtUkFTqoYyLRtXNhxA6pHzrKNVSYqZCSH8kGIPkVJmGtAiurJ+H66u3w+ZmRydJw+Dz3dzcPBaCp7euc86WpWkmLkqzs7OrCMQIjre61qKPUQqmekQtwFoStVI2nUMcnMFGru2Zx1HL1LMXN7OnTtZRyBEdFKpayn2EN4z04A2ALm5Ak6ThgIA8pIzGafRjxQzlxcSEsI6AiGik0pdS7GH8J6ZDnGLyO3DUegy/R0obK2gKSnF2fnf4XFCKgBgwMY5SD54BunHLgEABm1bhMTtR5ERfYVlZMHMDuMHoeMYH+22tu2aIutCIk7P2sAqrqBDhw5hyZIlrGNUi0ajgUwmYx2DcIz3uqa+Zzhcv4NWq9VYs2YNHBwcYGVlhW7duiE6OhpOTk5cXht49esD2O08CT+6Tkb6iXg079tFuy52+Ta4L3oPChsrtPXzRnFeAfMiBYQz/xV2AlGjVyBq9ApET18HVcEzxK/azTCtcUhMycWM4LOw67MTZt3/g1ZDwvDpd5fx4FEh62iEVBv1PcPhekBPmTIFwcHBmDZtGiIjIzFu3DiMHz8eycnJ8PDwYB2vSsVP8nF2wSa0HtwDbYZ6AgCKHuYhITQc3p9NRre5oxG7YjvbkOVUlllLJkP/jXNwOeQHKNMfsAloJI6eTUf3sYeweV8inihLoNEAGdkFCPo+Hu7jDiE5PY91REJqhPqe+Lgd0GFhYdi+fTsOHz6MhQsXYuDAgVi6dCl69+4NlUqFHj16sI4oqDhXietbjqDH4gnA/w5h3tpzCg06tkBCaASKc5WME1ZUWWYA6L5gLHITU5EWFccwnW7h4eGsIwh68KgQ/vOOoaSkFBpNxfWZOYXwn3MMmspWEpPFe12/ivqeuLgd0CEhIfD19YWPj0+Z5Z06dYK5uTnc3NwAAHfu3IGPjw8cHR3RtWtXnD59mkXcSiX8Oxw2Te3QaezLf8PTlPvI4+xU/leVz9yiX1e09OmGi8G7GCfTLTExkXUEQVsPJqGwqBTqKuavWq3B1b8e48zlrLoNRrjGe12XR31PPDINh3+up6eno02bNti6dSsmT55cZt348eORmJiI+Ph4AMDQoUMxYsQIzJw5E+fOncPYsWORkpICCwsLwe+h74k5HzceAGeLJjX7h1Si3/pZSNp9HNmxun/pEosf4ItHv9fo+4iR27qJHYbu/xTHJnyu9yGe2mQWos8DA0JDQ3VuFxoaKlak6uvwEVDPAZAJ/F2s0QAPIoCsg3WXizDDa11T3zNs39N37HJ5Fnd6ejoAoHnz5mWWFxYWIjo6GsOGDQPw/BZ4Z86cweHDhwEAffr0QcuWLXHy5EkMHTq0bkMbIbd5Y2BR3wb9NszWLnty+x7OL9rCMJWEyRTCwxkAoHm+HSGECZ76HpfvoG/dugUHBwesW7cOc+fO1S4PCgrCp59+io0bN2LmzJm4fPkyxowZg+TkZO0248aNw5AhQ0Q7yzty1HJmzxdt1tsFww6srNFrWeWuTWYhcXG6Pwfy8vJCbGys4Daenp6C6w3pH5+exn8OJlV5iPuF/6zsj4CRjnUTijDFa11T36seQ/U9Lv9Uf/311+Hm5oaQkBA0btwYrVq1wr59+xAREQEAXJ/BTdhZvHgx6wiCpo/rjNADSVWul8kAWxtzvDv09TpMRXjHe10Tw+HyJDG5XI69e/fC1dUVM2bMQEBAAOzt7TFr1iyYmZlpTxBr27YtsrKy8OzZM+1rU1JS0K5dO1bRCUP+/v6sIwjycLHH7PEula6TyZ5//Lzpn31gY83l382EEd7rmhgOlwMaABwdHXHy5Enk5+cjLS0NwcHB+PPPP+Hi4gJra2sAgL29Pfr27YutW7cCAM6dO4d79+5h4MCBLKMTRry8vFhH0GnDx72wak5PNGpQ9iTG11vXx4F1g/G34Z0YJSO8kkJdE8OQ1J/qFy9eRK9evcos+/777/H+++9j/fr1sLCwQFhYmM4zuAlhRS6X4ZMp3TBvYhdY9dwOADiz42306d6UbvlJCClDMgNaqVQiKSkJM2fOLLP89ddfx++/i39ZT3X0WDIBTT2dkR2XiLxbGega6I9zizYj6/wNuM54B219PaFMz8GZOd/CzMIcQ/csR96d+zg9+2tucwNAWz9veK8MwN6e06GwseImtzGwtDDT/ndf92YMkxBSc1X2kJgE9P8mEPVa2UNdrEL0jHUoLSph1kOqypkddxPDDq5Eo85tcXjIR3h65z5s2zRF/28CodFoUJDxEKcDv4FGrcbgnYth0dAGkSOW1Vlubg9xl2dra4vS0lIEBgayjlKGnVMbmNvaIMp/OSwb1YfCxgrXNh1G1vkbsHqtAVr07YLIEcvw+EYq2vp6QVVQhOjp61jHFsz9QvvhvZCf8RAAuMktpF+/fqwjECI6XutaqIc07tIe6mIVovyX49ZPJ/H6qP7MeohQTo2qFCcCViP1SIx2++K8fBybuApR/svx9G42Wg12BwAc//uqOs8umQHNq2benbU3f8+IvgqNWq1d91q3jrh/7joAIPP0VTTtyc+lM0K5AaDVIHdknP4TGo26spdzae3atawjECI6XutaqIcUZD6CzOz5eLFoUA/PHj9lkhHQ3euKcp6U+br4ST5KnhYAADQlpdCUsuuBNKBrycLOFu4Lx8F3fxC6zRkFCzvbl+sa1kOJ8vkTiorzCmDRoB6rmBUI5QaATuPeQPJ+th8dVNf8+fNZRyBEdLzWtVAPefboKcysLDDy9/VwmvQWUiMucJlTiHWzRmgxwI3p07ck8xk0r4pzlYj/8ifc/fUiWg/xQL2Wr2nXleQVoF6L51+b17dGcV4+q5gVCOVu3rcLHlxKgrpExTBh9Z05c4Z1BEJEx2tdC/WQlm90Q9HDPBwaMBfthvdCl+nv4Mq6fdzlrIrcQoH+G2bj3Eff0ztoKcu6kIBmvZ5f29q8jytk8pc/0pw/bqFZ7+frWvZ3w4NLVd+koq4J5W7k3BZt3uqJN3cvhZ1jG7h//B6rmIQQTgn1EADaJ1c9e/QU5g1s6jzfC7pyVqbPl9ORuD0KT5LSDR1PEA3oWsq9eRdqlQq++4OgVqmgKijSrit6mIesmAQM+zkYjV3bc/W4RqHcCVsjcHRsEH6b8Dlyk+4i/osfGSYlhPBIqIdkRF9BQ4dW8N0fhO6L3sXN7Ue5zAkAPpvno6VPN/TbMBtthnqiiYcj2vl5w+Ufb8N3fxDaDmN3HTod4hbB5ZDd2v9uN7wXugb6Iy8lE1nnb+DaxkO4tvGQdr3Cxgr9N36InD9uM0hallDuF15cUsBT7qroul8xIVLEc10L9ZBT//iqzLYse4hQzuhpFU/C+8FhYoVlg3cuRmFWriFjVkADWmSp4TFIDY+pcr2qoKhOr6PTl1Rzv+rgwYN0W0RidKRS11LpIbpyVoXFZVY0oHVo7NpBkt+bVW6WP69Vq1ZJopERUh0s6pr6Hh/flwa0Dt7BAawj1IhUcxNC2JNq/5Bq7qrQSWKEEEIIh2hAE6Px1Vdf6d6IEImhujZdNKCJ0XB2dmYdgRDRUV2bLhrQxGgMHz6cdQRCREd1bbpoQBNCCCEcogFNCCGEcIgGNDEa7u7urCMQIjqqa9NFA5oYjfj4eNYRCBEd1bXpogFNCCGEcIgGNCGEEMIhGtDEaOzYsYN1BEJER3VtumhAE0IIIRyiAU2MxqRJk1hHIER0VNemi55mpcOFZdvw6HoKk+/d2LVDjZ/Owip3bTITQvhAfa96DNX3aEDr8Oh6CrLO32Ado9qkmpsQwp5U+4dUc1eFDnETo/HBBx+wjkCI6KiuTRcNaGI0pk6dyjoCIaKjujZdNKCJ0fDz82MdQac8ZTH2/ZqCT9bHaZe9+9EJfLYlHr+eS0dJiZphOsIjKdQ1MQz6DFokvvuD0MTDEWqVCppSNZRp2biy4QBSj5xnHa1KUswsJCcnh3WEKt3LysdnW/7Af4/cQn6hqsy6PUdTsOfo8xNbWja1wfSxzljw966wsaZfT8J3XUuxh0gpM72DFtGV9fvwQ6eJCHMJwK09p+Dz3RzUb9+cdSxBUswsNTsP/wVX/wP4fm9iheFcXkZ2AZZvvIzu4w7i3B9ZdZSQkJqTYg+RSmYa0AagKVUjadcxyM0VaOzannUcvUgxc3nOzs6sI1Sw7NtLmPTP3/FEWVyt1/2VmoeBUyLw88lUAyUjUsFjXVdGij2E98w0oA1Abq6A06ShAIC85EzGafQjxczl7dy5k3WEMr7ZfR2fbfmjxq8vLlFj3MITOBtP76RNGW91XRUp9hDeM9OHXCJy+3AUukx/BwpbK2hKSnF2/nd4nPD8HdCAjXOQfPAM0o9dAgAM2rYIiduPIiP6CsvIgpkdxg9CxzE+2m1t2zVF1oVEnJ61gVVcQSEhIViyZAnrGACAxJRcfLQ2VnAbzdUpAACZ29YqtykuUeP9Zb/jyl5/+kzaRPFU15Whvmc4XL+DVqvVWLNmDRwcHGBlZYVu3bohOjoaTk5OXF56cPXrA9jtPAk/uk5G+ol4NO/bRbsudvk2uC96DwobK7T180ZxXgHzIgWEM/8VdgJRo1cgavQKRE9fB1XBM8Sv2s0wrbBDhw6xjqA154sYPCsW54zsW2l5WLPjT1H2VRtFz1TYdeQWPlkfh2XfXsKZy/eh0WhYxzJ6PNV1ZajvGQ7XA3rKlCkIDg7GtGnTEBkZiXHjxmH8+PFITk6Gh4cH63hVKn6Sj7MLNqH14B5oM9QTAFD0MA8JoeHw/mwyus0djdgV29mGLKeyzFoyGfpvnIPLIT9Amf6ATUAJuZmSi1/P3RN1n9/vTWR6CdbeX1PQYlAYJi6JxuptV/HZlj/Q//1w9Bh3CMnpecxyEX5Q3xMftwM6LCwM27dvx+HDh7Fw4UIMHDgQS5cuRe/evaFSqdCjRw/WEQUV5ypxfcsR9Fg8AZDJAAC39pxCg44tkBAageJcJeOEFVWWGQC6LxiL3MRUpEXFCbyavLAr/Lbo+8x8UIDjFzJE368+jkSn4d2PTmhPdHv1TfPVvx5jwPvhyH5YyCQb4Qv1PXFxO6BDQkLg6+sLHx+fMss7deoEc3NzuLm5AQCWL18OR0dHyOVy7Nu3j0XUKiX8Oxw2Te3QaezLf8PTlPvIu3OfYSph5TO36NcVLX264WLwLsbJdAsPD2cdAQAQe80wf20bar9CNBoNPvoqFjKUHcwvqNUa3MsuwMYfE+o8m6ngpa71RX1PPFyedZKeno5r165h3rx5FdalpaXB1dUVlpaWAABfX1+8//77mDx5cl3HLCNq9IoKy0qUhQhz4ffJTroyWzexg3fIFByb8DnUJcLX7/IgMTERTZo0YR0Df/712CD7vZr0yCD7FRL75wMk3nmic7st+xMRNIvvo1pSxUtdV4b6nmFxO6ABoHnzsheOFxYWIjo6GsOGDdMu69OnT42+h+yVQxlCPm48AM4WbH45oqOj4adnzvLEyO02bwws6tug34bZ2mVPbt/D+UVbqnxNbTIL0eeBAaGhoTq38/LyEitS1Vy+AcystV++OFu7KlWtL3929/6DRyBbN6T2+aqjgQfQbobOze7nFEImkwOgk8aqg9e6pr5n2L6n78mVXA5oe3t7AEBSUlKZ+9CuXr0amZmZXJ8gpsuZuRtZR9DbhSWhuLAklHUM6dEY6K9uQ+1XiFrPz5bVz0DDmVSF+l7NyDQcXiehVqvh7u6OzMxMrFmzBq1atcK+ffsQERGBtLQ0xMTEwNvbu8xr3njjDcyePRtjxowRNUvkqOXMni/arLcLhh1YWaPXsspdm8xC4uJ0n6jh5eWF2Fjha489PT0F14uhz8RfcP5Kts7t9LkO+lULJ3XFlwvq4AjAK54Vl6LVkDA8zH1W5TZyGRAw0hGhQf3rMJlx4LWuqe9Vj6H6Hpcnicnlcuzduxeurq6YMWMGAgICYG9vj1mzZsHMzEx7ghghr1q8eDHrCAAADxd7A+33NYPsV4ilhRk+er9rlevlMkChkGPu/7nWYSrTwktdk7rH5SFuAHB0dMTJkyfLLJs4cSJcXFxgbW1dxauIKfP392cdAQAwekh7fBsm7l/x1pZm8O3bWtR96uuj991w934+Nv6YALkMUL9yzM3C3Ax7vxqELg6NmWQzBbzUNal7XL6DrsrFixcrfP68bNkytG7dGufPn8e0adPQunVr3L4t/nWohH91cgKYHnx6Nkfn1+1E3ecEv46wa2Ap6j71JZfL8O2SPojd/Q4mjXDQLg+e3QMpUePwtk9bJrlMBS91Teoet++gy1MqlUhKSsLMmTPLLA8ODkZwcDCjVM/1WDIBTT2dkR2XiLxbGega6I9zizYjO+4mhh1ciUad2+LwkI/w9M59KGysMHTPcuTduY/Ts7/mMnfW+RuYcHMHHl27AwA4MeVLFOcqMXjnYlg0tEHkiGVMc/NOJpNhzQIvDJ/1qyj7q1/PHMunu4uyr9rw7NIEnl2aYNuhvwAA/5zKPhNhR6h/eH82GY06t8PT1CycW/g9NGo10/5RnR5t2bg+Bm//BGqVCsV5BYievg4yuZxJ35bMO2hbW1uUlpYiMDCQdZQy7JzawNzWBlH+y2HZqD4UNla4tukwss7fgEZVihMBq5F6JEa7vaqgCNHT1zFM/JxQbgB4nJimvR/ti7v/HP/7KpaRJcWvfxsEjHTQvaEe1i70RtsWtqLsixAxCPWP17p1hNxCgajRK5B78y5av/n8qCer/lHdHl2cm4+IEf9E1KgVeHg1Ga2HeDDr25IZ0Lxq5t1Ze/P3jOir0KjL3i+5KEf3TR5Y0JW7YafWGHYoGB5L/sYiXo3069ePdYQyvlvaB0N6taxyvcxtq84zuBcFdMWUUY5iRyMSwltdA8L9o367Znh8Iw0A8Oj6HTTt6cQk4wvV7dEatVp72zyZmRxPU9g9hpIGdC1Z2NnCfeE4+O4PQrc5o2BhJ413OrpyH+gbiMiRy2BhVw9t3urJKGX1rF27lnWEMqwsFfjlmzfx/ojqv5O2MJdjzQIv/Guup9431SHGibe6BoT7R97tDDTr7QIAaN6vCywa2rCKCaBmPdq+eye8HfUFWvTrgqdpui+ZNBQa0LVUnKtE/Jc/IWr0Clz+4kcubwZfGV25X3ydFhUHO2dpnAQ0f/581hEqsLJUYFvwABz59k04d2io12sGebXApR9HYMGkrjScCZd1LdQ/Hl2/g9zENAzd9yksbK1R+IDtUcSa9OicP27hiO/HSIuIhcP4QXWQsnI0oGsp60ICmvX631+LfVwhk0vjRyqUW2Ftqf26qacTnnJ8k/tXnTlzhnWEKg0f0BY3Do3GidBhmPluZ3h3bQK7+hawtjJDk0ZWGOzdEp9MccP1g6NwPNSPLlsiWjzWta6+d2XdPhwd8ymePVYi/dhlFhG1qtuj5eYvz50uVhaitKjYoPmESOYsbl7l3rwLtUoF3/1ByL50E6qCIsgUZtr1Ppvno5mXM+p3aIFr3/2Mu0f5eGSjUO4Gr7dA33UzUZJfBGVaNv74cg/jtMZBJpNhoFdLDPSq+nNpQqRAsO/JZPDd9ynUpWpknvkTOfF/8ZsVFXt0YdZj9Fw+ERq1Bs9ylTgdyO5qGxrQIrgcslv73+2G90LXQH/kpWQi6/wNRE8r+/mRwsYK/Td+iJw/2F+rLZT7l7cWVdh+8M7FKMzKrcOEhBBeCfWPyp4YxbJ/VKdHA0DUqLL5WfVtGtAiSw2PQWp4TJXrVQVFXF5HrCs3wP9lVrruV0yIFEmhrqXUP/TJWh6rvk0DWofGrh0k+b1Z5Wb58zp48CDdFpEYHRZ1TX2Pj+/L5dOsCCmP16f+mKIX127res410Y3qmgiRxinHhBBCiImhAU0IIYRwiAY0MRpfffUV6wiEiI7q2nTRgCZGw9nZmXUEQkRHdW26aEATozF8+HDWEQgRHdW16aIBTQghhHCIBjQxGu7u7qwjECI6qmvTRQOaGI34+HjWEQgRHdW16aIBTQghhHCIBjQhhBDCIRrQxGjs2LGDdQRCREd1bbpoQBNCCCEcogFNjMakSZNYRyBEdFTXposeN6nDhWXb8Oh6CpPv3di1A7yDA2r0Wla5a5OZEMIH6nvVY6i+RwNah0fXU5B1/gbrGNUm1dyEEPak2j+kmrsqdIibGI0PPviAdQRCREd1bbpoQBOjMXXqVNYRCBEd1bXpogFNjIafnx/rCEatpESNq0mPtF/fufcUGo2GYSLTQHVtuugzaGI0cnJyWEcwOkXPVNhzNAWhB24i9loOnhWXatd1GLYHr9lZ4q3erTDz3c7o694MMpmMYVrjRHVtumhAi8R3fxCaeDhCrVJBU6qGMi0bVzYcQOqR86yjVUmKmUndiTh9F1ODzuBedkGV2zzMfYawyGSERSZjsHdLbA3qh3Yt69dhSsKSFHuIlDLTIW4RXVm/Dz90mogwlwDc2nMKPt/NQf32zVnHEiTFzFWhB9uLQ63WYO4XMRg+61fB4Vze8QsZ6DLqAI5Epxkwnenhva6l2EOkkpkGtAFoStVI2nUMcnMFGru2Zx1HL1LMXN7OnTtZR5A8jUaDmZ+fw4Yfrtfo9coCFfznHkP47zSkxSKVupZiD+E9Mw1oA5CbK+A0aSgAIC85k3Ea/Ugxc3khISGsI0jef3+5hc17EwW30VydAs3VKVWuV5Vq8LdPonEvK1/seCZJKnUtxR7Ce2auP4NWq9VYu3YtNm/ejLt378LJyQlff/01pk6dCh8fH2zZsoV1xDLcPhyFLtPfgcLWCpqSUpyd/x0eJ6QCAAZsnIPkg2eQfuwSAGDQtkVI3H4UGdFXWEYWzOwwfhA6jvHRbmvbrimyLiTi9KwNrOIKOnToEJYsWcI6hmTdzynAnC9iRNnXE2UxpgWfxZFv3xJlf7X1NL8YBUWlaNzAEubm0npfwntdU98zHK4rdcqUKQgODsa0adMQGRmJcePGYfz48UhOToaHhwfreBVc/foAdjtPwo+uk5F+Ih7N+3bRrotdvg3ui96DwsYKbf28UZxXwLxIAeHMf4WdQNToFYgavQLR09dBVfAM8at2M0xLDGnTT4nIfVos2v7Cf7+LKzcfira/mog6k44h/4hEg97/RfOBu2E/YBfmfxmDzAf6f7ZOhFHfMxxuB3RYWBi2b9+Ow4cPY+HChRg4cCCWLl2K3r17Q6VSoUePHqwjVqn4ST7OLtiE1oN7oM1QTwBA0cM8JISGw/uzyeg2dzRiV2xnG7KcyjJryWTov3EOLof8AGX6AzYBiUGpVGps2S98aLsmvt8j/j71tWHXNQybeRSn4l4euszLL8H6XdfR872fcefeU2bZjBH1PfFxO6BDQkLg6+sLHx+fMss7deoEc3NzuLm54fHjx3j77bfh6OiIbt264a233sKtW7cYJS6rOFeJ61uOoMfiCcD/rg29tecUGnRsgYTQCBTnKhknrKiyzADQfcFY5CamIi0qjmE63cLDw1lHkKwbybm4n1Mo+n6Px2aIvk99xCfkYO7qCwCAUnXZm6loNMD9h4X4vyWnGCSrPinVNfU9cXE5oNPT03Ht2jWMHTu2wrq0tDS4urrC0tISMpkMc+fORVJSEq5cuYK3334bAQH8PEkp4d/hsGlqh05jX/6R8TTlPvLu3GeYSlj5zC36dUVLn264GLyLcTLdEhPZvVuTuks3DHMzjL9S8/BExMPm+vrupwTIBe6ZolZrcDY+m/kheH1Ira6p74lHpuHwXn0xMTHo3bs3wsPDy9zmrrCwEB07dsSwYcOwdevWCq+7ePEiRo4cifT0dJ3fQ987Hn3ceACcLZroH16HfutnIWn3cWTH6v6lSyx+gC8e/V6j7yNGbusmdhi6/1Mcm/C53od4apNZiD4PDAgNDdW5XWhoqFiRjEvTt4FmI8ssEjpTW4jMrdzvZtIy4FkdnyHr9C/Awl73dhlhwMPjhs9TBV7rmvqeYfuevmOXy7O47e2f/2IlJSWVGdCrV69GZmZmlSeIrV+/HiNHjqyLiCbBbd4YWNS3Qb8Ns7XLnty+h/OL+Dp7npCK9Dw4KOPyICJhiKe+x+U7aLVaDXd3d2RmZmLNmjVo1aoV9u3bh4iICKSlpSEmJgbe3t5lXhMUFITIyEicOHECNjY2omWJHLWc2fNFm/V2wbADK2v0Wla5a5NZSFyc7s+BvLy8EBsbK7iNp6en4HpT9f2eBMz47Jxe2754Z13hnXIVsk9NQJPG1jXOVhNjFxzHgeOpUKuF29up//jBp2eLOkpVEa91TX2vegzV97j881Eul2Pv3r1wdXXFjBkzEBAQAHt7e8yaNQtmZmZwc3Mrs/1nn32GI0eOICoqStThTKRl8eLFrCNIVo/OehwOroE2zevV+XAGgBnjOgsOZ7lcBuf2DTHAg7/bO5ZHdW26uDzEDQCOjo44efJkmWUTJ06Ei4sLrK1f/sIHBQUhIiICv/32G+zs7Oo4JeGJv78/6wiS5ebYCA1szZGnLBF1v/17sBmAA71aYNZ7nbHxxwTIALw6quVyGawszLAzxEcST9+iujZdXL6DrsrFixfLfP58/fp1fPrpp3j48CHeeOMNdO/eHd27d2cXkDDl5eXFOoJkWVkqMOkdB9H3O3WMk+j71IdMJsM3i3tjw8e90KpZvVeWA0P7tMK5/74Nzy7inQRlSFTXpovbd9DlKZVKJCUlYebMmdplrq6u9MB4QkQSON4Fm/cmorhELcr+erraMz2ELJPJ8OHfXDHrvc5QuG8DAKQefRdtmtsyy0RIdUhmQNva2qK0tFT3hgz0WDIBTT2dkR2XiLxbGega6I9zizYjPz0H/b8JhEajQUHGQ5wO/AYatRqDdy6GRUMbRI5YxmXurPM30GKAG9wC/SGTyxEXtANPbmVg6J7lyLtzH6dnf800NzEMh3YN8emMHljy9cVa78tcIcd/gvpzcQjZzOzlgUIazrVXVd/IjruJYQdXolHntjg85CM8/d91z64z3kFbX08o03NwZs63MLMwr/NeUt3MANDWzxveKwOwt+d0KGysmPQ/SR3i5pGdUxuY29ogyn85LBvVh8LGCtc2HUbW+RsozsvHsYmrEOW/HE/vZqPVYHcAwPG/r2KcWji3mZUFnCa+iV/fDUbU6BV4eDUZqoIiRE9fxzq2oH79+rGOIHkfvd8VQ/u0EtxG5rZV5xncGz7uha6OjcWMZrJ4qmuhvqFRleJEwGqkHnn5wBWr1xqgRd8uiByxDI9vpKKtr1ed95LqZn6h/fBeyM94fiMbVv2PBnQtNfPurL35e0b0VWjULw8PFj/JR8nT5zfl15SUQlMqzqFDMQjlbuLhCI1agyG7l6L/N4FQWFuyilkta9euZR1B8hQKOQ6sGwK//q1r9HqZDFi/yBsz3u0scjLTxVNdC/UNACjKeVLm69e6dcT9c8+fLZ55+iqa9nSsm6CvqG5mAGg1yB0Zp/+ERsO2Z9OAriULO1u4LxwH3/1B6DZnFCzsKh5Cs27WCC0GuHHxFJcXhHJbN2kIm2aNcGzC58iOuwnHv7/JMKn+5s+fzzqCUbCxVuDw129i7UfesLI00/t1ndo2wO/bhmPO/3XRvTHRG091rU+/K7N9w3ooUT6/x3txXgEsGtQT3N4QqpsZADqNewPJ+8W/I2J1SeYzaF4V5yoR/+VPuPvrRbQe4oF6LV8rs15uoUD/DbNx7qPvuXoHLZS7OK8AWbEJ0KjVyDx7DV1mvMMwqf7OnDnDOoLRMDOTY97ELhg9pD02701E6IGbyH5UVOm2PTq/hpnvdsYEv46wtqKWIjae6lpXvyuvJK8A9Vo838a8vjWK8/LrImYZ1c3cvG8XPLiUBHWJqo4SVo3eQddS1oUENOvlAgBo3scVMnnZH2mfL6cjcXsUniTpvj94XRLKnXPlNuwcnh/ibOzaHk/TsplkJOy1bWGLzz/sicwTE5D0yxj89OVAbPpnH2xe3hdRm4Yi+9QEXPppJKaMcqLhbAJ09bvycv64hWa9n2/fsr8bHlxKMnjG8qqbuZFzW7R5qyfe3L0Udo5t4P7xe3URs1L0G1VLuTfvQq1SwXd/ELIv3YSqoAgyxfPDgk08HNHOzxu2rZvA5R9v40ZoONIihW/ZV1eEcj97mIf752/A9+BKlBY+Q/TMDYzTEtbkchkc2jWEQ7uGrKMQhoT6BgD4bJ6PZl7OqN+hBa599zPuHo1DVkwChv0cjPz0HNz4d90/OrO6mRO2RiBhawQAYNjPwYj/4sc6z/wCDWgRXA7Zrf3vdsN7oWugP/JSMpF1/gZ+cJhYYfvBOxejMCu3DhNWTij3jS1HcGPLEe16hY0V+m/8EDl/3GYRVS+67ldMiBTxVtdCfSN6WsUT2q5tPIRrGw9pv2bRS6qb+YUXl8Ky6n80oEWWGh6D1PCKp+y/iofLrMrTlVtVUMT8um1dDh48SLdFJEaH57rWp9+Vx7qXSCkzDWgdGrt2kOT3ZpWb5c9r1apV3DYyQmqKRV1T3+Pj+9KA1sE7OIB1hBqRam5CCHtS7R9SzV0VOoubEEII4RANaGI0vvrqK9YRCBEd1bXpogFNjIazszPrCISIjuradNGAJkZj+PDhrCMQIjqqa9NFA5oQQgjhEA1oQgghhEM0oInRcHd3Zx2BENFRXZsuGtDEaMTHx7OOQIjoqK5NFw1oQgghhEM0oAkhhBAO0YAmRmPHjh2sIxAiOqpr00UDmhBCCOEQDWhiNCZNmsQ6AiGio7o2XfQ0Kx0uLNuGR9dTmHzvxq4davx0Fla5a5OZEMIH6nvVY6i+RwNah0fXU5B1/gbrGNUm1dyEEPak2j+kmrsqdIibGI0PPviAdQRCREd1bbpoQBOjMXXqVNYRCBEd1bXpogFNjIafnx/rCIQjSXeeYO3OPzHh45PaZQOnRGDuFzH4MfI2CgpVDNPpj+radNFn0CLx3R+EJh6OUKtU0JSqoUzLxpUNB5B65DzraFWSYmYhOTk5rCMQDly4mo1lGy/ht/MZFdadisvEqbhMbPgBaNTAAv8Y7YR/Tu2O+vUsGCTVD891LcUeIqXM9A5aRFfW78MPnSYizCUAt/acgs93c1C/fXPWsQRJMTMhlSkpUePjdbHo8/cjlQ7n8h7nFWP1tj/RdfRBRF/MrIOExkmKPUQqmWlAG4CmVI2kXccgN1egsWt71nH0IsXM5Tk7O7OOQBgpKVFj7MLjWL3tT6jVmmq9NjVDiTenRuHnk6kGSlc7UqlrKfYQ3jPTgDYAubkCTpOGAgDykqXxl7kUM5e3c+dO1hEII7NXncPPJ9Nq/PoSlRrjFp5A3LUHIqYSh1TqWoo9hPfM9Bm0iNw+HIUu09+BwtYKmpJSnJ3/HR4nPP+rfMDGOUg+eAbpxy4BAAZtW4TE7UeREX2FZWTBzA7jB6HjGB/ttrbtmiLrQiJOz9rAKq6gkJAQLFmyhHUMUsciTt/Fln03BbfRXJ0CAJC5ba1ym+ISNd5f9jsu/TgCVpb8tEbe65r6nuFw/Q5arVZjzZo1cHBwgJWVFbp164bo6Gg4OTlxeenB1a8PYLfzJPzoOhnpJ+LRvG8X7brY5dvgvug9KGys0NbPG8V5BcyLFBDO/FfYCUSNXoGo0SsQPX0dVAXPEL9qN8O0wg4dOsQ6AqljarUGH/5LvJN7btzOxaY9iaLtTwy81zX1PcPhekBPmTIFwcHBmDZtGiIjIzFu3DiMHz8eycnJ8PDwYB2vSsVP8nF2wSa0HtwDbYZ6AgCKHuYhITQc3p9NRre5oxG7YjvbkOVUlllLJkP/jXNwOeQHKNP5OwRITNev5+7h9t2nou5z056Ean+OTajvGQK3AzosLAzbt2/H4cOHsXDhQgwcOBBLly5F7969oVKp0KNHD9YRBRXnKnF9yxH0WDwBkMkAALf2nEKDji2QEBqB4lwl44QVVZYZALovGIvcxFSkRcUxTEdIRbsjbou+z79S83DxOr+XNvGM+p64uB3QISEh8PX1hY+PT5nlnTp1grm5Odzc3AAAI0eOhJubG9zd3eHl5YVjx46xiFuphH+Hw6apHTqNfflveJpyH3l37jNMJax85hb9uqKlTzdcDN7FOJlu4eHhrCOQOhZ33TDvbC4aaL81IbW6pr4nHn7OhHhFeno6rl27hnnz5lVYl5aWBldXV1haWgIAtm/fDjs7OwBAfHw83njjDTx69AhmZmZ1GRlRo1dUWFaiLESYC79PdtKV2bqJHbxDpuDYhM+hLuH/rkuJiYlo0qQJ6xikjqhUaty888Qg+75267FB9lsTPNc19T3D4nZAA0Dz5mUvHC8sLER0dDSGDRumXfZiOAPAkydPIJPJoNHo/vxI9sqhDCEfNx4AZws2vxzR0dHw0zNneWLkdps3Bhb1bdBvw2ztsie37+H8oi1VvqY2mYXo88CA0NBQndt5eXmJFYmwJrMAunxXZtGLs7WrUtX68md3b/p+KzYt61e7fHrgta6p7xm27+kzowBOB7S9vT0AICkpqcx9aFevXo3MzMwKJ4jNmjULkZGRePLkCfbv3w+Fgst/FgDgzNyNrCPo7cKSUFxYEso6BiGV05T+7381ZT47FGffJeLuz8RR36sZmUbfUV6H1Go13N3dkZmZiTVr1qBVq1bYt28fIiIikJaWhpiYGHh7e1d4XXR0NObNm4fff/8dtra2omSJHLWc2fNFm/V2wbADK2v0Wla5a5NZSFyc7hM1vLy8EBsbK7iNp6en4HoiLY7/by/+Ss3TuZ0+10G/6utPeiFwgmutsumD17qmvlc9hup7XJ4kJpfLsXfvXri6umLGjBkICAiAvb09Zs2aBTMzM+0JYuX5+PhALpfj7NmzdZyY8GDx4sWsI5A61tPF3jD7dTXMfmuC6tp0cXss2NHRESdPniyzbOLEiXBxcYG1tTUAQKlU4uHDh2jXrh2A5yeJ3b59G507d67zvIQ9f39/1hFIHXvX93WERSaLus92LW3h1YWfk7Kork0XtwO6MhcvXkSvXr20X+fn5+Pdd9+FUqmEQqGAlZUVdu3ahbZt2zJMSVjR51AgMS7D+7dB2xb1kJaZL9o+p491hpkZPwcXqa5Nl2QGtFKpRFJSEmbOnKld1qxZM8TExDBM9VyPJRPQ1NMZ2XGJyLuVga6B/ji3aDNyb97F4O2fQK1SoTivANHT10Eml2PonuXIu3Mfp2d/zWXukrwCeK18fslBvdb2SAiNwI1/h2PwzsWwaGiDyBHLmOYm5AWFQo61C70xZsEJUfbXsU19BI53EWVfxqqqvpEddxPDDq5Eo85tcXjIR3h65z5kCrMKyxQ2VnXeA6uTGQAm3NyBR9fuAABOTPkSxblKJv2Pnz8TdbC1tUVpaSkCAwNZRynDzqkNzG1tEOW/HJaN6kNhY4Vrmw4j6/wNFOfmI2LEPxE1agUeXk1G6yEeUBUUIXr6OtaxBXM/un5Hey/axzfScPe35ze6P/73VYxTE1LR6Dc7YIJfR8FtZG5bdZ4gZmYmw7aVA1DPxlzMeEZFqG9oVKU4EbAaqUdevmmqbFld98DqZgaAx4lp2h744u5nLPqfZAY0r5p5d9be/D0j+io0arV2nUatfn4JCACZmRxPU/h5nJlQ7hcU1pawbmqn/auSd/36Gf66VcKnrUH9MKRXyxq/Xi6XYcdnA9Dfo7nujesYT3Wtq28U5VS8cUxly+pSTTI37NQaww4Fw2PJ3+okY1Ukc4ibVxZ2tnCa+CZcp74NiwY2SPnlfJn/w+27d0Kvf/0Dpc+Kcf37XxgmLUtXbgBoNcgd907GM0pYfWvXrmUdgTBiZanAL9+8ibmrL2Dz3uo9jarZa9b4z8r+8OvfxkDpaoenutanb/CmJpkP9A1Eca4SvVdPRZu3euLurxfrKG1Z9A66lopzlYj/8idEjV6By1/8WOFm8Dl/3MIR34+RFhELh/GDGKWsSFduAGjr543UiAsM0tXM/PnzWUcgDFlZKvD9sr74bYsvvLvqPgvb2tIM08Y64/rBUdwOZ4Cvutanb/CmJplfbJMWFQc7Z3YnHdOArqWsCwlo1uv5SSXN+7hCJn/5I5WbvzxAUawsRGlRcZ3nq4pQbgCQKcxg59AKj2+ksohXI2fOnGEdgXBgSK9WiPnhHVz8cQRWTHfH8AFt4NS+IV5vXR/uzq8hYKQDNi/vi3vHx+P7ZX3xmp0V68iCeKprXX2DR9XNrLC21G7T1NOJ6Ud8dIi7lnJv3oVapYLv/iBkX7oJVUERZIrnD+po7NoePZdPhEatwbNcJU4Hsj1r+1VCuQGgRb8uyDxzjWFCQmrHw8UeHga6kYmp0tU3fDbPRzMvZ9Tv0ALXvvsZd4/GVbqM58z56Q/Qd91MlOQXQZmWjT++3FOneV9FA1oEl0N2a/+73fBe6Broj7yUTGSdv4GoUWWfnKKwsUL/jR8i5w/xn2NbXUK5M05dQcapK2W2H7xzMQqzcus4JSGEJ0J9I3paxc/Lyy9j0QOrm/mXtxZVWMai/9GAFllqeAxSw6u+NltVUMTldcS6cgP8X2ZFN3Mgxojnutanb5THugfWJDPApv/RgNahsWsHSX5vVrlZ/rwOHjxIt0UkRodFXVPf4+P7cvk0K0LK4/WpP4TUBtU1EcL/KXiEEEKICaIBTQghhHCIBjQxGl999RXrCISIjuradNGAJkbD2dmZdQRCREd1bbpoQBOjMXz4cNYRCBEd1bXpogFNCCGEcIiugyaSoM9lJCtWrKDLTYikUF0TIXQdNCGEEMIhOsRNCCGEcIgGNCGEEMIhGtCEEEIIh2hAE0IIIRyiAU0IIYRwiAY0IYQQwiEa0IQQQgiHaEATQgghHKIBTQghhHCIBjQhhBDCof8PSkKXX8dj3ggAAAAASUVORK5CYII=\n", + "text/plain": [ + "
" + ] + }, + "execution_count": 3, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "from qiskit.circuit.library import EfficientSU2\n", + "\n", + "# the rotation gates are chosen randomly, so we set a seed for reproducibility\n", + "ansatz = EfficientSU2(num_qubits, reps=1, entanglement='linear', insert_barriers=True) \n", + "ansatz.draw('mpl', style='iqx')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Solve with the VQE" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now that we have the problem and ansatz specified we can use the Variational Quantum Eigensolver (VQE) to solve for the minimal eigenvalue of our Hamiltonian.\n", + "\n", + "The VQE requires a classical optimization routine, along with an initial point, to calculate the parameter updates." + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "from qiskit.algorithms.optimizers import SPSA\n", + "\n", + "optimizer = SPSA(maxiter=50)\n", + "\n", + "np.random.seed(10) # seed for reproducibility\n", + "initial_point = np.random.random(ansatz.num_parameters)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "To access intermediate information we can pass a callback into the VQE. The callback is given the current number of function evaluations, the current parameters, function values and standard deviation in the expectation evaluation. " + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [], + "source": [ + "intermediate_info = {\n", + " 'nfev': [],\n", + " 'parameters': [],\n", + " 'energy': [],\n", + " 'stddev': []\n", + "}\n", + "\n", + "def callback(nfev, parameters, energy, stddev):\n", + " intermediate_info['nfev'].append(nfev)\n", + " intermediate_info['parameters'].append(parameters)\n", + " intermediate_info['energy'].append(energy)\n", + " intermediate_info['stddev'].append(stddev)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Local VQE\n", + "\n", + "Before running the Qiskit Runtime VQE program, let's first simulate this system locally using Qiskit's `VQE` class." + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [], + "source": [ + "from qiskit.providers.basicaer import QasmSimulatorPy # local simulator\n", + "from qiskit.algorithms import VQE\n", + "\n", + "local_vqe = VQE(ansatz=ansatz,\n", + " optimizer=optimizer,\n", + " initial_point=initial_point,\n", + " quantum_instance=QasmSimulatorPy(),\n", + " callback=callback)" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [], + "source": [ + "local_result = local_vqe.compute_minimum_eigenvalue(hamiltonian)" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Eigenvalue: -0.994140625\n", + "Target: -1\n" + ] + } + ], + "source": [ + "print('Eigenvalue:', local_result.eigenvalue)\n", + "print('Target:', target_energy)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "With the information from the callback, we can for instance compute the average of the estimation errors in the expectation evaluations. For an exact (statevector) simulation this would be 0, but for a shot-based readout as we have on real hardware, this will be a small finite error." + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Mean error: 0.014968352082564208\n" + ] + } + ], + "source": [ + "print('Mean error:', np.mean(intermediate_info['stddev']))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Runtime VQE: ``VQEProgram``\n", + "\n", + "To run the VQE using Qiskit Runtime, we only have to do very few changes from the local VQE run and mainly have to replace the `VQE` class by the `VQEProgram` class. Both follow the same `MinimumEigensolver` interface and thus share the `compute_minimum_eigenvalue` method to execute the algorithm and return the same type of result object. Merely the signature of the initializer differs sligthly.\n", + "\n", + "We start by choosing the provider with access to the Qiskit Runtime service and the backend to execute the circuits on. \n", + "\n", + "**Note:** To run this tutorial, replace the provider by your provider." + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [], + "source": [ + "from qiskit import IBMQ\n", + "\n", + "IBMQ.load_account()\n", + "provider = IBMQ.get_provider(project='qiskit-runtime') # replace by your designated provider" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We're using the IBM Montreal device, but you can replace this with another device that support the Qiskit Runtime." + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": {}, + "outputs": [], + "source": [ + "backend = provider.get_backend('ibmq_montreal')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Let's not forget to reset the callback!" + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "metadata": {}, + "outputs": [], + "source": [ + "intermediate_info = {\n", + " 'nfev': [],\n", + " 'parameters': [],\n", + " 'energy': [],\n", + " 'stddev': []\n", + "}\n", + "\n", + "def callback(nfev, parameters, energy, stddev):\n", + " intermediate_info['nfev'].append(nfev)\n", + " intermediate_info['parameters'].append(parameters)\n", + " intermediate_info['energy'].append(energy)\n", + " intermediate_info['stddev'].append(stddev)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Since the Qiskit Runtime VQE runs on real hardware, we might want to do readout error mitigation. For this purpose the `VQEProgram` supports a boolean flag `readout_error_mitigation` that can be set to `True` to enable error mitigation with a complete measurement fitter." + ] + }, + { + "cell_type": "code", + "execution_count": 13, + "metadata": {}, + "outputs": [], + "source": [ + "readout_error_mitigation = True" + ] + }, + { + "cell_type": "code", + "execution_count": 14, + "metadata": {}, + "outputs": [], + "source": [ + "from qiskit_runtime.vqe import VQEProgram\n", + "\n", + "# in this first release, the optimizer must be specified as dictionary\n", + "optimizer = {'name': 'SPSA',\n", + " 'maxiter': 50}\n", + "\n", + "runtime_vqe = VQEProgram(ansatz=ansatz,\n", + " optimizer=optimizer,\n", + " initial_point=initial_point,\n", + " provider=provider,\n", + " backend=backend,\n", + " shots=1024,\n", + " readout_error_mitigation=readout_error_mitigation,\n", + " callback=callback)" + ] + }, + { + "cell_type": "code", + "execution_count": 15, + "metadata": {}, + "outputs": [], + "source": [ + "runtime_result = runtime_vqe.compute_minimum_eigenvalue(hamiltonian)" + ] + }, + { + "cell_type": "code", + "execution_count": 16, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Eigenvalue: -0.947692038031636\n", + "Target: -1\n" + ] + } + ], + "source": [ + "print('Eigenvalue:', runtime_result.eigenvalue)\n", + "print('Target:', target_energy)" + ] + }, + { + "cell_type": "code", + "execution_count": 17, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Mean error: (0.009295459418135251+0j)\n" + ] + } + ], + "source": [ + "print('Mean error:', np.mean(intermediate_info['stddev']))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Additionally to final results, the `VQEProgram` also returns the history of the optimization, such that we can visualize the " + ] + }, + { + "cell_type": "code", + "execution_count": 18, + "metadata": {}, + "outputs": [], + "source": [ + "history = runtime_result.optimizer_history\n", + "\n", + "loss = history['loss']\n", + "timestamps = history['time']\n", + "\n", + "runtimes = np.concatenate(([0], np.diff(timestamps)))\n", + "runtimes_in_min = runtimes / 60" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "From this, we can extract the times per iteration and statistics on it:" + ] + }, + { + "cell_type": "code", + "execution_count": 19, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Total time taken: 24.77min\n", + "\n", + "Median time of the iterations: 26.76s\n", + "Average time per iteration: 29.72s\n", + "Standard deviation: 13.35s\n" + ] + } + ], + "source": [ + "print(f'Total time taken: {np.sum(runtimes_in_min):.2f}min\\n')\n", + "\n", + "# note that the median and average might differ, since the device get's calibrated every hour\n", + "# resulting in some iteration times that are much larger than the rest\n", + "print(f'Median time of the iterations: {np.median(runtimes):.2f}s')\n", + "print(f'Average time per iteration: {np.mean(runtimes):.2f}s')\n", + "print(f'Standard deviation: {np.std(runtimes):.2f}s')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "or plot the loss against the total time taken:" + ] + }, + { + "cell_type": "code", + "execution_count": 20, + "metadata": { + "scrolled": false + }, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAuYAAAH8CAYAAACdPGozAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAACGBUlEQVR4nOzdd3zV9fXH8dcne68bIECALLbK3qjgrHXUUScOrLittVq11l+tdjtqq3VbV7WO1lVXraOiQpCNgyFksOe9WWSP+/n9cW9iEkJIIHckeT8fj/u4937nuXwDnPvJ+Z6PsdYiIiIiIiKBFRLoAERERERERIm5iIiIiEhQUGIuIiIiIhIElJiLiIiIiAQBJeYiIiIiIkFAibmIiIiISBBQYi4i0gnGmAxjjDXGzO3kfvONMes6sN1GY8yzBxtfMDPGPGuM2RjoOEREgpUScxHp9Ywxg40xj3mT4hpjzG5jzBvGmGlBENsAY8ydxpixHdy+8YtD48NtjCkyxvzHH5/HGDPKG2+Gr891KIwxZ3r/fK5uZ5vx3m1+32r5YcaYfxhjthtjao0x24wxLxhjRrVxjFmtrkfrx1W++Hwi0j2FBToAEZFAMsbMAN7zvn0SWAv0B+YCC40x11prH222yyYgGqjzUUjDAXez9wOAXwEbgVWdOM4rwDtAKDAMuAb4xBgzyVr7dZdE2rZReOKdjyfm5i4neAaE3gVKgAuAR/ezzRzv8wuNC4wxZwIvAUXAU0AhkAlcBvzQGHOutfbfbRzrYeCLNpYvPpjgRaRnUmIuIr2WMSYZeBWoAmZYa/ObrfsT8F/gr8aYldbaLwCsZ7rkal/FZK2t6aJDrbLWNk8oPwf+A1yNJ0n3O2utr77MdJq1tsYY8y9gnjFmsLV2c/P1xpgQ4DxgpbV2jXdZNvA8UAAcaa11Ntv+L8DnwAvGmCOstYWtTrnAWvuy7z6RiPQEwTJyISISCFcCacDNzZNyAGttFXCJ9+0djcvbqjE3xsQZY+4zxhR6S2H2eGvKj2rv5MaYo4wxZcaYt40xkd5lTTXmxphZwFLv5s80K3+48yA+6+fe5+xm57/TGGPbiGuu9zwZzZZtNMa8b4yZaYxZYoypNsYUGGMubr4f8C/v20+axTvXu36fGnPv+seMMT80xqw2xlQZYxY3lu4YYy43xmzwnu9TY0xWG/FOMsa8Z4wp9e6/wBgzuwN/Ji8ABji/jXWz8fy24oVmy24GYoArmiflAN73VwJx3u1ERDpNibmI9Gan4hn9/mdbK72jnguA44wxUe0c51Hgx8AbeEaj7wb2AGP2t4Mx5ng8I9j/Bc7cz0j5Wr77UvAEcJH38Xo7sexPhve5+CD2bZSJ5zcMHwI3eY/1rDFmtHf9Z8CD3te/57t4PzvAcacDf8YzGn0nnnKed7z11z/F8+d7DzAVeLb5jsaYo/F86UgBfg3cCkQCH3i/2LTnczylSXPaWDcHT0nRS82WnQpstNZ+3sb2WGs/w1O+c2obq+ONMaltPEIPEKOI9CIqZRGR3mwU8O0Byke+BI4GhgL7q80+BXjSWntjR05qjDkVz8jyP4FLrbUNbW1nrd1ljPkPnoRzUfPSlA6IMcak4qkxHwrc713+aieO0dow4GhvAoox5p/AFuBS4GfW2gJvycz1wIfW2vkdPO4IYGTjby2MMcXA43g+91Brbal3eRhwmzEmx1qbZ4wx3u0WAMd7y4wwxjwGrMTz5WD6/k5qrbXGmBe9xxxtrV3t3T8KOAv4n7V2h3dZIp4R9Lbqx5v7CjjNGBNvrd3bbPkT3kdrI4EDdusRkd5BI+Yi0pvFA3sPsE3j+vh2tikFphhjBh7ohMaYHwKvAc8Bc/eXlHeBX+IZtd+JZ2R4BHCTtfZQEvP1jUk5gLV2D/AtsE95SSd90qqUqPGGyNcbk/JWyxvPNwbP6PqLgKNxFBpIwDOqP8UYE3OAczd+2Wk+an6K9xjNvwg1Xv+D/Xn5HXB8G4/NiIh4acRcRHqzvbSfcNNs/e52trkZT6K92RizEngfeN5a+22r7QYDLwNvWmuvPIh4O+Mp77migFnADcA+9eSd1FYSWQwkd/FxG5PxLftZ3ni+Yd7np9o5tgOo3N9Ka+0a7zU73xhzu3fUfQ6eG4Kblwx15Ata43oLOFst/8Za+9EB9hWRXk6JuYj0ZmuA8caYyHbKWY4AaoFt+zuItfZf3hKOHwAn4CnluMUYM9da+2KzTXfhqWk+0RgzzVq7qEs+RdvymiWC7xhj6oB7jTGfW2uXNYa+n333V/e8v9F9c7BBHuC4Bzpf4299fw4s38+2ezpw/heAPwHTjTFrgO8DbzQvRbHWlhpjtuP5eWjPEcBWa21tB84rItKCSllEpDd7G8+I8tltrfR2JTkSeMfbpWW/rLU7rbWPW2vPwjMyXgDc1WqzGjw3Bq4B3jMdmzToUEe5G/0BKAN+22xZMYAxJqnVtkMO4TxdFW9HNJa/7LXWfrSfR0daW76E50vABcAPgQhalrE0ehvINMbMbOsgxpgj8dxk+6+21ouIHIgScxHpzR7HU4N9jzEms/kKY0w08Ayezhz37O8AxphQ742BTay1JXgmnklqvb13FPZ7eMo3PjDGjDhAjBXe50MqF7HWlgGP4RmtH+dd3JjYNrV1NMbE8l2byIPRJfF20HIgD7jRGLNPiYkxpk9HDuK9wfN/wDl4PrsTT7ec1u7DUxbzuDHG0epcDjx/vqXAQ534DCIiTVTKIiK9lrW22Hsz5nvASmNM48yfaXhm/swErrXWtjc7YzywzRjzGp4OLmXADDzJd5sJmve8J+BpI/iRMebINiakaZSPZ2T7amNMOZ5a52+std907tMCnlaGNwG34UlCP8DzBeEpY8y9eEaNf4Sn/GPwQRwfPN1QGvB0OknCU6u9uJ3Pd9CstW5jzGV4avrXGGOeBrbi6Z5yNJ6Sl470MwfPCPnxQCrwSFuTIXk7wVyMZ4T9a2PM3/C0R8wA5gGJwLn7+awzvV1lWltjrV3RwRhFpIdTYi4ivZq1dqEx5gjgF3iS1QF4/m20wDEdaPlXiWe69eOB04BwPKPlPwMeaOe8u4wxx+HpmPKxNznfp47dWltnjLkITynKw97j3wV0OjG31u40xvwDuMQYM8xau94YcwbwCPAbPL89+AueLwLPdPb4zT7X5Xj+PJ/EU69+KZ4/ky5nrf3MGDMVTxeaa/B0U9mJZ2Kmv3XiUK/j6ZceQ9tlLI3ne80YMx7Pl5vLgb54fvtcBYy31u6v9eG13kdrfwKUmIsIAMbb9lVERLyMMcfiGUX/DPh+ME0lL8HHO4r+LPB3a+3cwEYjIt2ZasxFRFqx1n6Mp9b4WOBp70Q2Im2y1v4dzwj6JcaY3wU6HhHpvjRiLiIiIiISBDRiLiIiIiISBII2MTfGRBpj/mqMcRpjKowxbxlj0g+wz23GmKXGmDJjzB5jzNvGmMP8FbOIiIiIyMEK2sQcT2eAs4Dz8UzwkYBn9rr9zUgHnmmnHwGmA8cA9XhakaX4NFIRERERkUMUlDXm3sk69gCXWmv/4V02CM9U1idZa9ua+KGt48ThmezhdGvt2+1tm5qaajMyMg4pbhERERGRA1m+fLnTWrvPJGjB2sd8Ap5evR80LrDWbjHGrMUzGt6hxBzPxB8heKedbk9GRgbLli07iFBFRERERDrOGLOpreXBWsqShmfmOGer5bu86zrqAWAVsKitlcaYK4wxy4wxy/bs2XMwcYqIiIiIdAm/JubGmN8aY+wBHrO66Fz3AzOBs6y1DW1tY619wlo70Vo7sU+ffX6bICIiIiLiN/4uZfkL7Ux17LUZmIpnGudUPLXmjfrhmb66XcaYPwPnAbOttQUHFamIiIiIiB/5NTG31jrZtzxlH8aY5UAdcDzwondZOjASyD3Avg8A5+JJytcdaswiIiIiIv4QlDXm1tpS4CngHmPMccaYccDzwFfAR43bGWPWGWOua/b+YeBS4AKg2BiT5n3E+fcTiIiIiIh0TrB2ZQG4AU8f8leAaOBj4OJW9eLD8ZS7NLrG+/xxq2PdBdzpkyhFRESkw8rKyti9ezd1dXWBDkXEJ8LDw+nbty8JCQmd3jdoE3NrbQ3wY+9jf9uY9t4LfLOtlD998C2PXTSByLD25mYSERHxrbKyMnbt2sXAgQOJjo7GGP23LT2LtZaqqiq2bdsG0OnkPChLWaTr/HPZFj75dg9biqoCHYqIiPRyu3fvZuDAgcTExCgplx7JGENMTAwDBw5k9+7dnd5fiXkPtyDPc69tUUVtgCMREZHerq6ujujo6ECHIeJz0dHRB1WupcS8B9teUkXBngoAiipqAhyNiIgIGimXXuFgf86VmPdgC/O+60xZVKGbbERERESCmRLzHmxhnpPkmHBAI+YiIiL+MHfuXE455ZQObfvss88SF7f/js4bN27EGMOyZcu6Kjy/6cyfg3xHiXkPZa1lQZ6Lo4b1IS4yTCPmIiIih2Dbtm1cccUVpKenExERwcCBA7n88svZunVri+0eeOABXnjhQJOcd8ygQYPYsWMHY8eOBWD+/PkYY3A625+rsTGhb3wkJiYydepU3n777S6Jq7n9xdSVfw4dVVtbS2pqKr/5zW/aXP/oo48SExNDaWkpAG63m7/+9a+MHTuW6OhoEhISOOaYY/jPf/7TYr/Gz9jWY926rp3LUol5D7V+VznO8hpm5KSSHBuuEXMREZGDVFhYyMSJE/nmm2947rnnyMvL44UXXmD16tVMmjSJjRs3Nm2bmJhIUlJSl5w3NDSUtLQ0wsIOrrv1+++/z44dO1i8eDGTJ0/mrLPO4ptvvumS2A6kK/8cOioiIoKLLrqIZ599FmvtPuufeuopfvjDH5KYmAjABRdcwC9/+UuuvPJKVq9ezaJFi5g0aRKnnHIKjz766D77r169mh07drR4DB06tEs/gxLzHqqxG8uMnFRSYiNxqSuLiIjIQbn22msJCQnho48+4thjj2Xw4MHMnj2bjz76iJCQEK699tqmbVuXcHz22WdMnTqVuLg4EhMTmTx58n6T4+LiYmbMmMGJJ55IRUVFi1KWjRs3Mnv2bAD69OmDMYa5c+e2G7fD4SAtLY0RI0bwu9/9jrq6Oj755BNg/2UyxhheffXVFtu89tprHH/88cTExDBq1Cg+/PDDpvX7i6n1n8OsWbO4+uqruemmm0hJSaFPnz488MAD1NTUcO2115KUlMTgwYN5/vnnW8Szbds2zjvvPJKTk0lOTubkk09mw4YN+/3M8+bNo6CggPnz57dY/uWXX7J8+XLmzZsHwD//+U9eeeUVnn32Wa6++mqysrIYPXo0d999Nz/+8Y+54YYb9vltSN++fUlLS2vxCA3t2jlilJj3UAvznGSlxjIwKZqUmHCKK5WYi4iIdFZRURHvv/8+1157LTExMS3WxcTEcM011/Cf//yH4uLiffatr6/nBz/4ATNnzuTLL79k8eLF3HDDDW0mc9u3b+eoo44iPT2dt99+m9jY2BbrBw0axGuvvQZ8N3L7wAMPdOgz1NXV8eSTTwKeWSk76/bbb+f666/nyy+/ZNKkSZx33nmUl5d3OqZ//OMfxMfHs3jxYn7+859zww03cPrppzNs2DCWLVvGJZdcwrx589ixYwcAlZWVzJ49m6ioKD799FMWLVpE//79Oe6446isrGzzHKNHj2bKlCk8/fTTLZY/9dRTDB06lKOOOqoplqFDh3L66afvc4ybb76Z2traps/mT0E786ccvLoGN18UuDhrfDoAKbGRfLtzb4CjEhER2dddb69mzfYyv55z1IAEfnXq6A5tu2HDBqy1jBw5su1jjRqFtZYNGzYwefLkFuvKysooKSnh1FNPJTs7G4ARI0bsc4y8vDyOP/54TjzxRB555BFCQvYdNw0NDSUlJQXwjNympqYeMPajjjqKkJAQqqqqcLvdZGZmcs455xxwv9Z++tOfcuqppwLw+9//nr///e+sWrWKmTNndiqm0aNHc+eddwJw44038sc//pHw8HB+8pOfAHDHHXdw9913s3DhQn74wx/y8ssvY63lmWeeaWo/+Pjjj9O3b1/eeeed/X6WefPmcf311/PQQw+RmJhITU0N//jHP7jllluatlm/fv1+r+nAgQNJSEjg22+/bbE8IyOjxfukpKR9RtUPlUbMe6BVW0qorG1gRo7nL4gjLgJXRW2b9VYiIiJy6CIiIvZZlpKSwty5cznxxBM5+eSTuf/++9m8eXOLbWpra5k5cybf//73eeyxx9pMyg/Wiy++yMqVK3nrrbcYOnQoTz/9dFMi3RlHHHFE0+sBAwYAHNSsls2PY4yhb9++HH744U3LwsPDSU5Objr28uXLKSwsJD4+nri4uKZyoOLiYvLz8/d7nvPOO4/Q0FBeeuklAN58803Kysq45JJLOhVv62v6ySefsGrVqqbH559/3qnjdYRGzHugBRuchBiYluUAIDkmgpp6N1V1DcRE6JKLiEjw6OjIdaAMHToUYwxr1qzhjDPO2Gf9mjVrCAsLIzMzs839n3nmGW644Qbef/993nrrLW6//XbefPNNTjzxRMCTjJ5wwgm89957bNq0iSFDhnRZ7Onp6QwdOpShQ4cSFxfH2WefzZo1a0hNTW36AtB80G5/M1U2L39pHLl2u92djqd1GY0xps1ljcd2u92MHTuWl19+eZ9jtfcFIy4ujnPOOYenn36aq666iqeeeoqTTz6ZtLS0pm2GDRvG2rVr29x/27ZtlJWVMWzYsBbLMzMzO/SbikOhEfMAWrO9jIueWsyO0qouPe7CPCeHpyeR6O1h7oj1fONzlavOXEREpDNSUlL43ve+xyOPPLJPXXNlZSUPP/wwZ5xxRlOnj7aMGTOGW2+9lfnz5zNr1iyee+65pnXGGJ599llmzpzJ7Nmz9xlRb65xBLehoaHTn+Poo49m1KhR/PrXvwY8N2sCTfXcAKtWrer0cQ8lpgMZP348eXl5pKamkpOT0+JxoJH/efPmsXTpUt555x0+/vjjpps+G82ZM4cNGzbw5ptv7rPvPffcQ1RUFOeee25XfpwOUWIeQLGRoXy+wckbK7d12TH3VtexcksJM3McTcuSvYm5bgAVERHpvIcffpiGhgaOO+44/ve//7Flyxbmz5/P8ccfT3h4OA8++GCb+xUWFvLzn/+c3NxcNm3axCeffMJXX33FqFGjWmwXEhLCc889x/Tp05k1a9Z+k/MhQ4ZgjOHdd99lz549lJeXd+pz3HTTTTzxxBNs2bKF6Ohopk6dyt13383q1avJzc3lZz/7WaeO1xUxtWfOnDn069ePH/zgB3z66acUFhby2WefcdNNN7XbmQVg2rRpjBo1iosvvpi0tDROOumkFuvPPvtszj33XObOnctjjz1GYWEha9as4ec//zkPPfQQf/vb33A4HC322b17Nzt37mzxqK3t2txKiXkADXHEMikjmdeWb+2y+u8lhUU0uG1TfTlASuOIuVomioiIdFpmZibLli1j9OjRXHTRRWRkZDB79mxCQkJYtWpVixKJ5mJiYli/fj1nn302w4YN45JLLmHOnDnceuut+2zbPDnf38j5wIEDueuuu7j99tvp168f1113Xac+xymnnEJGRkbTBDyNnUsmTZrElVdeyW9/+9tOHa8rYmpPTEwMn332GVlZWZx99tmMGDGCSy65hOLiYpKTkw+4/2WXXUZxcTFz587dpxOOMYYXX3yRX//61zz66KOMHDmS0aNH85e//IVPP/2UOXPm7HO80aNH079//xaPzz77rMs+L4DRDYEeEydOtIGY8valJZu57fWv+fe1MxgzKOmQj3fX26t5aclmVt1xAlHhnh/Cjc4KZt03nz+dPYazJqQf8jlEREQOxtq1a/fbCaO7eeihh7jpppv417/+xWmnnRbocKQLrF+/ntmzZ3P00UfzwgsvHPKNuO39vBtjlltrJ7ZerhHzADv5iP5EhoXw2oquabezMM/JpIyUpqQcVMoiIiLS1a677jqef/55Vq9eTVVV194rJoExbNgwPvnkE4YNG3ZQ9fZdQS06AiwhKpwTRqfx1pfbuf3kkUSGHfwMUrvKqlm/q7ypf/l35wgjPNSolEVERKQLHUxPcAluw4YNa+q1HggaMQ8CZ44fSEllHZ+s23NIx1mY5wRoUV8Onjqq5JgIipWYi4iIiAQtJeZB4MicVPrERx5yOcuCPCfJMeGM6p+wz7qU2AiNmIuIiIgEMSXmQSAsNIQzxg3kk3W7cZXXHNQxrLUszHMyPSeVkBCzz/qU2AiKlJiLiIiIBC0l5kHizPEDqXdb3v5y+0Htn7+nnF1lNczMaXtGquRYlbKIiIiIBDMl5kFiRFoCowck8NqKg5tsaMEGT335/hJzh0pZRERERIKaEvMgctb4dL7eVsr6XXs7ve+CPBeDU2IYlBLT5vqU2AhKq+qoa3AfapgiIiIi4gNKzIPIaWMHEBZiOn0TaH2Dmy8KXPt0Y2mucfbPksq6Q4pRRERERHxDiXkQSY2LZNbwPry5chsN7o7PyPrl1lLKa+r3W8YC3yXmugFUREREJDgpMQ8yZ41PZ1dZDQu8Pck7YmGeE2NgWrZjv9ukxCgxFxER6SxjTLuPuXPnBiy2jIwM7rvvvoCdX7qeZv4MMseM7EtidDivr9jK0cP6dGifBXlORg9IaBoVb0tKnBJzERGRztqxY0fT63feeYfLL7+8xbLo6OhOHa+2tpaIiP3/f92T9KbP2lU0Yh5kIsNCOXVMf/67eid7qw9cD15RU8/KzcXt1pdD81KWg+uTLiIi0hulpaU1PZKSklosq6io4OKLLyYtLY3Y2FjGjx/PO++802L/jIwM7rzzTn70ox+RlJTEnDlzAHj66acZPHgwMTExnHrqqTzyyCMY03IekrfffpsJEyYQFRVFZmYmt99+O7W1ngG2WbNmsWnTJm6++eam0fv9qa2t5dZbbyU9PZ2YmBgmTZrEf//736b18+fPxxjDxx9/zJQpU4iJiWHixImsWLGixXFyc3M5+uijiYmJYeDAgVx99dWUlZU1rZ81axZXX301P/vZz+jTpw8zZswA4N1332X48OFERUVx1FFH8fLLL2OMYePGjVRUVJCQkMCrr77a4lwffvgh4eHh7Nq1qyOXqcdQYh6EzhqfTnWdm/e+3nHAbZdsLKKuwbZbXw6Q3FTKops/RUREukJ5eTknnXQSH374IV9++SVnnXUWZ555JuvWrWux3f3338+IESNYtmwZv//971m0aBHz5s3j2muvZdWqVZx22mn86le/arHPf//7X+bMmcN1113H6tWrefrpp3n11Vf5xS9+AcDrr79Oeno6d9xxBzt27Ggxit/apZdeyqeffsqLL77IN998wyWXXMKpp57Kl19+2WK72267jT/+8Y+sWLECh8PBnDlzsNZzz9vXX3/NCSecwGmnncaXX37J66+/zqpVq/jRj37U4hgvvPAC1lo+//xz/v73v7N582bOPPNMTj75ZL788kuuv/56brnllqbtY2NjOf/883n66adbHOfpp5/mlFNOoV+/fh28Gj2DSlmC0NhBSWSlxvLaim2cO2lwu9su3OAkIiyESRkp7W4XHhpCQlSYRsxFRCTobLro4n2WxZ/0PVIuuAB3VRVbrrhyn/WJZ5xB0plnUF9czLbrf7LP+uTzzyPh+9+nbscOtt9ya4t1Q57/e5fEPWbMGMaMGdP0/vbbb+ftt9/m1Vdf5f/+7/+alh999NEtktE77riDE044gVtv9cQ1bNgwli5dypNPPtm0ze9+9ztuvvlmLr30UgCys7O5++67ufDCC7n33ntJSUkhNDSU+Ph40tLS9htjfn4+L730Ehs3bmTwYE9Ocd111/HRRx/x+OOP88gjjzRt+5vf/IbZs2c3xThz5ky2bdtGeno69957L+eeey433XRT0/aPPvoo48aNY/fu3fTt2xeAzMxM/vSnPzVtc9ttt5GVlcX9998PwPDhw1m/fj2333570zaXX345U6dOZdu2bQwcOJDi4mLefPNN/vWvfx3wGvQ0GjEPQsYYzpqQzpLCIrYUVba77YI8JxOHJBMVHnrA4zriIilSu0QREZEuUVFRwS233MKoUaNITk4mLi6OZcuWsXnz5hbbTZw4scX7devWMXny5BbLpkyZ0uL98uXL+d3vfkdcXFzT44ILLqCiooKdO3d2OMYVK1ZgrWXUqFEtjvXuu++Sn5/fYtsjjjii6fWAAQMA2L17d1M8L7zwQotjNJaqND/OhAkT9vmskyZNavezTpw4kcMPP5znnnsOgBdffJGUlBROOumkDn/OnkIj5kHq9HEDue+Db3ltxVZuOG5Ym9vs2VvDup17ufnE4R06ZnJMuEbMRUQk6LQ3gh0SHd3u+rDk5HbXh/fv32Uj5K397Gc/4/333+e+++5j6NChxMTEcPHFFzfVgTeKjY3t9LHdbje/+tWvOPvss/dZ16dPx5pDNB7HGMPSpUsJDw9vsa71javN1zfWrLvd7qbnefPm8dOf/nSfcwwcOLDp9cF8VoB58+bxwAMP8Itf/IKnn36aSy65hNDQAw869jRKzIPUwKRopmU5eH3FNn5y7NA2b+rIzfe0VDxQfXmjlNhItha3PwIvIiIiHbNgwQIuvvhizjrrLACqq6vJz89n2LC2B9QajRgxgqVLl7ZYtmTJkhbvx48fz7p168jJydnvcSIiImhoaGj3XOPGjcNay86dO5vKVA7G+PHjWb16dbvxtGXEiBH8+9//brGs9WcFmDNnDjfffDMPPfQQK1as4OWXXz7oWLszlbIEsbPGp7O5qJJlm4rbXL8wz0lCVBiHDUzs0PFSYsMprlS7RBERka4wbNgw3njjDVasWMHXX3/NhRdeSHV19QH3u/766/nggw+499572bBhA0899RRvvPFGi23uuOMOXnzxRe644w6++eYb1q1bx6uvvtqiVj0jI4PPP/+cbdu24XS2Pf/JsGHDmDNnDnPnzuXVV1+loKCAZcuWcd999/H66693+LPeeuutLFmyhKuuuoqVK1eSl5fHO++8w5VX7lv/39xVV11Ffn4+P/vZz/j22295/fXXefzxxwFaDDomJSVx9tlnc9NNN3HUUUcxdOjQDsfWkwRtYm6MiTTG/NUY4zTGVBhj3jLGpHdi/9uMMdYY85Av4/Sl7x2WRkxEKK8t37rPOmstCzY4mZ6dSmjI/lskNZcSG0lRRW3THdYiIiJy8O6//3769u3LkUceyUknncTUqVM58sgjD7jftGnTePLJJ3nwwQc54ogjePPNN7n11luJiopq2ubEE0/k3Xff5ZNPPmHy5MlMnjyZP/7xj003cAL8+te/ZsuWLWRnZ7db3vLMM89w6aWXcssttzBixAhOOeUUPvvsM4YMGdLhz3rEEUfw2WefsXHjRo4++mjGjBnDbbfddsCuKUOGDOG1117jrbfeYsyYMfz5z39u6kDT/PMCXHbZZdTW1nLZZZd1OK6exgRrkmaMeRT4AXAJ4ALuB5KACdbadn9vY4yZCrwElAGfW2uvO9D5Jk6caJctW3aoYXe5G/+5ig9X72Lp/x3X4gbPQmcFs++bz29OP4yLpnbsL9aTnxXwu/fW8tWdJ5AQFX7gHURERLrQ2rVrGTlyZKDDCEo//elP+eijj/j6668DHYrPPfDAA9xxxx2UlJS0GDV/5ZVXuPLKK9m+fTsxMTEBjLBrtPfzboxZbq2d2Hp5UI6YG2MSgcuAm621H1prVwAXAUcAx3Vg338APwLargHpRn44Pp29NfV8sKZlg/0FeZ2rLwdI9k4yVKzZP0VERALq3nvvZdWqVeTl5fHYY4/x2GOP7dMTvKd4+OGHWbJkCYWFhbz00kv85je/Ye7cuU1JeWVlJQUFBfz+97/n8ssv7xFJ+cEKysQcmACEAx80LrDWbgHWAtMPsO8TwKvW2k98F57/TM1yMCAxap9yloUbnAxMiibD0fEfXoc3MXcpMRcREQmoZcuWceKJJ3LYYYfxwAMP8Ic//IEbbrgh0GH5RF5eHmeccQYjR47kl7/8JVdddRX33ntv0/p77rmH4cOHk5KSwi9/+csARhp4wdqVJQ1oAFrfybDLu65NxpjLgRzgQt+F5l8hIYYzxg/k0fn57C6rpm9CFA1uS26+k+8dltbuFLytacRcREQkOLzyyiuBDsFv/vznP/PnP/95v+vvvPNO7rzzTv8FFMT8OmJujPmt94bM9h6zDvLYw4HfAxdYazs0i44x5gpjzDJjzLI9e/YczGn94szx6bgtvLlqGwDfbCulrLqeGZ0oYwGNmIuIiIgEM3+XsvwFGHmAxxJgJxAKtM48+3nXtWWad/vVxph6Y0w9cDRwjfd9ZOsdrLVPWGsnWmsndqZZv79l94lj3OAkXlu+zdONxVtfPj27c4l5ijcxL1JiLiIiARKsTSdEutLB/pz7tZTFWutk3/KUfRhjlgN1wPHAi95l6XgS99z97PYm0LqtyjPABjwj6d06Gz1zfDq/fPMbVm8vY2GekxFp8fSJ3+e7RrtiIkKJCAtRKYuIiAREeHg4VVVVvfrmPukdqqqq9plptSOC8uZPa20p8BRwjzHmOGPMOOB54Cvgo8btjDHrjDHXefcpsdZ+0/wBVABF3vfd+iv6qUf0JyI0hH8s3syyjcWd6sbSyBiDIzZCpSwiIhIQffv2Zdu2bVRWVmrkXHokay2VlZVs27aNvn37dnr/YL35E+AGoB54BYgGPgYubtXDfDj7lrv0SEkxERw3qi8vL92MtTBj6MF97JTYCI2Yi4hIQCQkJACwfft26uo6dDuYSLcTHh5Ov379mn7eOyNoE3NrbQ3wY+9jf9u025LEWjuri8MKqDPHpfPe1zsJDzVMzkg5qGOkaMRcREQCKCEh4aASFpHeIGgTc9nX0cP7kBoXQXafOGIjD+7SpcRGsMlV2cWRiYiIiMihUmLejYSHhvDcjyYTd5BJOUByjEpZRERERIKREvNuZvSAxEPa3xEbwd6aemrqG4gMC+2iqERERETkUAVlVxbxnZS4xtk/ddONiIiISDBRYt7LpMRokiERERGRYKTEvJfR7J8iIiIiwUmJeS/TlJhXKjEXERERCSZKzHuZpsS8vCbAkYiIiIhIc0rMe5mkmAiMUSmLiIiISLDpdGJujDnJGPOOMWaNMWaQd9k8Y8yxXR+edLXQEENSdLhKWURERESCTKcSc2PMHOCfwAYgEwj3rgoFbuna0MRXUmIjNGIuIiIiEmQ6O2J+C3C5tfanQH2z5V8AY7sqKPEtR2ykEnMRERGRINPZxHwosKiN5eVAwqGHI/6QHBuuxFxEREQkyHQ2Md8ODGtj+VFA/qGHI/6QohFzERERkaDT2cT8CeBBY8wM7/tBxphLgHuAR7s0MvGZlNhwiivrcLttoEMREREREa+wzmxsrb3HGJMIfAhEAZ8ANcB91tqHfRCf+EBKbCQNbktZdR1JMRGBDkdERERE6GRiDmCtvd0Y8ztgFJ4R9zXW2vIuj0x8xuGdZMhVUavEXERERCRIdDoxB7DWVgLLujgW8ZNkb2JeXFELfQIcjIiIiIgAnUzMjTFvtbfeWnvaoYUj/tB8xFxEREREgkNnR8xdrd6HA2OAQcDrXRKR+FxK8xFzEREREQkKnb3589K2lhtj/gSUdUlE4nMpGjEXERERCTqdbZe4P48D13bRscTHosJDiYkIVS9zERERkSDSVYn58C46jvhJckyESllEREREgkhnb/58sPUioD9wEvB0VwUlvueIi1Api4iIiEgQ6ezNn4e3eu8G9gA/RYl5t5ISG4GrXIm5iIiISLDo7M2fs30ViPhXSkwEG3ZpXigRERGRYNFVNebSzaTERujmTxEREZEgcsAR8wNNKtScJhjqPpJjI6iqa6CqtoHoiNBAhyMiIiLS63WklKX1pELSAzTO/llUWcvAiOgARyMiIiIiB0zM9zepkHRvjZMMFZXXMjBJibmIiIhIoKnGvJdKaTZiLiIiIiKB19l2iRhjZgPnA4OBiObrrLXHdFFc4mNNiXlFTYAjERERERHo5Ii5MWYu8B8gHpiFp4d5MjAeWNPFsYkPOWIjASiqqAtwJCIiIiICnS9l+RlwnbX2fKAOuM1aOw54AVBT7G4kPiqM0BCjEXMRERGRINHZxDwL+Mj7ugaI875+CJjbRTGJH4SEGJJj1MtcREREJFh0NjF34SljAdgGHOZ97QDU2qObSYkNV2IuIiIiEiQ6m5h/Dpzgff1P4EFjzDPAS8CHXRmYMSbSGPNXY4zTGFNhjHnLGJPegf36G2OeM8bsMcZUG2PWGGOO7srYegrN/ikiIiISPDrbleU6IMr7+g9APTADT5L+2y6MC+AvwA/wdIBxAfcD7xhjJlhrG9rawRiTBCwEFgAn47k5NQvY3cWx9QiO2EjW7iwLdBgiIiIiQicTc2ttUbPXbuDuLo8IMMYkApcBl1prP/QuuwjYBBwH/Hc/u94C7LDWXtxsWaEvYuwJkmPDKdaIuYiIiEhQ6Gy7xFXGmJuMMQN8FZDXBCAc+KBxgbV2C7AWmN7OfqcDi40xrxhjdnvjvc4YY3wabTeVEhtJSVUdDW4b6FBEREREer3O1pi/h6ecZZMx5iNjzCXGmLgD7XQQ0oAGwNlq+S7vuv3JAq4BCoATgQeAPwLXtrWxMeYKY8wyY8yyPXv2HHLQ3U1KTDjWQolm/xQREREJuE4l5tbaX1hrM4HZwHrgPmCXMeZlY8zJB9rfGPNbY4w9wGPWwXwQrxBghbX2NmvtSmvtM8CD7Ccxt9Y+Ya2daK2d2KdPn0M4bfeUEtc4yZAScxEREZFA6+yIOQDW2gXW2muA/sC5wHDgrQ7s+hdg5AEeS4CdQCiQ2mr/ft51+7ODfWcgXQsM7kBsvY4jNgIAlxJzERERkYDrbFeWJsaYQcAFwBxgNJ5OKO2y1jrZtzylrWMvxzOz6PHAi95l6XgS99x2dl2I50tCc8Pw3DQqrSTHeBJz3QAqIiIiEnidvfkz2VuX/SmebicX4+lhnmmt7bJe4dbaUuAp4B5jzHHGmHHA88BXfDfzKMaYdcaY65rt+mdgqjHmdmNMjjHmbOB64OGuiq0nccRpxFxEREQkWHR2xHwnnt7grwA3WGtXdn1ITW7A0yf9FTyzin4MXNyqh/lwmpW7WGuXGmNOB34P/BLY7H1+xIdxdlsaMRcREREJHp1NzE8BPvb2MPcpa20N8GPvY3/b7NMG0Vr7LvCuD0PrMSLCQoiPDNOIuYiIiEgQ6OwEQx/6KhAJjJS4CHVlEREREQkCnUrMjTEpwO+AY4G+tKpRt9YmdF1o4g/JMREUq4+5iIiISMB1tpTlKWAc8ASwHdCUkd2cIzaCHaXVgQ5DREREpNfrbGJ+LHC8tXaxL4IR/0uJjWD19rJAhyEiIiLS63V2gqHdQLkvApHASImNoKiyFmv1yw8RERGRQOpsYn478GtjTJwvghH/S4mNoLbeTUVtw4E3FhERERGf6Wwpy/8BGcBuY8wmPLNzNrHWHtFFcYmfJMd+18s8LvKgJ4IVERERkUPU2UzsVZ9EIQHjiP1u9s9BKTEBjkZERESk9+psH/O7fBWIBEaKNzEvqqgJcCQiIiIivVtna8wxxkQZY35ojLnVGJPkXZbt7XEu3cx3iXndAbYUEREREV/q7ARDOcBHQByQBPwLKAGu9r6f16XRic9pxFxEREQkOHR2xPwvwAdAP6Cq2fK3gNldFJP4UVxkGBGhIRoxFxEREQmwzt78OR2Yaq1tMMY0X74ZGNBlUYnfGGNIjg3XiLmIiIhIgHW6xhwIb2PZYKD0EGORAEmJjaSoojbQYYiIiIj0ap1NzD8Abmz23hpjEoC7gHe7LCrxq5TYcCXmIiIiIgHW2cT8RmCmMeZbIAp4BdgIpAE/79rQxF80Yi4iIiISeJ3tY77dGDMWOB8YjyexfwL4h7W2qr19JXg5YiNwKTEXERERCahOz8HuTcCf9j6kB0iOiWBvdT11DW7CQw/mtgMREREROVSd7WN+8X5WWaAayLPWrjzkqMSvUuI8vcyLK2rpmxAV4GhEREREeqfOjpg/DETg6czi9i4LARqbYIcbY1YC37PW7umaEMXXUmK8kwxVKjEXERERCZTO1i2cA6wEZuC5+TPK+3o5cAYwDjDA/V0Yo/hY0+yf5aozFxEREQmUzo6Y3w/MtdYubrZskTHmRuAZa+1IY8xNwPNdFqH4nMNbyqIbQEVEREQCp7Mj5hlAZRvLK73rAAqB5IMPSfwt2VvKUlypxFxEREQkUDqbmC8B7jfGpDUu8L6+D2gcRR8KbO2a8MQfkmM8k7m6VMoiIiIiEjCdTcznAQOAzcaYjcaYjcBm77J53m1igd92WYTic2GhISTFaPZPERERkUDq7ARDG4wxhwEnAMO9i9cBH1prrXebN7s0QvGLlJgIilTKIiIiIhIwBzPBkAX+631ID5ESG6GuLCIiIiIBdMDE3Ntx5RFrbbX39X5Za9UmsZtKjo1gS1Fb9/WKiIiIiD90ZMT8x8BzeGb2/HE721nUv7zbcsRGsGpLSaDDEBEREem1DpiYW2sz23otPUtKbATFFbVYazHGBDocERERkV6nQ11ZjDG5xpikZu//YIxJafY+1Riz2QfxiZ+kxEZQ77aUVdcHOhQRERGRXqmj7RKnAhHN3l8LJDV7Hwqkd1FMEgApsZ7Lq5aJIiIiIoHR2T7mjdqqdbCHEogElhJzERERkcA62MRcehgl5iIiIiKB1dHE3LLviLhGyHuQ7xLzmgBHIiIiItI7dXSCIQO8YIxpzNqigCeNMY2NryO7PDLxq+8S87oARyIiIiLSO3V0xPw5YDvg8j5eALY0e78d+HtXBmaMiTTG/NUY4zTGVBhj3jLGtHuDqTEm1BjzG2NMoTGm2vv8W2NMp2c47W1iIsKICg/RiLmIiIhIgHQoYbXWXurrQNrwF+AHwPl4kv/7gXeMMROstQ372edWPB1jLgG+Bo7A86WiBviNrwPu7hyxkbhUYy4iIiISEEE5kmyMSQQuAy611n7oXXYRsAk4DvjvfnadDrxtrX3b+36jMeYtYIqPQ+4RkmPDKVZiLiIiIhIQwdqVZQIQDnzQuMBauwVYiyf53p8FwGxjzAgAY8wo4BjgPd+F2nOkxEaqK4uIiIhIgATliDmQBjQAzlbLd3nX7c/dQDywxhjTgOfz/c5a+0hbGxtjrgCuABg8ePChxtztpcSEU+gsD3QYIiIiIr2SX0fMvTdi2gM8Zh3CKc4FLgYuAMZ7X19jjLmsrY2ttU9Yaydaayf26dPnEE7bM6TERlJUrhFzERERkUDw94j5X/B0dGnPZmAqEAqkAnuaresHfN7OvvcC91lrX/a+/9oYMwS4DXjqYALuTRxxEVTUNlBd10BUeGigwxERERHpVfyamFtrnexbnrIPY8xyoA44HnjRuywdGAnktrNrDJ4SmOYaCN5a+qCSHOPpZV5cWUv/xOgARyMiIiLSuwRlwmqtLcUzwn2PMeY4Y8w44HngK+Cjxu2MMeuMMdc12/Vt4OfGmJONMRnGmDOAG4E3/Bh+t9U4yZBL5SwiIiIifhesN38C3ADUA68A0cDHwMWtepgPx1Pu0ujHePqVPwL0BXYATwK/9kO83Z4j7rsRcxERERHxr6BNzK21NXgS7R+3s41p9X4vnoT+Bl/G1lM1lrKoZaKIiIiI/wVlKYsEhkOlLCIiIiIBo8RcmiRGhxNiVMoiIiIiEghKzKVJSIghOSYCl0pZRERERPxOibm0kBIboUmGRERERAJAibm0kBwbQZFKWURERET8Tom5tOCIjVBXFhEREZEAUGIuLSTHRlCsxFxERETE75SYSwuO2AiKK2txu22gQxERERHpVZSYSwspsRG4LZRU1QU6FBEREZFeRYm5tJASq9k/RURERAJBibm0oMRcREREJDCUmEsLSsxFREREAkOJubSgxFxEREQkMJSYSwvJMZ7EfNnGIkp1A6iIiIiI34QFOgAJLlHhoUzJTOH1ldt45+sdHD+yH2eOH8hRw/oQHqrvcSIiIiK+osRc9vHyFVP5ZlsZr63Yyltfbufdr3fgiI3g1DEDOGt8OocNTMAYE+gwRURERHoUY60mkgGYOHGiXbZsWaDDCDp1DW4+/XYPr6/cykdrdlPb4GZo3zjOGD+Q08cOZEBSdKBDFBEREelWjDHLrbUT91muxNxDifmBlVbW8e7XO3h9xVaWbSrGGJie7eCMcel877A04iL1CxgRERGRA1FifgBKzDtno7OCN1Zu4/WVW9lSVEW/hEjevHYG/RM1gi4iIiLSnv0l5rqbTw5KRmosPz1+GJ/dPJt/zJtCeXU9Vz6/nOq6hkCHJiIiItItKTGXQ2KMYUZOKn85bxxfbS3ltte/Rr+FEREREek8JebSJY4f1Y8bjx/GGyu38dSCwkCHIyIiItLtKDGXLnPd7BxOOiyN37+3ls/W7wl0OCIiIiLdihJz6TIhIYb7zh7DsH7xXPfiCjY6KwIdkoiIiLTBWssmVwWvLN3Mq8u3UtfgDnRIgrqyNFFXlq6zpaiSUx9aQGpcJG9cM534qPBDPuaa7WXUNbgZMyjp0AMUERHphbYUVbKowMUXBS6+yHexvbS6aV1maiy3fm84J45O0ySCfqB2iQegxLxr5eY5uejpJcwe3pcnLppASMjB/SWva3Dz1//l8fAneYQaw1NzJ3Lk0D5dHK2IiEjPs72kikX5rqZkfGtxFQApsRFMzUphWpaDadkONrkq+cN/1pG3u5yJQ5K57fsjmTAkOcDR92xKzA9AiXnXe2ZhIXe9vYbrjx3KjccP6/T++XvKufGVVXy5tZQzxw1kzY4yNrkqeWHeZCYMSfFBxCIiIt1XXYObD1bv4rP1e/ii0MUmVyUASTHhTMlsTMRTGdo3bp8Bs/oGN/9avpX7P1zPnr01fP/wNG45cQQZqbGB+Cg9nhLzA1Bi3vWstdzy6lf8a/lWHp0znpMO79/h/V5YvJnfvbuGqPBQfn/G4Xz/8P7s2VvDOY8vwllew8tXTGX0gEQffwIREZHgV1FTz8tLt/DU5wVsL60mPiqMKZme0fBpWQ5GpMV3+DfXFTX1PPl5AU98VkBtvZsLpw7h+mOHkhIb4eNP0bsoMT8AJea+UVPfwHlPfMG3O/fy2tXTGdk/od3td5dVc8trXzH/2z0cNawP9/7wCPolRDWt31ZSxdmP5lJT7+afV00ju0+crz9Cl6moqeeBjzcwJj2J7x+uGj4RETk0RRW1PJu7kb8v2khJZR2TM1K4alYWRw/rS+hBlpA22l1WzZ8/2sArSzcTGxHG1bOz+dGMTKLCQ7so+t5NifkBKDH3nV1l1Zz61wVEhofw1rUzSd7Pt+73v9nBba9/TWVtA7efPJKLpg5pM3kt2FPOOY8vIjw0hH9dNY305Bhff4RDtqWoksv/vox1O/cCcNSwPvz6tNH6FaGIBIXymnoWbHDyybrdfLp+D6VVdcRGhhITEUZMRCgxEaHERnpex0aEEd3qfURYCHUNbuoarPfZTW2Dm7r6Vu8bLHX1burdlqlZKZw9cRCJ0YfeIKC32VJUyd8+L+CVZVuornNz/Kh+XHV0tk/qwjfs2svd76/jo7W76Z8YxU0nDOfMcQMP+t4x8VBifgBKzH1r5eZizn38CyZmJPP3H00mLPS7Tp17q+u46+01vLp8K4cPTOTP544lp2/7I+Frtpdx3hOLSImN4J9XTaNvfFS72wfSksIirnphOXUNbh48fxybnBX86YP11DS4ufrobK6eld0lIxDVdQ28uXIb7369g5k5qVw4dQixkWFd8AlEJFjUN7ixQHjooXc7LnRW8L91u/lk3W4WF7qoa7DER4Zx1LA+DEiKorK2gcraBipq6j3PtfVUeZ8razzP1XVtt9gzBiJCQ4gIDSE8LITwUEN44/vQEOrdbvL3VBAdHsqZ4wcyd3oGQ/vFH/Jn8qfqugYK9njaAoeFGkJDDGEhjc8hhIYYwkNbvg8LMYeU0K7dUcbjn+bz9lc7CDHwg7EDufKoLL/82S3Kd/GH/6zlq62ljEiL53uHpTF2UBJjByWRFNO9ylyc5TWs2V7G1uIqLpgyOCAxKDE/ACXmvvfq8q387F9fcumMDH516mjAk7Te+M9VbC+p4trZOVx/7NAO/4ezYnMxF/5tMYNTYnj5iqlB+Q/DK0s3839vfsOg5BievGRiU+nN7rJqfvvuWt76cjtDHDHcddpoZg3ve1DnKK2q4x+LN/HMwo3s2VtD/8QodpRWkxwTzrwjs7h42pAuaVkpIr5V3+Bm994adpRWsaO0mh0l1Z7n0iq2l1azs7SKPXtrABiQFE2GI5bBjhgyHDEMTolliCOGIY4YYiLa/kJeW+9mSWGRJxn/djeF3rkmcvrGccyIvswe3peJGcmdSvob3JaqugZq6hoID/su8e5IGcU320p5Lncj//5yO7X1bmbkOLhkWgbHjux3yGUYvlBaVceKTcUs2VjE0sIivtpaSu1B9P6OjwpjQGI0/ZOi6J8YzYDEKPonfffcPzGqxWCNtZbFhUU89mk+87/dQ2xEKOdPHsxlR2bSPzG6Kz/iAbndlre/2s7jnxawdmcZjSlkVmosYwcnMW5QEuMGJzM8Lb5LvjweKrfbsrmokjU7yli9vZQ128tYs6OMXWWev0chBlbf9T2iI/xfnqPE/ACUmPvHr99ew9MLC/ndGYextbiKxz7NZ3BKDPefM/agfgW3MM/Jpc8sZdSABF6YN4W4QxghLqqoJTIspEtGmesb3PzuvbU8s3AjRw5N5aHzx5MYs29yvDDPyS/f/IYCZwXfPzyNX54yqsP/0O4srebphYW8uHgz5TX1HDk0lauOzmZ6toOVW0r468cb+OTbPSRGh/OjGZnMnZGhXxmL+JG1lvKaekoq6yiqqKWospbiilqKK+so9r4vKq9lZ1k1O0ur2b23Gner/5JjIkLpnxjFgKRo0hI8iZtnYphKNhVVsslVQUllXYt9+sRHNiXrGY4YEqLDyc13smCDk4raBiLCQpiW5eCYEX05ZkRfBqUEthywqKKWl5Zs5oUvNrGjtJr05GgunjaEcycObvPfTX/ZXVbdlIQv2VjMOm8iGhZiOGxgIpMzUzgiPZGwkBAa3JZ6t9v7bL97bnC3eF/fYCmqqGG79wvXjpJqXBW1+5w7JTaC/omexH1PeQ1fbinBERvBpTMyuGhqRkD/XBqV19Tz1dYSVm72PFZtKcZZ7vksUeEhHD4wkXGDkxk7KIlxg5P2+b+tMf+0FmyzZd+9hsZ3zVPVxtcW2+w1uK1ls6uyRQK+dsdeymvqAQgNMeT0iWP0gARGDUhgVH/Pc6AG9ZSYH4ASc/+ob3BzyTNLWJjnAuC8SYP45SmjDikZ/mD1Tq7+xwomZ6TwzKWTOlUW0uC2zP92Ny8u3swn3+4mNjKMC6cO4dIZGQddHlNaVcd1L67g8w1OfjQjk198f0SL0p3WauobePKzAv76vzzCQgw/PX4Yl0zP2O9ow4Zde3n8swL+vWobbgsnH96fK4/OarNLzVdbS/jr//L4cM0u4iPDmDsjgx/NyNxvnb+IHJziilpeXLKZBRucFFfWUlRRS3FlLXUNbf8fG2IgOSaCpJhw0rwJWGMi5hlJ9bxOiAo74I3ipVV1bHZVstFVweaiSjY6K5qS9saRwf6JUcwe0Zdjhvdleo5jv6PqgVTf4OaDNbt4NncjSwqLiA4P5fRxnjKX4Wm+LdWorK1no7OSr7eVsKSwmKUbi9hc5Gk1GB0eyvghSUzKSGFyRgpjByd16Z9fdV2D9zcknt+MND17E3eAC6cO5uyJg4L6xktrLVuLq1i1xZusbylm9baypt8qGNMywfal2IhQRvb/LgEfPSCRof3igurPT4n5ASgx95/iilruens1Jx8xgONH9euSY765chs//ecqjh3Rl0cvnHDAX6FtL6nin8u28MrSLeworaZPfCQ/nJDOZlcl732zg/DQEH44IZ0rjszq1A2aBXvKmffcMrYUV/Lb0w/j3Ekdr13b7KrkV299wyff7mFEWjy/Pf0wJmZ4+rVba1m2qZjH5ufz8brdRIWHcN6kwVw2M7NDo12rt5fy0P/y+M83O4mNCOWiaRnMOzKT1LjIDsfX0zW4LcWVtTjLa3Du9T6X17Cn2fvSqjrSk6MZ1i+eYf3iGNovniEpMe1+8ZL21TW4+Wz9HhbmuRjiiOGI9ERG9k/o8v9A3d6Si66+72Kjs4KnFhTyr+Wem/DGpCfSLyGKlNgIkmIiSIkNJzkmotn7CFJiIoiPCvPLzXNVtQ0UVdYyIDGqW3WCWrO9jOdyN/Lmqm3U1LuZluXgB2MH0Cc+kqSYCJJjwkmKiSAxOrzDZS8VNfVsdFWwyfslZqOzgo0uzxeZ3d4SIYDkmHAmepPwSZkpjB6QEBRlGd1RTX0Da3fsZeXmYoobfzPg/Tk0373EYJq9bra82c/sAbfFMCApmlEDEhiSEhP0N6d2u8TcGHMFcD4wDkgEMq21Gzuw31nAb4BsIB+43Vr7xoH2U2Le/b3wxSb+781vOG3MAP587th9/rFuPTpugSOH9uGCyYM5dmTfpn94C50VPPFZAa8t30q9281Jh/fn6qOzOWxg+33TP9+wh2v/sYKw0BAeu3ACkzM7PwmStZYP1uzirrdWs720mnMmpnPUsD48vaCQFZtLSImN4JJpGVw0bchB9ZT9dudeHvokj3e+2k5kWAgXThnCFUdl0Teh478dqG9w47YQEda9/qOqqW9gs6uSQmcFhc4KNroq2Frsqdl1ltdSVFGzTxkBQHioITUukj7xkcRHhbG5qJItRVVN6yPCQsjuE8ewfnEM6xfP0L6e50EpMS1+Bq217K2pZ3dZDbv3VrNnb03T693NXlfXuZme7eDYkf04cmhqj7yB11rLyi0lvLlyG+98tYOiilrCQ03T6HJ4qGF4WjxHpCcxJj2RMYOSyOkT16EvQG63ZUtxJRt2lbN+917yGp93l1NT72b84GROHN2PE0enMcRxcF2RrLUs3VjM3z4v4MO1uwgPCeEHYwcw78gsn4/s9jbFFbW8vHQLzy/a2GL6+EbGQEJUeFOi3vicFBNOfGQYO0qr2eSqpNBV0VSf3yg1LpLM1BiGODwlPxmpsQzvF092n30n3xHpat0xMb8BiAaqgD/TgcTcGDMN+Bz4FfA6cCZwFzDDWru4vX2VmPcMj87P5+7313H+5MH8/ozDMMawvaSKV5Zu4Z/LvhsdP2diOudNGtzuaPPusmqeXriRf3yxib2taribf4u31vJs7kZ+++5ahvaN48mLJx5yzWZFTT0P/m8DT31eSL3bMiglmsuPzOLsCYO65CaVvN3lPPJJHm+u2kZYaAjjByfR4LbU1rupbbDU1jdQ1+B539jmrPF1Y/KaHBPuuXEp6btfvw9IjGaA9+altMQov48yNbgt24qrKHCWs9GbgBd4k/BtxVUtEu+U2AgGJUfTJz6S1LjGRwSpzd73iYskIXrfUoKKmnrydpezftdeNjQ+7ypnW8l3CXtUeAg5feOIDAv1JOHepLu1yLAQ+iZE0jc+ir7xkRgDCzY4KauuJyIshBnZDo4b1Y9jR/QjLTF4uw91RKGzgjdXbuPNVdvY5KokMiyE40f144xxAzlyaB+c5TV8tbWEL7eW8tXWEr7aWsreak99aHR4KIcNTOCI9CSOSE9kTHoSxsD6XeVs2O3581+/ay/5e8pb/Dn3T4xiqPcLU2xEKB+v283q7WUATZ0lThydxoi0+AOOKNc3uPnPNzv52+cFfLm1lKSYcC6cMoSLpw8J6s5QPUF9g5ttJVWe+vzKWkoqaymuqKOkqs7zurLx2bu8spaK2oamevsMRywZqbFkODw3yWakxh7SPUkih6rbJeaNjDETgaV0LDF/BUix1h7fbNlHwB5r7fnt7avEvOe45/11PDI/nx9OSKe4orbd0fGOKKuu48XFm3lqQSF79tZw+MBErp6VzYmj02hwW3711je8tGQLx4/qx1/OHdulI5z5e8rZXFTJkTmpPimX2Ois4PHP8tmwq5yIME83hQhvZ4UIb4uz1ssb/+x2lXk6Rmwv8XSQKK1qeQOaMdAnLrKp20BidDiRYSFEhocS6T1WZHgIkWGh3uXfvY4ICyHEmKZWbeU19VTU1FPhfd/8deO68pp6tpdUtajpjYsMIzM1lsxUz3/KWd7nTEesT26e2ltdx4bd5WzYtZf13kSxvsF6E29v8p3gGX3vGx9Fn/jINmuI6xrcLN1YxEdrdvPh2p1NI/RHpCdy3Mh+HDeyHyP7HziR7Kzquga2FFWy0eWpT278lf+2kipSYiK8HUAau394RhkPdOOUq7yGd7/ewRsrt7FycwnGwPRsB6ePHcj3Dktrt2OQ223Z6Krgq62lrNpSwldbS1i9vYya+n2/4AxIjCKnXzzD+sYx1FtmNLRvXJvH31JUyX9X7+SD1btYuqkIa2FwSow3Se/HuEHJLUZMy2vqeXnJZp5ZuJFtJVVkOGK4bGYmZ01ID8pabfFocNug7O4iAr0nMd8M/NVae2+zZTcD11lrh7S3rxLznsNay6/eWs3fF23q8Oh4R1TXNfDGym088VkBhc4KMlNjSYwOZ9WWEq6bncONxw/r1b/+rKip97R1K6lu8dyYvJfX1FNT76amzk1NfUObZSMdEeud2MTz8Exu0vh+YFL0d8l3aiypcRHdqq62LdZaNuwu56O1u/hozS5WbinBWk8ietwoT5I+uJM/2+U19WwuqvR09nB5fqOw2VXJjrLqFjdnxUd5vtgMTIqmqKKWzUWV7GhVTpAYHd4iUR/iTdx3llbz5sptfLp+D/Vuy8j+CZwxbgCnjRl4SCP/dQ1u1u/ay9dbSwkxhqH94sjZTwLeEXv21vDR2l28/81OcvOd1DVY+sRHcsKofswe3pclG4t4afFm9tbUMzkjhXlHZgZtOz8R6T56S2JeC8yz1v692bKLgSette3e5XaEw2Hf/v7JLZbFn/Q9Ui64AHdVFVuuuHKffRLPOIOkM8+gvriYbdf/ZJ/1yeefR8L3v0/djh1sv+XWfdanXHop8cfMpqagkJ2/+tU+61OvvorY6dOpXruWXb//wz7r+/z0p8SMH0flipXs+fOf91nf7xe3ETVyJBW5uTgffWyf9Wl33UVkViZ7//cJRc88s8/6AffcTXj//pS99x7FL728z/qBDz5AWHIyJa+/Qekb+5bxD3ricUKioyl68UX2/uf9fdYPed5zmVxPPU35/Pkt1pmoKAY/+QQAex55hMpFX7RYH5qURPpfHwRg95/up2rVqpbr+/Wj+Ib/Y0T/eFx3/5GatetarI/IyKD/b34NwI5f3kHtxo0t1keOHEHaL34BwLabb6F+586mdRbYMTCH3w88hvW79vL85tfpZ1vWLsZMm0qfa64BYPPlV2CrWyYzcbNm4bjsRwBsuujiff5sesPPXuiQIRR/9D9K//4cbmuxbovbelpeVd/8S+odfYnL/YTQd94g1DtpR4gxGIL7Zy8sLY2B994DwM7f/75Lf/YA7KjDWXLceXy4dhfTnruX2OryFutX9RnKSyM8vzT8de6TRDa0/E3GkrRRvDZ0FgB3f/4I4aEhRIWHEhkeQlRYKHVHzSbp/AsYHG0o++mPaZ1+xpz2A4qPOoEthduJ/f0vqa5roLq+gWrvF653M6bzWfpYUitLuP3Ll5vKgmK8JVjB8LO3v3/36t2WdZfeyDs73dR99AHHb1iAMaapdV1cZFiv/tmLHjuWvjfdCMDWH19PQ0lJi/X6d6/3/p+rn73O/+xlvPB8m4m5X38HZ4z5LXD7ATabba2d74dwGm8wvQJgRLxu2OlJjDEcnt7+zZoHfWwgu08cb103g5p6N84b/0NDSc0B95OWwrwJYWWIAQw0K50flp5EeH8HZetiKA6i9lbBICYilHMmDeKcSYPYtKQvJTtDqG82yUn/w9KYfOYYANILEwipbfmzmTZ2AMedPp7Bjhiit75MaKvfKMSnJ5EyKAl3VRV72zh/eGgIw/rFkxXRl22tRr7dFiacOootYycRW+Kkb/F/90nsg1lYiOGE0WmcfFx/nIk72PbcV0R7S69ERPzBryPmxphUIPUAm2221lY220elLCIiIiLSY+yvlMWvI+bWWifg9OEpFgHHA/c2W3Y8kOvDc4qIiIiIHLKgvZ3cGJMGpAHDvItGGWOS8IyoF3m3+RhYYq29zbvNA8BnxpifA28CZwCzgZl+DF1EREREpNOCuXDuKmAl8A/v+3e9709rtk020L/xjbU2FzgPmAt8BVwMnHugHuYiIiIiIoEW9F1Z/EU15iIiIiLiD/urMQ/mEXMRERERkV5DibmIiIiISBBQYi4iIiIiEgRUY+5ljNkDbArQ6VPxbRtJCS663r2Lrnfvouvdu+h69z5ddc2HWGv7tF6oxDwIGGOWtXUDgPRMut69i65376Lr3bvoevc+vr7mKmUREREREQkCSsxFRERERIKAEvPg8ESgAxC/0vXuXXS9exdd795F17v38ek1V425iIiIiEgQ0Ii5iIiIiEgQUGIuIiIiIhIElJiLiIiIiAQBJeYiIiIiIkFAibmIiIiISBBQYi4iIiIiEgSUmIuIiIiIBAEl5iIiIiIiQUCJuYiIiIhIEFBiLiIiIiISBJSYi4iIiIgEASXmIiIiIiJBQIm5iIiIiEgQUGIuIiIiIhIElJiLiIiIiAQBJeYiIiIiIkFAibmIiIiISBBQYi4iIiIiEgSUmIuIiIiIBAEl5iIiIiIiQUCJuYiIiIhIEFBiLiIiIiISBJSYi4iIiIgEASXmIiIiIiJBQIm5iIiIiEgQUGIuIiIiIhIElJiLiIiIiAQBJeYiIiIiIkFAibmIiIiISBBQYi4iIiIiEgSUmIuIiIiIBAEl5iIiIiIiQUCJuYiIiIhIEAgLdADBIjU11WZkZATk3BUVFcTGxgbk3OJ/ut69i65376Lr3bvoevc+XXXNly9f7rTW9mm9XIm5V0ZGBsuWLQvIuefPn8+sWbMCcm7xP13v3kXXu3fR9e5ddL17n6665saYTW0tVymLiIiIiEgQ6DaJuTHmNmPMUmNMmTFmjzHmbWPMYa22edYYY1s9vghUzCIiIiIiHdVtEnNgFvAIMB04BqgHPjLGpLTa7iOgf7PH9/0Yo4iIiIjIQek2NebW2hObvzfGXASUAjOAt5utqrHW7vRnbCIiIiIih6o7jZi3Fo8n/uJWy2caY3YbY9YbY540xvQNQGwiIiIiIp1irLWBjuGgGGP+CQwFJlprG7zLzgMqgUIgA/gtEApMsNbWtHGMK4ArAPr16zfh5Zdf9k/wrZSXlxMXFxeQc4v/6Xr3LrrevYuud++i6937dNU1nz179nJr7cTWy7tlYm6MuR84D5hprS1oZ7sBwCbgXGvt6+0dc+LEiVbtEsUfdL17F13v3kXXu3fR9e5+6hvcfLWtlEX5LvL3lHP/OWM7tX8XtktsMzHvNjXmjYwxf8aTlM9uLykHsNZuN8ZsxTOyLiIiIiK9iNttWbOjjEX5LnLznSzdWEx5TT0AI9LiKa+pJy4yeNJhn0VijHnwIHa701pb1M4xHwDOxZOUr+tADKnAQGDHQcQiIiIiIt2ItZa83eXkehPxxYVFlFTWAZDVJ5bTxw1gWlYqU7NScMRFBjjaffnyK8J1wCKgtoPbzwT+ArSZmBtjHgYuAk4Hio0xad5V5dbacmNMHHAn8BqeRDwD+AOwG3jjYD6AiIiIiAQvay0bXZUsynexqMDFonwXznLPbYXpydGcMKof07NTmZrlIC0xKsDRHpivx+7PsNbu7siGxpi9B9jkGu/zx62W34UnIW8ADgcuBpLwJOefAOdYaw90bBEREREJctZathRVsajAyaJ8F18UFLGzrBqAvvGRzMxxMD07lWnZDgalxAQ42s7zZWJ+KZ4+4x11JbBrfyuttaa9na21VcCJ7W0jIiIiIt3L1uLvRsQXFxSxraQKgNS4CKZkOZiW5WBatoOs1FiMaTddDHo+S8yttc91cvsXfRWLiIiIiHQPNfUN/OfrnSzMc/JFoYstRZ5EPDkmnKlZDq48OotpWQ5y+sZ1+0S8Nb/fhmqMSaLVxEbt3fApIiIiIr3DjtIqrnphBV9uKSExOpwpmSn8aEYm07IdDOsbT0hIz0rEW/NLYm6MGQI8BswCIpqvAiyeSYBEREREpJdalO/iuhdXUF3XwEMXjOOkw/oT2sMT8db8NWL+DJ4bMi8DtuNJxkVERESkl7PW8tSCQv7wn3VkOGJ4/KJp5PTtnTOq+isxnwxMtdZ+46fziYiIiEiQq6ip59bXvuKdr3bwvdFp3HfOmKCa8Mff/PXJC4Hg6+IuIiIiIgFR6KzgyueXkbe7nFu/N4Krjs7qcTdzdpa/EvOfAH8wxlxjrc3z0zlFREREJAh9tGYXP31lFWGhhr//aAozh6YGOqSg4K/E/N94Rsy/NcbUAPXNV1prE/wUh4iIiIgESIPb8sBH63nwf3kcNjCBxy6cQHpy95sIyFf8lZhf56fziIiIiEgQKqms5YZXVjH/2z38cEI6vz39MKLC1ZivOb8k5p2dbEhEREREeo4128u48oVl7Cyt5renH8acKYN7fT15W3yWmBtjUhonDjLGpLS3rSYYEhEREelZrLXk7S7n43W7+ctH60mMDueVK6cxfnByoEMLWr4cMd9jjOlvrd0NOGm7d7kmGBIRERHpIbYUVZKb7yQ330Vuvos9e2sAmJbl4MHzx9EnXk362uPLxPwYoHEkfLYPzyMiIiIiAbB7bzWL8l3k5rnILXCypagKgNS4SKZnO5iR42B6diqDUnSDZ0f4LDG31n7a1msRERER6Z5Kq+r4osBFbp5nVHzD7nIAEqLCmJrlYN7MLKZnO8jpG6ca8oPg16mVvLXmfYGQ5suttWv8GYeIiIiIHFh1XQMrNhWzMN/JgjwXX28twW0hOjyUSZkp/HBCOtOzUxk1IIHQECXih8ovibkxZhzwDHB44yI8teWqMRcREREJEg1uy5rtZSzIc5Kb72RJYRE19W5CQwzjBiXx42OGMiMnlbGDkogICznwAaVT/DVi/jSwDc8MoLto+0ZQEREREfEjay2bXJVNiXhuvouSyjoAhveLZ86UIcwc6mBypoO4SL8WWvRK/voTHgqcba3N89P5RERERKQN20qqWJTv8j6cbC+tBmBAYhTHj+zHzKGpTMt20Dc+KsCR9j7+SswXACMBJeYiIiIifrS7rJpFBd5EvMDFJlclACmxEUzNSuHq7FRm5qSS4YjRDZsB5q/E/DLgb8aYLOAboK75SmvtZ36KQ0RERKRHc5XX8EVBEYsKnCzKd5G/pwKA+KgwpmQ6uGRaBtOyHQzvF0+IbtgMKv4sZRkHnNjGOt38KSIiInKQSivr+KLQMyL+RYGLdTv3AhATEcrkzBTOmThInVO6CX8l5o8DHwF/QDd/ioiIiBy0vdV1LN1Y1FSasnp7GdZCZFgIEzOS+dkJw5iWncoR6YmEh6pzSnfir8Q8Hfi+tTbfT+cTERER6REqa+tZtrG4qU78622lNLgtEaEhjBucxE+OHcq0LAdjBycRGaYihO7MX4n5h8AEQIm5iIiISDsaJ/X5osAzIr5qSwl1DZawEMOYQUlcfXQ207IdTBiSTFS4EvGexF+J+fvAn4wxRwBfs+/Nn6/7KQ4RERGRoFJT38DKzSVNNeIrN5dQ2+AmxMDh6UlcNjOLadkOJg5JJla9xHs0f13dR7zPv2hjXYdu/jTG3AacCQwHaoAvgNustd8028YAvwKuAJKBxcC11trVhxS9iIiISBeprXfz5daSpl7iKzYXU1PvxhgYPSCBS6YPYVq2g0kZKcRHhQc6XPEjvyTm1tquuPNgFp4EfylggF8DHxljRllri7zb3ALcBMwFvgXuAD40xgy31u7tghhEREREOqWuwc1XW0s9pSn5LpZtKqK6zg3AyP4JzJniScQnZ6aQGK1EvDfrNr8Psda2aLVojLkIKAVmAG97R8tvAP5orX3Nu80lwG7gAjydYURERER8qnki/kWBi2Ubi6mqawA809yfN2kwU7McTMlMITk2IsDRSjDxWWJujDkHeNNaW9vB7c8A/mutrezgKeKBEKDY+z4TSAM+aNzAWltljPkMmI4ScxEREfGBugY3X29rTMSLWLaxiMra7xLxcyamMzXLMyLuiIsMcLQSzHw5Yv4SnkR5Twe3fw4YCxR0cPsHgFXAIu/7NO/zrlbb7QIGdvCYIiIiIu2qb3DzzfYy3i2o5ZmCJSzbWESFNxEf1i+OH074LhFPVSIuneDLxNwALxhjajq4fVSHD2zM/cBMYKa1tuFggvMe5wo8N4rSr18/5s+ff7CHOiTl5eUBO7f4n65376Lr3bvoevdMbmvZutfN2iI3a10NfFvcQFW9Z92AOCdT0kIZmRLJ8ORQEiIt4ASXk29cAQ1bfMDXf8d9mZg/18nt/wGUHWgjY8yfgfOA2dba5qPrO73P/YDNzZb3a7auBWvtE8ATABMnTrSzZs3qZMhdY/78+QTq3OJ/ut69i65376Lr3TNYa8nfU+Gd0MfJonwXxZWeTs8ZjhhOH5/K9GwH7p3f8oMTZwc4WvEnX/8d91libq29tKuPaYx5ADgXT1K+rtXqQjwJ+PF4OrdgjIkCjgRu7upYREREpOfYUlTJonwXuflOcvNd7N7r+YX/gMQojhnRj+nZDqZlOxiQFN20z/zi9YEKV3qobtOVxRjzMHARcDpQbIxprCkvt9aWW2utMeYvwC+MMeuA9cD/AeXAiwEIWURERILUrrLqFon41uIqAFLjIpiW7RkRn57tYHBKDJ7GbyK+120Sc+Aa7/PHrZbfBdzpfX0PEA08zHcTDJ2gHuYiIiK9W1FFLV8UfJeIF+ypACAxOpypWSlcfqRnds2hfeOUiEvAdJvE3Fp7wL8l1lqLJ0m/09fxiIiISPAqq65jSUERud5R8XU7PWN0sRGhTM5M4fxJg5mW7WBk/wRCQ5SIS3DoNom5iIiIyP5U1zWwbGMxuflOFua7+HprCW4LkWEhTMpI4eYTBzAt28HhAxMJD+2KCclFup4ScxEREel2GmfXzM3zlKYs31xMbb2bsBDD2EFJXDc7h+k5qYwbnERkWGigwxXpEL8l5saYk4BrgSzgRGvtFmPMPKDQWtu6blxERESkidttWbdzb1ON+JLCIsprPM3ER/VP4JJpQ5iek8rkjBRiIzXuKN2TX35yjTFzgMeAvwHHAuHeVaHALex7Q6eIiIj0cluKKlmQ52RBnqeXeFFFLQCZqbH8YOwAZuSkMjXLQUpsRIAjFeka/vpKeQtwubX2Ze8oeaMvgF/7KQYREREJYkUVtZ4acW8yvqXI08KwX0Iks4b1YXqOp41h817iIj2JvxLzocCiNpaXAwl+ikFERESCSGVtPUsKPZ1TFmxwsmaHZwLw+MgwpmY7uGxGJjOHppLdRy0MpXfwV2K+HRgGbGq1/Cgg308xiIiISAB5btgsYWGei4V5TlZsLqauwRIRGsL4IUn87IRhTM9J5YiBiYSpc4r0Qv5KzJ8AHmxWxjLIGHMkngmB7vRTDCIiIuJHrW/YXFzgoqK2AYDRAxL40YxMpuekMikjmZgI3bAp4pe/Bdbae4wxicCHQBTwCVAD3GetfdgfMYiIiIhvWWvZ5KpkYb6T3DwXiwq+u2EzKzWW08cN1A2bIu3w29dTa+3txpjfAaOAEGCNtbbcX+cXERGRrrerrNp7w6aLRfkutpV4bthMS4hi1vA+TM/WDZsiHeXX3xtZayuBZf48p4iIiHSdPXtr+KLAMxr+Rb6LAmcFAEkx4UzLcnDVrGxmZDvITI3VDZsineSvPuaRwDXAbKAvnhHzJtbayf6IQ0RERDqnuKKWxYUucvM9I+Ibdnt+2R0XGcbkzBQumDKYqVkORvVPICREibjIofDXiPmTwCnAv4E1gPXTeUVERKQTyqrrWFxQxKJ8z6j4up1lWAvR4aFMykzhzPHpTMt2cNiABHVOEeli/krMTwN+YK391E/nExERkQ6ob3Dz5dYSPlvv5PMNe1i1pQS3hYiwECYOSebG44YxLdvBEelJRIQpERfxJX8l5rsBp5/OJSIiIu3Y5Krgsw1OPl+/h0X5LvbW1GMMHJGexLWzc5iencq4wUlEhYcGOlSRXsVfifkvgN8bY+Zaa4v9dE4REREBSqvqWJTv5LMNThZscLK5qBKAgUnRnDKmPzNz+jAjx0FSjFoYigSSvxLzD4Argd3GmJ1AXfOV1tosP8UhIiLS41XXNbBiUzG5+S5y851N5SmxEaFMy05l3pGZzMxJVecUkSDjr8T873j6l/8F2IVu/hQREekytfWeqe4bE/EVm0uorXcTGmI4fGAi183OYebQPowbnES4btgUCVr+SsyPB46x1i720/lERER6rAa3ZfX2Um8i7mLZxiIqaxswBkb1T+DiqUOYnuNgUkYK8VHhgQ5XRDrIX4n5ZqDGT+cSERHpUay1fLtrL7l5nkR8caGLvdX1AAztG8fZEzwtDKdkOkjWVPci3Za/EvOfAvcYY66x1ub56ZwiIiLdkrWWTa7KptKURfkuXBW1AAxOieGUI/ozNcvBtGwHfeOjAhytiHQVfyXm/wIigW+NMTVAffOV1toEP8UhIiISlHaWVpOb72yaYXNbSRUA/RIiOWpYH6ZnexLx9OSYAEcqIr7ir8T8Oj+dR0REpFsoqaxtGhHPzXNR4KwAIDkmnGnZDq6alc30bAdZ6pwi0mv4JTG31j7nj/OIiIgEq9p6Nys3F/P5Bief5zn5amsJ1tvCcEqWgwumDGZatoORaQmEhCgRF+mNfJaYG2NSrLVFja/b27ZxOxERkZ7CWkv+ngoWbNjD5xucfFHgoqK2gdAQw9hBSfzk2KEcOTSVI9LVwlBEPHw5Yr7HGNPfWrsbcNJ273LjXa45f0VEpNvbW2t5+8vtLNjg5PMNe9heWg1AhiOGM8YP5MihfZiW7SBBLQxFpA2+TMyPAYqavT7kSYWMMUcBPwMmAAOAS621zzZb/yxwSavdFltrpx7quUVERFqrrK1nSWERi/JdLMx3snpbJZaVJESFMSMnlWuPSeXInD4MduiGTRE5MJ8l5tbaT5u9nt9Fh40DvsEzk+jf97PNR8BFzd7XdtG5RUSkl6utd7NqS0nTDZsrtxRT12AJDzWMG5zM6TnhXHzCJA4fmEiYylNEpJP8cvOnMaYBaCxrab7cAey21naolMVa+x7wnnffZ/ezWY21duchhCsiIgKA221Zs6OM3HwnC/NcLG02w+bhAxP50cxMZmSnMikjheiIUObPn8+4wcmBDltEuil/tUvc3+3lkXT9iPZMY8xuoAT4FLi99RcCERGRtlhrKXBWeNoY5jlZVOCipLIOgJymGTZTmZblIDFGdeIi0rWMtYdc+r3/gxtzo/flvcBdQHmz1aHAkcAga+24gzh2OXBdqxrz84BKoBDIAH7rPc8Ea21NG8e4ArgCoF+/fhNefvnlzobRJcrLy4mLiwvIucX/dL17F13v4FdU7Watq4E1LjdrXA0U13j+X0yJMox2hDLSEcrIlBCSow5cmqLr3bvoevc+XXXNZ8+evdxaO7H1cl8n5oXel0OArUBDs9W1wEbgDmvt4oM49j6JeRvbDAA2Aedaa19v73gTJ060y5Yt62wYXWL+/PnMmjUrIOcW/9P17l10vYNPcUUtiwr2ndgnJTaCadkOpmc7mJGdyhBHTKcn9tH17l10vXufrrrmxpg2E3OflrJYazO9J/8EONNaW+zL87Vx/u3GmK3AUH+eV0REgkdj55TcfBcL85ys2VG2z8Q+07NTGZEWr4l9RCSg/DXz52x/nKc1Y0wqMBDYEYjzi4iI/9U1eDqnLMxr2TklIjSE8UOS+Olxw5iR49DEPiISdPx18yfGmHOBY4G+QIt/Ca21p3XwGHFAjvdtCDDYGDMWT7/0IuBO4DU8iXgG8AdgN/DGocYvIiLBye22rNu519s5xcniwpadUy6bmcWMHAcTh3g6p4iIBCt/tUu8F7gB+ATYzsFPNjTRe4xGd3kfzwFXA4cDFwNJeJLzT4BzrLV7D/J8IiISZKy1bC6qZGGeZ1KfRfkuiio8Db6y+sRy1vh0ZuQ4mJrlICkmIsDRioh0nL9GzC8GzrfWvnooB/FOVNReAeCJh3J8EREJTttKqliU77lh84t8V9NU9/0SIpk1rA/Tc1KZkeOgf2J0gCMVETl4/krMQ4BVfjqXiIh0c7vKqlmU7/I8ClxsLqoEIDkmnGnZDq7OcjAt20F2n7hOd04REQlW/krMnwAuxFMDLiIi0oKrvIZFBd8l4gV7PC0ME6LCmJLlYO70DKZlOxjeT51TRKTn8ldingRcYIw5HvgKqGu+0lp7vZ/iEBGRILJyczFPfl7A+9/sxO1tYTg5M4XzJg1iWlYqowYkEKpEXER6CX8l5qP4rpRlRKt1vpvhSEREgo7bbflo7S6e/LyApRuLSYwO54qjsjlxdD8OH5hImFoYikgv1aP7mIuISPCormvgtRVb+dvnhRQ6K0hPjubOU0dx9sRBxEb6rXuviEjQ0r+EIiLiU67yGp7/YhN/X7SJoopajkhP5KELxvG90WkaHRcRacZffczfam99RycYEhGR7qNgTzlPLSjk1eVbqal3c9zIvlx+ZBaTM1PUSUVEpA3+GjF3tXofDowBBgGv+ykGERHxsfKaehYXuHhl6RY+XLuL8JAQzhw/kHlHZpLTNz7Q4YmIBDV/1Zhf2tZyY8yfgDJ/xCAiIl2vrsHNqi0lLNjgZGGek1VbSqh3W5Jiwrludg4XTRtC3/ioQIcpItItBLrG/HFgAXBXgOMQEZEOsNaybudeFuZ5EvHFhUVU1jZgDBwxMJHLj8piZk4qE4YkExUeGuhwRUS6lUAn5sMDfH4RETmArcWV5Oa5WJDnJDffibO8FoCs1FjOGp/OjBwHU7McJMVEBDhSEZHuzV83fz7YehHQHzgJeNofMYiISMfsra7ji4IiFmzYw+cbnBQ4PbNwpsZFMjMnlRnex4Ck6ABHKiLSs/hrxPzwVu/dwB7gpygxFxEJqAa35autJXy+wcmCDU5WbC6m3m2JDg9lalYKc6YO4cihqQztG6duKiIiPhTwCYaMMbFAvT/iEBERjy1FlZ5EPG8PC/NclFbVYQwcNiCRK47KYuZQT514ZJjqxEVE/CVgNebGmCjgOuAWoG+g4hAR6Q1KK+tYVOBiYZ6TBXlOCr3lKf0TozhxdD+OHNqHGTmppMSqTlxEJFB8mpgbYyKAXwEnAHXAPdbaN40xFwN/BCzwZ1/GICLSG9XUN7BiUwkL8vawIM/F11tLcFuIiQhlapaDi6cN4cihfcjuE6vyFBGRIOHrEfM7gWuBD4EZwL+MMU8CxwK3AS9aa+t8HIOISI/ndnvaGDYm4ksKXVTXuQkNMYwdlMR1xwzlyKGpjElPIiIsJNDhiohIG3ydmJ8DzLXWvmGMGQOsBJKB0dZa1ZWLiBwkay2biypZlO9iYb6L3DwnrgpPG8OcvnGcN2kwM3JSmZqVQnxUeICjFRGRjvB1Yj4IWApgrf3SGFML3K2kXESk87YUVbKowMUX+S6+KHCxvbQagD7xkRw1rE9TK8O0RM20KSLSHfk6MQ8Hapq9rwNKfXxOEZEeYXtJFYvyXZ5kvMDF1uIqAFJiI5ialcLVWZ6JfXLUxlBEpEfwR1eWPxhjKr2vI4A7jTEtknNr7fV+iENEJKjtLqsmN9/FonwXXxS62OTy/NOZFBPOlMwU5s3MZFq2p594SIgScRGRnsbXiflnQHaz97nA4FbbWB/HICISlBpbGC7Kd7Iw30Xe7nIA4qPCmJLp4OJpGUzLcjAiLV6JuIhIL+DTxNxaO8uXxxcR6U6qahtYurGI3HwXuflOvtlWittCdHgokzJTOHtCOtOzUxk1IIFQJeIiIr1OwCYYEhHp6eoa3GwobuCrjzewMM/Jys0l1Da4CQsxjBucxI+PGcqMnFTGDlILQxERUWIuItJlGtyWb7aVestTXCzdWERlbQPGrGdU/wTmzshgeraDSRkpxEbqn18REWmpW/3PYIw5CvgZMAEYAFxqrX222XqDZ6bRK/D0S18MXGutXe3/aEWkp2uc1Cc338kXBS4WFxaxt9rTDTanbxxnjU8noXon8049imRNdS8iIgfQrRJzIA74Bvi799HaLcBNwFzgW+AO4ENjzHBr7V5/BSkiPZO1lrzd5SwqcJGb52JxoYviSs/kxRmOGE45oj9TsxxMy3LQN8HTS3z+fKeSchER6ZBulZhba98D3gMwxjzbfJ13tPwG4I/W2te8yy4BdgMXAI/7M1YR6f6stWxyeSb1aWxj6Cz3TM0wMCmaY0f2Y1qWg2nZDgYkRQc4WhER6e78kpgbYxqA/tba3a2WO4Dd1trQLjhNJpAGfNC4wFpbZYz5DJiOEnMR6YDGSX08ibizaXbNvvGRzMhxMD3bwbSsVAalRGtSHxER6VLGWt+3ETfGuIG0NhLzAUC+tbbTQ03GmHLgusYac2PMdGAhMMRau7nZdk8DA621J7ZxjCvw1KPTr1+/CS+//HJnw+gS5eXlxMXFBeTc4n+63sGltMayrqiBta4G1hY1sKvS829iXDiMSAlllCOUkSmhpMWag0rEdb17F13v3kXXu/fpqms+e/bs5dbaia2X+3TE3Bhzo/elBa7yJtONQoEjgXW+jKE91tongCcAJk6caGfNmhWQOObPn0+gzi3+p+sdWKVVdSz2lqbk5jtZv6sCgPjIMKZk9eGK7NQundRH17t30fXuXXS9ex9fX3Nfl7L82PtsgHlAQ7N1tcBG4KouOtdO73M/YHOz5f2arRORXqaytp6lG4vJzXeyKN/VNKlPVHgIkzJSOGNcOtOyHRw2IIGwUPUSFxGRwPH1zJ+ZAMaYT4AzrbXFPjxdIZ4E/Hhgqfe8UXhG5W/24XlFJIjU1rtZubm46WbNlVuKqWuwhIcaxg1KbprUZ8ygRCLDuuL2FhERka7hl5s/rbWzu+I4xpg4IMf7NgQYbIwZCxRZazcbY/4C/MIYsw5YD/wfUA682BXnF5Hg0+C2rN5eysI8T2nK0o1FVNe5CTFw+MBELpuZxfRsBxMzkomJ6FaNqEREpJfx2f9SxpgHgdustRXe1/tlrb2+g4edCHzS7P1d3sdzeHqX3wNEAw/z3QRDJ6iHuUjPYa2lwFlBbp6TBXlOvigoorTK00t8eL94zps0mOnZDqZkOUiMDg9wtCIiIh3ny+Gjw4HwZq8PmbV2Pp569f2tt8Cd3oeI9BA7S6tZmOdkYb6T3DwXO8s8LQwHJkXzvdFpTM9xMD07lT7xkQGOVERE5OD5LDFvXr7SVaUsItI7lFbV8UWBy5OM5znJ3+PpnJISG8G0bAczslOZkeNgcEqMeomLiEiP4a8JhiZba5fsZ92F1toX/BGHiASn6roGVmwqZmG+kwV5Lr7eWoLbQkxEKJMzUzzlKTkORqYldEkLQxERkWDkrzuh3jHGHGWtbdGz3BhzEfAYoMRcpBdpcFvWbC9jQZ6T3HwnSwqLqKl3ExpiGDcoqalzythBSUSEqYWhiIj0Dv5KzP8EfGCMmW6t3QpgjLkYeBQ4108xiEiAWGvZ5Kpkgbc0ZVGBi5LK727YnDNlCDOHOpic6SAuUp1TRESkd/JXu8S7jTF9gI+MMTOBk/Ek5Wdba9/1RwwiEhh5u8v5ycsrWb29DIABiVEcP7IfM4emMi3bQd/4qABHKCIiEhz8NjRlrf2ZMcaBp4VhGvBDa+17/jq/iPjfW19u57bXviIyPJQ7Tx3FUcP6kJkaqxs2RURE2uDLPuZntrH4PeBY4CUgqnEba+3rvopDRPyvpr6B3727lr8v2sSEIck8dME4+idGBzosERGRoObLEfNX21n3I+8DwAKaF1ukh9hSVMl1L67gy62lXH5kJrd8bwThobqBU0RE5EB82cdc/xOL9DIfr93Fjf/8Erfb8tiFE/jeYWmBDklERKTbCFj7A2NMuLW2LlDnF5GuU9/g5k8frufR+fmMHpDAI3PGM8QRG+iwREREuhW/jGobY643xpzV7P3TQJUx5ltjzHB/xCAivrG7rJo5f1vMo/PzOX/yYF67erqSchERkYPgr3KT64E9AMaYo4CzgQuAVXh6nItIN5Sb7+T7Dy7gq62l3H/OGP5w5uFEheuWERERkYPhr1KWgUCh9/WpwL+stf80xnwNfO6nGESki7jdlkc/zedPH3xLZmosL14+hWH94gMdloiISLfmr8S8DOgLbAGOB+71Lq8DNLuISDdRXdfAys0lPPFZPp98u4dTxwzgj2ceTqxm6xQRETlk/vrf9APgSWPMCiAH+I93+Wi+G0kXkSBT3+Dm622l5Oa7yM13smxjMTX1biLCQvjND0Zz4dQhmixIRESki/grMb8W+B0wGM+Mn0Xe5ePxTDYkIkHA7bZ8u2uvJxHPc7K4sIjymnoARqTFM2fKEKZnO5iclUJCVHiAoxUREelZ/JKYW2vLgB+3sfxX/ji/iLTNWkuhs4JFBS5y81wsKnBRVFELQGZqLKeNHcD0bAdTsxykxkUGOFoREZGezWeJuTEmpXFk3BiT0t62zUbQRcTHthZXkpvv4ot8F7n5LnaWVQOQlhDFrOF9mJ6dyrRsBwOTogMcqYiISO/iyxHzPcaY/tba3YATsG1sY7zL1V9NxEd2l1W3GBHfXFQJgCM2gqnZDqZlOZie7SAzNVb14iIiIgHky8T8GKBxJHy2D88jIs0UV9TyRYGr6YbN/D0VAMRHhTE1y8GlMzKYlu1gWN94QkKUiIuIiAQLnyXm1tpP23otIl1rb3UdSzcWkZvnScbX7CgDICYilEkZKZwzcRDTs1MZNSCBUCXiIiIiQcsvN38aYwYAZwHDgVrgW+Cf1tpif5xfpCeprmtgxaZicvNdLMx38tXWUhrcloiwECYMTubG44cxI8fBEelJhIf6a3JfEREROVQ+T8yNMVcCfwEi8Uw0BJAA3G+Mucxa+7J3u7HW2lW+jkeku6lrcPPV1lJy85zk5rtYvrmY2no3oSGGMemJXH10NtOzHYwfkkxUuG7XEBER6a58mpgbY74HPAz8FbjXWrvdu3wAcCvwvDFmM3ANsB5Y5ct4RLqTugY3f/pgPc8v2khFbQMAo/oncPHUIczISWVSZgpxmnFTRESkx/D1/+o3A/dZa3/efKE3Qf+JMaYK+AhwAXf6OBaRbmNHaRXXvbiS5ZuKOW3MAE46LI0pWQ5SYiMCHZqIiIj4iK8T84nAde2sfw64BZhtrc071JMZY+4EWk9atMtam3aoxxbxl8/W7+GGV1ZRU9fAg+eP47QxAwIdkoiIiPiBrxPzcKCqnfVVQHVXJOXNfAvMava+oQuPLeIzDW7LAx+t56+f5DG8XzwPzxlPdp+4QIclIiIifuLrxPxb4ATgif2sP9G7TVeqt9bu7OJjivjU7r3V/OSlVSwqcHH2hHR+/YPDiI7QjZwiIiK9ia8T86eBPxpjdlhr326+whhzGvAH4I4uPmeWMWY7UAMsBn5hrS3o4nOIdJlF+S6uf3kle6vruOeHR3DOxEGBDklEREQCwNeJ+UPADODfxpj1wFrv8pHAUDy9zB/qwvMtBuYC64C+wP8BucaY0dZaVxeeR+SQud2WRz/N508ffEtGaizPXzaZEWkJgQ5LREREAsRYa31/EmPOAS4AhnkXrQf+Ya39l4/PGwcUAH+01t7fxvorgCsA+vXrN+Hll1/2ZTj7VV5eTlycaol7i/LycoiI5YmvavjK2cCUtFDmHhZJdJhm5eyJ9Pe7d9H17l10vXufrrrms2fPXm6tndh6uV8S80AyxnwCrLPWXt3edhMnTrTLli3zU1QtzZ8/n1mzZgXk3OJ/f3vjY55eB87yWn556igunDIYY5SU91T6+9276Hr3LrrevU9XXXNjTJuJeY+encQYEwWMAD4JdCzSu1lr2VxUydtfbufPS6oZkBzNa1dP5/D0xECHJiIiIkGiRyXmxpj7gLeBzXhqzH8JxOLply7iVztKq8jNc5Gb7+KLAhfbSjydQ8f1DeXZq44kMSY8wBGKiIhIMOlRiTmQDrwEpAJ7gC+AqdbaTQGNSnoFZ3kNi/K/S8QLnRUAJMWEMy3LwZVHZzEty8HWNcuUlIuIiMg+elRibq09L9AxSO9RWlnHF4UuFuV7Ht/u2gtAXGQYUzJTmDNlMNOyHYxMSyAk5Lsa8m1rVU8uIiIi++pRibmIL1XXNbBiUzEL8pwszHfx9dYS3BaiwkOYlJHCD8YNYHp2KocNSCAsNCTQ4YqIiEg347PE3BjzYEe3tdZe76s4RA5Wg9uyZnsZC/Kc5OY7WVJYRE29m9AQw7hBSVx3zFBmZDsYOziJyDDN0ikiIiKHxpcj5od3cLue3a9Rug1rLZtclU2JeG6+i5LKOgCG94vngimDmZmTyuTMFOKjVCMuIiIiXctnibm1dravji3SVXbvrWZRvouFeU4W5n3XOWVAYhTHj+zHjJxUpmc76JsQFeBIRUREpKdTjbn0KmXVdSwuKGKhd1R8/a5yABKjPZ1Trjo6ixk5qWSmxmrSHxEREfErvyXmxphhwA+BwUBE83XW2h/5Kw7pXarrGlixuZjcPBcL8518tbWUBrdtumHzjHHpzMhxMHpAIqEhSsRFREQkcPySmBtjTgZeA1YCE4ClQDYQCXzujxikd6hvcPPN9jJy850syne1uGHziPRErpmVzfTsVMYP0Q2bIiIiElz8NWL+a+Aua+0fjDF7gYuA7cDzwCI/xSA9kNttWbdzb4tEfG9NPQDD+sVxwZTBzMhOZUqWbtgUERGR4OavxHw48Ir3dR0QY62tNsb8GngXuN9PcUg3Z60lf08Fi7xdU74ocFHs7ZySmRrLqWMHMD3bwdQsB6lxkQGOVkRERKTj/JWY7wUa21rsAHKAb7znT/ZTDNJNbSupYuEGJwu9yfievTUADEyK5tiR/Zie7WBatoP+idEBjlRERETk4PkrMV8MzATW4Bkh/5MxZgxwBiplkVbKquuaWhgu2OCkwFkBQJ/4SKZlOZie7WB6diqDUqLVOUVERER6DH8l5jcCcd7XdwLxwFnAeu866cVq692s2lLCgg17WJDn5Etv55SYiFCmZKYwZ+oQjhyaytC+cUrERUREpMfyS2JurS1o9roSuNof55XgZK0lb3c5n29wsiDPyeICFxW1DYQYGDMoiWtmZTMzJ5Vxg5OJCAsJdLgiIiIifuGvdokFwCRrravV8iRghbU2yx9xSOBsLa5s6iXevE48MzWWM8enMyMnlWnZDhKj1TlFREREeid/lbJkAG01jY4EBvopBvEjZ3kNi/Jd5OZ7prrfXFQJQGpcJNOzHczIcTAjJ5X05JgARyoiIiISHHyamBtjzmz29mRjTGmz96HAscBGX8Yg/rG3uo4lhUUszPMk4+t27gUgPiqMqVkOLp2RwYwc1YmLiIiI7I+vR8xf9T5b4KlW6+rwJOU3+TgG8YHmU93n5n93w2ZkmGeq+5tPHMCMnFQOG5BAWKjqxEVEREQOxKeJubU2BMAYU4inxtzpy/OJ79Q3uPl6Wym53vKUZRuLm6a6H5OeyNVHZzM9x8H4wclEhWuqexEREZHO8ldXlkx/nEe6jrWW9bvKWZjnuVlzcYGraar7EWnxXDh1CNOzHUzO1FT3IiIiIl3BXzd/Yow5GbgVGIWntGUNcLe19j1/xSDt21pc6ZnUJ8/FonwnzvJaAIY4YjhlzICmGTY11b2IiIhI1/NXu8R5wCPAP4DnvIuPBN4wxlxtrX3aH3FIS8UVtSwqcLEgz0lunpONru86p8zISWVGdirTcxzqnCIiIiLiB/4aMb8VuNFa+1CzZU8ZY5YDPweUmPtBdV0DSzd6OqcszHPyzfZSrIW4yDCmZKZw8bQMZmqGTREREZGA8FdiPhh4v43l/wHu81MMvU59g5tvtpexMM/JwjwnyzYVU1vvJjzUMG5wMjccO4yZQx0ckZ5EuDqniIiIiASUvxLzzcDxQF6r5ScAm/wUQ4/XeMNm46Q+zW/YHNk/gUumDWF6TiqTM1KIjfTb7QUiIiIi0gG+nmDoaeAneEbF/2qMGQ/kelfPAC4CfuzLGHq6LUWVTZ1TcvNdOMs9U93rhk0RERGR7sXXw6aXAD+31j5ujNmNZzKhxtlA1wLnWGv/7eMYepQ9e2vIzXd6JvYpcLKlqAqAPvGRzMxxMF03bIqIiIh0S75OzJvuILTWvgG84ePzeU5qzDXAzUB/YDVwg7X2c3+cu6uVVdexuKDIOyruZP2ucsAz1f20LAeXzchkRk4qObphU0RERKRb80ehsfXDOZoYY84FHgCuARZ4n/9jjBllrd3sz1gORnVdA8s2FnvqxPNdfL21BLeFqHDPVPdnjEtneraDwwYmEhqiRFxERESkp/BHYr7zQCO51tqunMP9RuBZa+2T3vc/NsZ8D7gauK0Lz9Ml6hvc5BU38M3/NrAwz8XyzZ7OKaEhhrGDkrhudg7Tc1IZNziJyDBNdS8iIiLSU/kjMb8CKPHDeTDGRAAT2LcF4wfAdH/E0Blf/H97dx4rV1mHcfz7FMSlBYJUCtWAIIooaIViQJYWEhTBRHCJ4gYmgAa3GFEDAawmGKOIQIAg/mOjYl1BqzEqSoOCCy1BqArIorLIXqhXa0vh5x/nFIex1Sp37kzPfD/JZGbe98w5b/tkbn9973vOufUBjlu4lInVa4Gb2H2HrXjHvjux/64z2WfnZzLDK6dIkiSNjVQNbqVJkseA7avq3oEd5InHmw3cCcyrqit62k8H3lpVu/VtfwLNfxyYNWvW3osWLZqKYT7u4dXFJTevYZfpj/Cy2dPZcguXpoyDiYkJZsyYMexhaIqY93gx7/Fi3uNnsjI/+OCDl1XV3P72QU/JTun68v9VVV0EXAQwd+7cmj9//pSP4bWvgiVLljCMY2s4zHu8mPd4Me/xYt7jZ9CZD/p2j1M9BXw/8Cgwq699FnD3FI9FkiRJ2mgDLcyratpULWNpj7cGWEZzl9Feh/KvGxtJkiRJI6eLZxeeBXwpya+BK4F3A7OBC4c6KkmSJOk/6FxhXlVfS7ItcCrNDYaWA4dX1Z+GOzJJkiRpwzpXmANU1QXABcMehyRJkrSxBn3ypyRJkqSNYGEuSZIkjYCB3mBoU5LkPmBY69Bn0lzqUePBvMeLeY8X8x4v5j1+JivznarqWf2NFuYjIMnS9d39Sd1k3uPFvMeLeY8X8x4/g87cpSySJEnSCLAwlyRJkkaAhflouGjYA9CUMu/xYt7jxbzHi3mPn4Fm7hpzSZIkaQQ4Yy5JkiSNAAtzSZIkaQRYmA9RkhOT3JbkH0mWJTlw2GPSk5fkoCTfTXJnkkpybF9/kixIcleSVUmWJHnxkIarJynJyUmuTrIyyX1JFifZo28bM++IJO9Jcl2b98okv0hyRE+/WXdY+32vJOf1tJl5h7RZVt/j7p7+geZtYT4kSd4EnAN8EngZcBXwgyQ7DnVgmgwzgOXAB4BV6+n/CPAh4H3APsC9wI+TbDllI9Rkmg9cALwCOARYC1yW5Jk925h5d9wBfBTYC5gL/BS4NMlL2n6z7qgk+wInANf1dZl599wI7NDz2LOnb6B5e/LnkCT5FXBdVR3f0/YH4JtVdfLwRqbJlGQCeG9VfbF9H+Au4LyqOqNtezrNF/ukqvr8sMaqyZFkBvAwcGRVLTbz7kvyIHAyzdUazLqDkmwNXAMcB3wMWF5V7/X73T1JFgBvqKo91tM38LydMR+CJFsAewM/6uv6Ec2sm7prZ2B7erKvqlXAFZh9V2xJ87N1RfvezDsqyWZJ3kzzW7KrMOsuu4hm4uzyvnYz76Zd2qUqtyVZlGSXtn3geVuYD8dMYDPgnr72e2gCV3ety9fsu+sc4FrgF+17M++YJHu2vw1bDVwIHFVV12PWnZTkeGBX4NT1dJt59/wKOBY4DDieJserkmzLFOS9+WTsRJIESc4CDgAOqKpHhz0eDcyNwBxga+ANwMIk84c4Hg1Ikt1ozgU7oKoeGfZ4NHhV9YPe90l+CdwKHAP8ctDHd8Z8OO4HHgVm9bXPAu7+983VIevyNfuOSfI54GjgkKq6tafLzDumqtZU1c1Vtaw9J+ha4IOYdRftR/Nb7t8mWZtkLTAPOLF9/UC7nZl3VFVNAL8Fns8UfMctzIegqtYAy4BD+7oOpVmnqO66jebL+3j2SZ4GHIjZb7KSnMO/ivIb+rrNvPumAU/FrLvoUporcszpeSwFFrWvb8LMO63N84XAX5iC77hLWYbnLOBLSX4NXAm8G5hNs15Rm7D2qhy7tm+nATsmmQM8WFV/TnI2cEqSG2h+qJ8KTAAXD2G4epKSnA+8HTgSWJFk3TrDiaqaqKoy8+5I8ing+8DtNCf6voXmkplHmHX3VNVDwEO9bUn+RvPzfHn7/mzMvDOSnAksBv4MbAecBkwHFk7Fd9zCfEiq6mvtiQSn0lwjczlweFX9abgj0ySYC/Seuf/x9rGQ5oSSTwNPB84HtqE50eSVVfXXqR2mJsmJ7fNP+to/DixoX5t5d2wPfLl9fpjmmtavrqoftv1mPX7MvFueA3yVZgnTfTTryvftqc8GmrfXMZckSZJGgGvMJUmSpBFgYS5JkiSNAAtzSZIkaQRYmEuSJEkjwMJckiRJGgEW5pIkSdIIsDCXJA1MkgVJlg97HJK0KfA65pLUEUm+CMysqtf0vp6iYz+X5nbV+1TV0p72GcBTq+qBqRiHJG3KvPOnJGmDkmwOPFr/5yxOVU3Q3K5akvRfuJRFkjomyQLgGOCIJNU+5rd9z06yKMmK9vH9JM/v/WyS5UmOTXILsBqYnuSwJD9rP/Ngkh8m2b3nsLe1z1e3x1vSu7+e/U9LclqS25OsTnJ9ktf29D+3/fzrk/w4yd+T/C7JoYP525Kk0WFhLkndcybwdeAyYIf2cVWSZwCXA/8A5gH7AX8BLmv71tkZeAvwRuCl7fbTgbOBlwPzgYeBxUm2aD/z8vb5sPZ4r9vA2D4AfBj4KLAncAnw7SRz+rY7Azi3Pf7VwKJ2WYwkdZZLWSSpY6pqIskqYHVV3b2uPcnbgADvXLc0Jcm7gHuB19AU8wBbAG+vqnt6dvut3mMkeSewkqYg/zlwX9v1QO8x1+Mk4Myqurh9f3qSg9r2t/Vs97mqWtwe6xTgHcCc9liS1EnOmEvS+NibZjb8r0kmkkzQzHxvAzyvZ7s7+opykjwvycVJbkmyEriH5t+QHTf24Em2AmYDV/Z1/Rx4UV/bdT2v72qft9vYY0nSpsgZc0kaH9OAa4E3r6fvwZ7Xf1tP//eAO4B3AXcCa4Hf0cyuT4b+k0sfebyjqpKAk0mSOs7CXJK6aQ2wWV/bNcDRwP1V9dDG7ijJtsALgROr6vK2bS+e+G/Imva5/5iPq6qVSe4C9gd+0tN1AE2RL0ljzdkHSeqmPwJ7JNktycwkTwG+QrME5TtJ5iXZOclBST7be2WW9VgB3A8cn2TXJPOAC2lmzde5F1gFvCrJrCRbb2BfnwFOSnJ0khck+QRwIM0Jq5I01izMJambvgD8HlhKc2Lm/lX1d+Ag4FbgG8ANwEKaNeYrNrSjqnoMeBPwEmA5cD5wGs2lFNdtsxZ4P3AczZrw72xgd+fSFOefbvd1FPD6qvrN//nnlKTO8M6fkiRJ0ghwxlySJEkaARbmkiRJ0giwMJckSZJGgIW5JEmSNAIszCVJkqQRYGEuSZIkjQALc0mSJGkEWJhLkiRJI8DCXJIkSRoB/wRiBw/HHU28pQAAAABJRU5ErkJggg==\n", + "text/plain": [ + "
" + ] + }, + "metadata": { + "needs_background": "light" + }, + "output_type": "display_data" + } + ], + "source": [ + "import matplotlib.pyplot as plt \n", + "\n", + "plt.rcParams['font.size'] = 14\n", + "\n", + "fig, (ax1, ax2) = plt.subplots(2, 1, sharex=True, figsize=(12,8))\n", + "\n", + "# plot loss and reference value\n", + "ax1.plot(loss, label='Qiskit Runtime VQE')\n", + "ax1.axhline(y=target_energy, color='tab:red', ls='--', label='Target energy')\n", + "\n", + "# plot time taken\n", + "ax2.plot(np.cumsum(runtimes_in_min))\n", + "\n", + "# settings\n", + "ax1.set_title('Qiskit Runtime VQE')\n", + "ax1.set_ylabel('Eigenvalue')\n", + "ax1.legend(loc='best')\n", + "ax2.set_ylabel('Total Qiskit Runtime [min]')\n", + "ax2.set_xlabel('Iteration')\n", + "ax2.grid();\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Qiskit Runtime VQE: Direct call\n", + "\n", + "Instead of interacting with the Qiskit Runtime via the `VQEProgram` we can also directly call the the Qiskit Program in the cloud -- like running the Quantum Kernel Alignment (QKA) or Circuit Runner. Under the hood, `VQEProgram` is nothing but a proxy that takes care of creating the dictionaries interaction with the Qiskit Program running in the Qiskit Runtime on the cloud." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "To call the VQE code on the Qiskit Runtime on the cloud, we interact with the `run` method of the provider\n", + "```\n", + "provider.run(progam_id, inputs, options, callback)\n", + "```\n", + "where \n", + "* `program_id` determines the Qiskit Runtime program run, i.e. `'vqe'` for the VQE\n", + "* `inputs` contains the input for the VQE algorithm\n", + "* `options` specifies the backend\n", + "* `callback` is the callback used inside the VQE" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Let's gather all the inputs to the VQE in a dictionary." + ] + }, + { + "cell_type": "code", + "execution_count": 21, + "metadata": {}, + "outputs": [], + "source": [ + "vqe_inputs = {\n", + " 'ansatz': ansatz,\n", + " 'operator': hamiltonian,\n", + " 'optimizer': {'name': 'SPSA', 'maxiter': 5}, # let's only do a few iterations!\n", + " 'initial_point': initial_point,\n", + " 'readout_error_mitigation': True,\n", + " 'shots': 1024\n", + "}" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The backend options only need to contain the name of the backend." + ] + }, + { + "cell_type": "code", + "execution_count": 22, + "metadata": {}, + "outputs": [], + "source": [ + "backend_options = {\n", + " 'backend_name': backend.name()\n", + "}" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The callback acts a little different than before. Instead of only being passed the VQE callback arguments, it is passed a tuple with the job ID as first argument, followed by the VQE arguments." + ] + }, + { + "cell_type": "code", + "execution_count": 23, + "metadata": {}, + "outputs": [], + "source": [ + "intermediate_info = {\n", + " 'nfev': [],\n", + " 'parameters': [],\n", + " 'energy': [],\n", + " 'stddev': []\n", + "}\n", + "\n", + "def raw_callback(*args):\n", + " job_id, (nfev, parameters, energy, stddev) = args\n", + " intermediate_info['nfev'].append(nfev)\n", + " intermediate_info['parameters'].append(parameters)\n", + " intermediate_info['energy'].append(energy)\n", + " intermediate_info['stddev'].append(stddev)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now we can call the Qiskit Runtime VQE directly, without going through the convenience of the `VQEProgram`. \n", + "Note, that the result will not be of the same type as `VQE` or `VQEProgram`, but a plain dictionary." + ] + }, + { + "cell_type": "code", + "execution_count": 24, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Job ID: c2aga1418tu3tem8ghn0\n" + ] + } + ], + "source": [ + "job = provider.runtime.run(\n", + " program_id='vqe',\n", + " inputs=vqe_inputs,\n", + " options=backend_options,\n", + " callback=raw_callback\n", + ")\n", + "print('Job ID:', job.job_id())" + ] + }, + { + "cell_type": "code", + "execution_count": 25, + "metadata": {}, + "outputs": [], + "source": [ + "result = job.result()" + ] + }, + { + "cell_type": "code", + "execution_count": 30, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Reached -0.41066745294814644 after 16 evaluations.\n" + ] + } + ], + "source": [ + "print(f'Reached {result[\"optimal_value\"]} after {result[\"optimizer_evals\"]} evaluations.')" + ] + }, + { + "cell_type": "code", + "execution_count": 31, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Available keys: ['optimizer_evals', 'optimizer_time', 'optimal_value', 'optimal_point', 'optimal_parameters', 'cost_function_evals', 'eigenstate', 'eigenvalue', 'aux_operator_eigenvalues', 'optimizer_history']\n" + ] + } + ], + "source": [ + "print('Available keys:', list(result.keys()))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.7.9" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +}