diff --git a/README.md b/README.md index aef3a01..c01e689 100644 --- a/README.md +++ b/README.md @@ -5,8 +5,8 @@ * [Installation](#installation) * [Quick start](#quick-start) - + [Pure state simulation](#pure-state-simulations) - + [Mixed state simulation](#mixed-state-simulations) + + [Pure state simulation](#pure-state-simulation) + + [Mixed state simulation](#mixed-state-simulation) * [Converting from TKET](#converting-from-tket) * [Examples](#examples) * [Contributing](#contributing) @@ -131,13 +131,15 @@ An example of this can be found in the [`pytket-qujax_heisenberg_vqe.ipynb`](htt Below are some use-case notebooks. These both illustrate the flexibility of qujax and the power of directly interfacing with JAX and its package ecosystem. -- [`heisenberg_vqe.ipynb`](https://github.com/CQCL/qujax/blob/main/examples/heisenberg_vqe.ipynb) - an implementation of the variational quantum eigensolver to find the ground state of a quantum Hamiltonian. -- [`maxcut_vqe.ipynb`](https://github.com/CQCL/qujax/blob/main/examples/maxcut_vqe.ipynb) - an implementation of the variational quantum eigensolver to solve a MaxCut problem. Trains with Adam via [`optax`](https://github.com/deepmind/optax) and uses more realistic stochastic parameter shift gradients. -- [`noise_channel.ipynb`](https://github.com/CQCL/qujax/blob/main/examples/noise_channel.ipynb) - uses the densitytensor simulator to fit the parameters of a depolarising noise channel. -- [`qaoa.ipynb`](https://github.com/CQCL/qujax/blob/main/examples/qaoa.ipynb) - uses a problem-inspired QAOA ansatz to find the ground state of a quantum Hamiltonian. Demonstrates how to encode more sophisticated parameters that control multiple gates. -- [`variational_inference.ipynb`](https://github.com/CQCL/qujax/blob/main/examples/variational_inference.ipynb) - uses a parameterised quantum circuit as a variational distribution to fit to a target probability mass function. Uses Adam via [`optax`](https://github.com/deepmind/optax) to minimise the KL divergence between circuit and target distributions. -- [`classification.ipynb`](https://github.com/CQCL/qujax/blob/main/examples/classification.ipynb) - train a quantum circuit for binary classification using data re-uploading. -- [`generative_modelling.ipynb`](https://github.com/CQCL/qujax/blob/main/examples/generative_modelling.ipynb) - uses a parameterised quantum circuit as a generative model for a real life dataset. Trains via stochastic gradient Langevin dynamics on the maximum mean discrepancy between statetensor and dataset. +- [`heisenberg_vqe.ipynb`](https://github.com/CQCL/qujax/blob/develop/examples/heisenberg_vqe.ipynb) - an implementation of the variational quantum eigensolver to find the ground state of a quantum Hamiltonian. +- [`maxcut_vqe.ipynb`](https://github.com/CQCL/qujax/blob/develop/examples/maxcut_vqe.ipynb) - an implementation of the variational quantum eigensolver to solve a MaxCut problem. Trains with Adam via [`optax`](https://github.com/deepmind/optax) and uses more realistic stochastic parameter shift gradients. +- [`noise_channel.ipynb`](https://github.com/CQCL/qujax/blob/develop/examples/noise_channel.ipynb) - uses the densitytensor simulator to fit the parameters of a depolarising noise channel. +- [`qaoa.ipynb`](https://github.com/CQCL/qujax/blob/develop/examples/qaoa.ipynb) - uses a problem-inspired QAOA ansatz to find the ground state of a quantum Hamiltonian. Demonstrates how to encode more sophisticated parameters that control multiple gates. +- [`barren_plateaus.ipynb`](https://github.com/CQCL/qujax/blob/develop/examples/barren_plateaus.ipynb) - illustrates how to sample gradients of a cost function to identify the presence of barren plateaus. Uses batched/vectorized evaluation to speed up computation. +- [`reducing_jit_compilation_time.ipynb`](https://github.com/CQCL/qujax/blob/develop/examples/reducing_jit_compilation_time.ipynb) - explains how JAX compilation works and how that can lead to excessive compilation times when executing quantum circuits. Presents a solution for the case of circuits with a repeating structure. +- [`variational_inference.ipynb`](https://github.com/CQCL/qujax/blob/develop/examples/variational_inference.ipynb) - uses a parameterised quantum circuit as a variational distribution to fit to a target probability mass function. Uses Adam via [`optax`](https://github.com/deepmind/optax) to minimise the KL divergence between circuit and target distributions. +- [`classification.ipynb`](https://github.com/CQCL/qujax/blob/develop/examples/classification.ipynb) - train a quantum circuit for binary classification using data re-uploading. +- [`generative_modelling.ipynb`](https://github.com/CQCL/qujax/blob/develop/examples/generative_modelling.ipynb) - uses a parameterised quantum circuit as a generative model for a real life dataset. Trains via stochastic gradient Langevin dynamics on the maximum mean discrepancy between statetensor and dataset. The [`pytket`](https://github.com/CQCL/pytket) repository also contains `tk_to_qujax` implementations for some of the above at [`pytket-qujax_classification.ipynb`](https://github.com/CQCL/pytket/blob/main/examples/pytket-qujax-classification.ipynb), [`pytket-qujax_heisenberg_vqe.ipynb`](https://github.com/CQCL/pytket/blob/main/examples/pytket-qujax_heisenberg_vqe.ipynb) diff --git a/examples/barren_plateaus.ipynb b/examples/barren_plateaus.ipynb new file mode 100644 index 0000000..9961d90 --- /dev/null +++ b/examples/barren_plateaus.ipynb @@ -0,0 +1,412 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "d2fb25e7", + "metadata": {}, + "source": [ + "# Barren plateaus and vanishing gradients in parameterized quantum circuits" + ] + }, + { + "cell_type": "markdown", + "id": "05a9bb56", + "metadata": {}, + "source": [ + "## Definition" + ] + }, + { + "cell_type": "markdown", + "id": "1f5dd510", + "metadata": {}, + "source": [ + "In this notebook, we will numerically study the phenomenon of _barren plateaus_ in parameterized quantum circuits (PQCs).\n", + "\n", + "Given a parameterized quantum circuit $U(\\theta)$ and an initial state $\\ket{\\psi(0)}$, define \n", + "\n", + "$$\\ket{\\psi(\\theta)} := U(\\theta)\\ket{\\psi(0)}.$$ \n", + "\n", + "If $c$ is some cost function, $c(\\ket{\\psi(\\theta)})$ has a barren plateau if for all indices $j$ and some $k > 0$,\n", + "$$\n", + "\\frac{\\partial c}{\\partial \\theta_j} = O(\\exp(-k L))\n", + "$$\n", + "where $L$ is the number of qubits. \n", + "\n", + "In practical terms, this means that training of this parameterised circuit does not scale favourably with the number of qubits. \n", + "\n", + "Barren plateaus were first described [here](https://doi.org/10.1038/s41467-018-07090-4) and have since been extensively studied in the literature (see e.g. Section 6.1 of [this manuscript](https://doi.org/10.1016/j.physrep.2022.08.003) and the citations therein for a more in-depth discussion)." + ] + }, + { + "cell_type": "markdown", + "id": "09de3bc1", + "metadata": {}, + "source": [ + "## Numerical characterization of vainishing gradients" + ] + }, + { + "cell_type": "markdown", + "id": "c0aa4342", + "metadata": {}, + "source": [ + "It follows from a direct application of [Chebyshev's inequality](https://en.wikipedia.org/wiki/Chebyshev%27s_inequality) that\n", + "$$\n", + "P\\left(\\left| \\frac{\\partial c}{\\partial \\theta_j} \\right| \\geq \\epsilon \\right) \\leq \\frac{\\text{Var}\\left(\\frac{\\partial c}{\\partial \\theta_j}\\right)} {\\epsilon^2} \n", + "$$\n", + "Thus, if it is the case that $\\text{Var}\\left(\\frac{\\partial c}{\\partial \\theta_j}\\right) = O(\\exp(-kL))$ for some $k>0$, we are in the presence of a barren plateau. It has been empirically observed that this scaling is independent of the index $j$ chosen; here, we will choose the first angle in the vector of parameters $\\theta$." + ] + }, + { + "cell_type": "markdown", + "id": "6fd3a862", + "metadata": {}, + "source": [ + "## Imports" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "id": "16f5bf05-00b5-42a6-934e-2ab4402c6a49", + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "No GPU/TPU found, falling back to CPU. (Set TF_CPP_MIN_LOG_LEVEL=0 and rerun for more info.)\n" + ] + } + ], + "source": [ + "from typing import List, Tuple\n", + "\n", + "import jax\n", + "import jax.numpy as jnp\n", + "from jax.random import PRNGKey\n", + "\n", + "import qujax\n", + "from qujax import all_zeros_statetensor, print_circuit, repeat_circuit" + ] + }, + { + "cell_type": "markdown", + "id": "17c06a73", + "metadata": {}, + "source": [ + "## Define and check ansatz" + ] + }, + { + "cell_type": "markdown", + "id": "5389a6b8", + "metadata": {}, + "source": [ + "Here, the PQC will be a hardware efficient ansatz, which is known to have barren plateaus due to its expressibility" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "id": "397640bd", + "metadata": {}, + "outputs": [], + "source": [ + "def hardware_efficient_ansatz(\n", + " n_qubits: int,\n", + ") -> Tuple[List[str], List[List[int]], List[List[int]], int]:\n", + " \"\"\"\n", + " Builds and returns the circuit specification for a hardware efficient ansatz\n", + " \"\"\"\n", + "\n", + " gates: List[str] = []\n", + " qubit_inds: List[List[int]] = []\n", + " param_inds: List[List[int]] = []\n", + "\n", + " parameter_index = 0\n", + "\n", + " # Ry layer\n", + " for i in range(n_qubits):\n", + " gates.append(\"Ry\")\n", + " qubit_inds.append([i])\n", + " param_inds.append([parameter_index])\n", + " parameter_index += 1\n", + "\n", + " # Rx layer\n", + " for i in range(0, n_qubits):\n", + " gates.append(\"Rx\")\n", + " qubit_inds.append([i])\n", + " param_inds.append([parameter_index])\n", + " parameter_index += 1\n", + "\n", + " # CRz layer\n", + " for i in range(n_qubits - 1):\n", + " gates.append(\"CZ\")\n", + " qubit_inds.append([i, i + 1])\n", + " param_inds.append([])\n", + "\n", + " return gates, qubit_inds, param_inds, parameter_index" + ] + }, + { + "cell_type": "markdown", + "id": "1111e958", + "metadata": {}, + "source": [ + "Print one repetition of the circuit to visually check for correctness" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "id": "95e6e9bc-8ab0-4411-8e34-7e5cbdeb77d9", + "metadata": {}, + "outputs": [], + "source": [ + "gates, qubit_inds, param_inds, nr_of_parameters = hardware_efficient_ansatz(4)" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "id": "7d0f4aee-aad5-4b62-8483-5073746374e0", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "q0: ---Ry[0]---Rx[4]-----◯-------------------\n", + " | \n", + "q1: ---Ry[1]---Rx[5]-----CZ------◯-----------\n", + " | \n", + "q2: ---Ry[2]---Rx[6]-------------CZ------◯---\n", + " | \n", + "q3: ---Ry[3]---Rx[7]---------------------CZ--\n" + ] + } + ], + "source": [ + "print_circuit(gates, qubit_inds, param_inds);" + ] + }, + { + "cell_type": "markdown", + "id": "a7d28ba0", + "metadata": {}, + "source": [ + "## Define cost function" + ] + }, + { + "cell_type": "markdown", + "id": "f8ae9870", + "metadata": {}, + "source": [ + "We will be working with the cost function\n", + "$$\n", + " c(\\ket{\\psi(\\theta)}) = \\langle \\psi (\\theta) | X_1 | \\psi (\\theta) \\rangle.\n", + "$$\n", + "Expectation values are known to present with a barren plateau for sufficiently expressible circuits (such as the hardware efficient ansatz we employ) and sufficiently deep circuits." + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "id": "42bfbc11", + "metadata": {}, + "outputs": [], + "source": [ + "observables = [[\"X\"]]\n", + "qubits_to_measure = [[1]]\n", + "coefficients = [1.0]\n", + "\n", + "# Get function that computes expectation value from quantum state\n", + "statetensor_to_expectation_func = qujax.get_statetensor_to_expectation_func(\n", + " observables, qubits_to_measure, coefficients\n", + ")" + ] + }, + { + "cell_type": "markdown", + "id": "06c44f5c", + "metadata": {}, + "source": [ + "## Measure gradients" + ] + }, + { + "cell_type": "markdown", + "id": "0a7f77cd", + "metadata": {}, + "source": [ + "Define parameters" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "id": "063a6fd8", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Total number of gradient samples: 5000\n" + ] + } + ], + "source": [ + "n_qubits = [2, 4, 6, 8, 10]\n", + "circuit_depth = 100\n", + "rng_seed = 0\n", + "batch_size = 100\n", + "n_batches = 50\n", + "\n", + "print(\"Total number of gradient samples:\", n_batches * batch_size)" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "id": "c1752b93", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Computing variance for 2 qubits\n", + "Variance for 2 qubits: 1.3082\n", + "Computing variance for 4 qubits\n", + "Variance for 4 qubits: 0.3036\n", + "Computing variance for 6 qubits\n", + "Variance for 6 qubits: 0.0788\n", + "Computing variance for 8 qubits\n", + "Variance for 8 qubits: 0.0195\n", + "Computing variance for 10 qubits\n", + "Variance for 10 qubits: 0.0047\n" + ] + } + ], + "source": [ + "results_list = []\n", + "rng_key = PRNGKey(rng_seed)\n", + "\n", + "for q in n_qubits:\n", + " print(f\"Computing variance for {q} qubits\")\n", + " # qujax specification of hardware efficient ansatz\n", + " gates, qubit_inds, param_inds, nr_of_parameters = hardware_efficient_ansatz(q)\n", + "\n", + " # Get function that returns one application of the circuit\n", + " params_to_statetensor = qujax.get_params_to_statetensor_func(\n", + " gates, qubit_inds, param_inds\n", + " )\n", + "\n", + " # Allow for an arbitrary number of circuit repetitions while avoiding compilation overhead\n", + " repeated_circuit = repeat_circuit(params_to_statetensor, nr_of_parameters)\n", + "\n", + " def expectation_func(angles: jax.Array, statetensor_in: jax.Array):\n", + " return statetensor_to_expectation_func(repeated_circuit(angles, statetensor_in))\n", + "\n", + " gradient_func = jax.grad(expectation_func)\n", + "\n", + " def gradient_of_first_angle(\n", + " angles: jax.Array, statetensor_in: jax.Array\n", + " ) -> jax.Array:\n", + " return gradient_func(angles, statetensor_in)[0]\n", + "\n", + " # Batched gradient sampling\n", + " vectorized_gradient_of_first_angle = jax.jit(\n", + " jax.vmap(gradient_of_first_angle, (0, None))\n", + " )\n", + "\n", + " initial_state = all_zeros_statetensor(q)\n", + "\n", + " sample_list = []\n", + " for i in range(n_batches):\n", + " rng_key, parameters_rng = jax.random.split(rng_key)\n", + "\n", + " random_angles = jax.random.uniform(\n", + " parameters_rng, (batch_size, circuit_depth * nr_of_parameters)\n", + " )\n", + " samples = vectorized_gradient_of_first_angle(random_angles, initial_state)\n", + "\n", + " sample_list.append(samples)\n", + "\n", + " variance = jnp.var(jnp.stack(sample_list))\n", + "\n", + " print(f\"Variance for {q} qubits: {variance:.4f}\")\n", + "\n", + " results_list.append(variance)" + ] + }, + { + "cell_type": "markdown", + "id": "660d86d3", + "metadata": {}, + "source": [ + "# Plot Results" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "id": "b7980e0b", + "metadata": {}, + "outputs": [], + "source": [ + "import matplotlib.pyplot as plt" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "id": "86ed5a5c", + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAk4AAAGwCAYAAABfKeoBAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABLjUlEQVR4nO3deVhU9f4H8PeZAYYdBJRFQFERxQUUUFk0zaVw19wzt0wtFc3bXr+rbVp200RwL5cyt1RcS8lEAfcFxQ13QUFR2fdlzu8PlESwBhg4M8P79Tw8z+UwzHlP3Zy35/s53xFEURRBRERERP9KJnUAIiIiIm3B4kRERESkIhYnIiIiIhWxOBERERGpiMWJiIiISEUsTkREREQqYnEiIiIiUpGe1AF0iVKpRGJiIszMzCAIgtRxiIiISAWiKCIzMxMODg6Qyf75mhKLkxolJibCyclJ6hhERERUBQkJCXB0dPzHx7A4qZGZmRmAkn/w5ubmEqchIiIiVWRkZMDJyan0ffyfsDip0dPlOXNzcxYnIiIiLaPKmA2Hw4mIiIhUxOJEREREpCIWJyIiIiIVsTgRERERqYjFiYiIiEhFLE5EREREKmJxIiIiIlIRixMRERGRiliciIiIiFTE4qTBMvIKkZSeW+HPktJzkZFXWMuJiIiI6jYWJw2VkVeIsT+dwPDlx5CYVrY8JablYvjyYxj70wmWJyIiolrE4qShsvOL8DirAPEpORix4u/ylJiWixErjiE+JQePswqQnV8kcVIiIqK6g8VJQ9lbGGHjpE5wtjJGfEoOBi85gv0Xk0pLk7OVMTZO6gR7CyOpoxIREdUZLE4azMGypDw1tDTC/Yw8TPr5TJnS5GDJ0kRERFSbWJw0nIOlET7r27LMsYBm1rAzN5QoERERUd3F4qThEtNyMW/vlTLHfj2RgDE/nUB6DgfDiYiIahOLkwZ7dhDc2coYW9/2hZWJAQAg6vojBAYfxsXEdIlTEhER1R0sThoqKT233CC4VyMr7J4eULpMl5iWh0FLjmDzqQSJ0xIREdUNLE4aykShB2tTg3KD4A6WRtj2jh8aWhrB3FAPBUVKfPDbeXy09TzyCoslTk1ERKTbBFEURalD6IqMjAxYWFggPT0d5ubm1X++vEJk5xdVuOVAUnoujPTl+PnoHSz48ypEEWjd0BxLX/eCk5Vxtc9NRERUV1Tm/ZtXnDSYuaH+C/dpsrcwgqWxAaZ3d8W6CR1Qz1gfF+5loO/iKBy8klzLSYmIiOoGFicd0Nm1PnYHdYaHkyXScwsxfs1JLNgfh2IlLyYSERGpE4uTjmhoaYTNkzvhjU6NAADBf13HuNUnkJJdIHEyIiIi3cHipEMUenJ8ObA1Fg73gKG+DJHXHqFvcCTOxqdKHY2IiEgnsDjpoEHtHBE21R8uNiZITM/DsOVH8fPR2+B9AERERNXD4vSc3bt3w83NDa6urli1apXUcaqshZ05dk7zx6ut7FBYLOL/dlzEu5tikFNQJHU0IiIircXtCJ5RVFQEd3d3HDx4EObm5mjfvj2OHz8OKysrlX5f3dsRqIMoilgVeQvf/HEFxUoRzW1NsWy0F5rUN5U6GhERkUbgdgRVdOLECbRq1QoNGzaEmZkZevfujX379kkdq1oEQcBbXZrg14kdUd9MgasPstA/JBq/xyZJHY2IiEjr6FRxOnz4MPr16wcHBwcIgoCwsLByj1myZAlcXFxgaGgILy8vREZGlv4sMTERDRs2LP3e0dER9+7dq43oNa5jE2vsmR6ADo2tkJVfhLfXn8HXey6hqFgpdTQiIiKtoVPFKTs7Gx4eHggJCanw55s2bcLMmTPx6aef4uzZs+jcuTMCAwMRHx8PABUOTwuC8MLz5efnIyMjo8yXJmtgboj1b3XEpC5NAAArI29h1KrjSM7IkzgZERGRdtCp4hQYGIivvvoKgwcPrvDnCxYswJtvvomJEyeiZcuW+OGHH+Dk5ISlS5cCABo2bFjmCtPdu3dhb2//wvPNmzcPFhYWpV9OTk7qfUE1QF8uwye9W2Lp6+1hqtDDiVsp6LM4CsdvPpY6GhERkcbTqeL0TwoKCnD69Gn06tWrzPFevXrhyJEjAIAOHTrgwoULuHfvHjIzM7F371688sorL3zOjz/+GOnp6aVfCQkJNfoa1CmwjT12TPNHc1tTPMzMx6hVx7Hy8E1uWUBERPQP6kxxevToEYqLi2Fra1vmuK2tLe7fvw8A0NPTw/fff49u3bqhXbt2eP/992Ftbf3C51QoFDA3Ny/zpU2a1jdF2FR/DPB0QLFSxNd7L+Od9WeQmVcodTQiIiKNpCd1gNr2/MySKIpljvXv3x/9+/ev7ViSMTbQww/DPeHdqB6+2H0Jv1+4j7j7mVg62gtudmZSxyMiItIodeaKk42NDeRyeenVpaeSk5PLXYWqawRBwBu+jbF5si/sLQxx81E2BoZGI+ysbtxRSEREpC51pjgZGBjAy8sL4eHhZY6Hh4fDz89PolSapZ1zPeyeHoDOrjbILSzGzE0x+O+OC8gvKpY6GhERkUbQqeKUlZWFmJgYxMTEAABu3bqFmJiY0u0GZs2ahVWrVuGnn37C5cuX8e677yI+Ph5Tpkyp1nlDQ0Ph7u4OHx+f6r4EyVmbKrBmfAdMf7kZAGDd0TsYvvwYEtNyJU5GREQkPZ36yJWIiAh069at3PGxY8dizZo1AEo2wJw/fz6SkpLQunVrLFy4EF26dFHL+TXxI1eq468rDzBzYwwy8opQz1gfwSPbobNrfaljERERqVVl3r91qjhJTdeKEwAkpORgyi+ncTExA4IA/Kdnc7zTtRlkshdvDEpERKRN+Fl1pDZOVsbY+rYfRvg4QRSB/+2/ionrTiE9h1sWEBFR3cPiRP/KUF+Ob15ri/mvtYVCT4a/riSjz+JIXLiXLnU0IiKiWsXiRCob5uOErW/7wcnKCHdTczF46RFsOhkvdSwiIqJaw+JEldK6oQV2T+uM7i0aoKBIiQ+3xuKD384hr5BbFhARke5jcVIDXdqOQBUWxvpYOcYb77/iBpkAbD51F4OXHEH84xypoxEREdUo3lWnRrp4V92/ibr2CEEbzyIluwDmhnpYMMwTPdzr9k7sRESkXXhXHdWaAFcb7AkKQDtnS2TkFWHiulP4bt8VFCvZx4mISPewOFG12VsYYdMkX4zzawwACD14A2N+Oo7HWfnSBiMiIlIzFidSCwM9Geb0b4VFIzxhpC9H9PXH6Ls4CmfiU6WORkREpDYsTqRWAzwbYsc0fzSpb4Kk9DwMX34Ua4/cBkfpiIhIF7A4kdo1tzXDzmkB6N3GDoXFImbvvIgZG2OQnV8kdTQiIqJqYXGiGmGq0EPoqPb4rE9LyGUCdp5LxMDQaFxPzpI6GhERUZWxOKlBXdvHSVWCIGBi5ybY8FYnNDBT4FpyFgaERGFvbJLU0YiIiKqE+zipUV3cx0lVyZl5mP7rWRy/lQIAmBjggg8DW0Bfzu5ORETS4j5OpHEamBli/cSOmPxSEwDAqqhbGLXyGB5k5EmcjIiISHUsTlRr9OQyfBzYEstGe8FMoYeTt1PRJzgKx24+ljoaERGRSlicqNa92toOO6cHwM3WDI+y8vH6quNYfugGtywgIiKNx+JEknCxMcH2qX4Y1K4hipUi5v1+BVN+OY2MvEKpoxEREb0QixNJxthADwuGeeCrga1hIJdh38UH6L84CpeTMqSORkREVCEWJ5KUIAgY3akRNk/xRUNLI9x+nINBS6Kx7cxdqaMRERGVw+JEGsHTyRK7pgegs6sN8gqVmLX5HD7dHov8omKpoxEREZVicVIDboCpHlYmBlgzvgOCursCANYfj8ewZUdxLy1X4mREREQluAGmGnEDTPU5eCUZMzfFID23EPWM9fHDiHZ4qXl9qWMREZEO4gaYpPW6tWiA3dMD0KahBVJzCjFu9Qks+vMalEr2fCIikg6LE2ksJytjbJnii5EdnCGKwMI/r2LC2pNIzS6QOhoREdVRLE6k0Qz15Zg3uA2+G9IWCj0ZIuIeou/iKJy/myZ1NCIiqoNYnEgrDPV2wrZ3/NDI2hj30nIxZOlR/Ho8nruNExFRrWJxIq3RysECO6cFoEdLWxQUK/HJ9li8/9t55BZwywIiIqodLE6kVSyM9LHiDS98+GoLyATgt9N3MXjpEdx+lC11NCIiqgNYnEjryGQC3u7aFL9M7AgbUwNcTspAv5Ao7L94X+poRESk41icSGv5NbXB7umd0d7ZEpl5RZj082l8+8cVFBUrpY5GREQ6isWJtJqdhSE2TvLFeP/GAIClETcw5qcTeJSVL20wIiLSSSxOasCPXJGWgZ4Ms/u1wuKR7WBsIMeRG4/RJzgSp++kSB2NiIh0DD9yRY34kSvSu56cick/n8aNh9nQkwn4tE9LjPNrDEEQpI5GREQaih+5QnVWswZm2DEtAH3a2qNIKeLzXZcwfcNZZOcXSR2NiIh0AIsT6RxThR5CRrbDf/u6Q08mYPf5JAwIjcb15EypoxERkZZjcSKdJAgCJgS4YOOkTrA1V+B6chb6h0Rj9/lEqaMREZEWY3Einebd2Aq7p3eGbxNr5BQUY9qvZ/H5rosoKOKWBUREVHksTqTz6psp8PObHfB216YAgNXRtzFy5THcT8+TOBkREWkbFieqE/TkMnz4aguseMMLZgo9nL6Tir6LI3Hk+iOpoxERkRZhcaI6pVcrO+yaHoAWdmZ4lFWA0T8ex9KIG+CuHEREpAoWJ6pzGtuYYPs7/hjcviGUIvDtH1cw6efTSM8tlDoaERFpOBYnqpOMDOT4fqgH5g5qAwO5DOGXHqB/SBQuJWZIHY2IiDQYixPVWYIgYFRHZ/z2ti8aWhrhzuMcDFoSjd9O35U6GhERaSgWJ6rz2jpaYvf0ALzUvD7yi5R4b8s5fLwtFnmFxVJHIyIiDcPiRASgnokBVo/zwbs9mkMQgA0n4jFs+VEkpORIHY2IiDQIi5MahIaGwt3dHT4+PlJHoWqQyQTM6OGKNeM7wNJYH+fvpqNfSBQi4pKljkZERBpCEHkfttpU5tOVSbPdTc3B1PVncO5uOgQBCHrZFUHdXSGXCVJHIyIiNavM+zevOBFVwLGeMTZP8cXrHZ0hisCiA9cwfs1JpGQXSB2NiIgkxOJE9AIKPTm+HtQG3w/1gKG+DIevPkS/xVE4l5AmdTQiIpIIixPRv3jNyxHb3/FHY2tj3EvLxdBlR7H++B3uNk5EVAexOBGpoKW9OXZOD0Avd1sUFCvx6fYL+M+Wc8gt4JYFRER1CYsTkYrMDfWx/A0vfBzYAjIB2HbmHgYticatR9lSRyMiolrC4kRUCYIgYPJLTbF+YifYmBrgyv1M9F8chX0X70sdjYiIagGLE1EV+Da1xp6gzvBuVA+Z+UWY/PNpzPv9MoqKlVJHIyKiGsTiRFRFtuaG2DCpE94McAEALD90E6N/PI7kzDyJkxERUU1hcSKqBn25DP/X1x2ho9rDxECOYzdT0Dc4Cidvp0gdjYiIagCLE5Ea9Glrjx3TAuDawBTJmfkYseIYVkXeRHpuAZLScyv8naT0XGTkFdZyUiIiqg5+5Ioa8SNXKDu/CB9ti8Wuc4kAAEtjfZgq9LB5si8cLI1KH5eYlosRK47B2tQAayd0gLmhvlSRiYjqPH7kCpFETBR6CB7hiTn93CGXCUjLKcTd1Fy8tvQIEtNKrjw9LU3xKTl4nFWA7PwiiVMTEZGqWJyI1EwQBIzzd8Hmyb6ob6YAACSl56Hf4iicvpNSWpqcrYyxcVIn2FsY/cszEhGRpuBSnRpxqY6e9ygrH5N/Po3Td1LLHH9amp5dviMiImlwqY5IQ9iYKrB5si8Gt29Y5viMHs1YmoiItBCLkxqEhobC3d0dPj4+UkchDfQgIw+nbpe94vSfzeexMPwqPyiYiEjLcKlOjbhUR897dhDc2coYc/q74531Z5BXWLLDeKcmVlg8sn3pLBQREdU+LtURaYCk9Nxyg+Avt7DFgVkvoZ5xyfYDx26moNfCQ/jrygOJ0xIRkSpYnIhqiIlCD9amBuUGwRvWM8aeoM6wMzeEob4MqTmFmLDmFP674wLyCoslTk1ERP+ES3VqxKU6el5GXiGy84sq3HIgKT0XejIBSyJuYHX0bQBAc1tTBI9shxZ2/P8PEVFtqcz7N4uTGrE4UVVFxCXjvS3n8SgrHwZ6Mnwc2ALj/BpDEASpoxER6TzOOBFpma5uDfDHzM7o5lYfBUVKfL7rEsavOYmHmflSRyMiomewOBFpCBtTBX4a54PP+7eCgZ4MEXEPEbjoMA5eSZY6GhERPcHiRKRBBEHAWL/G2DUtAC3szPAoqwDj15zEnJ0XOThORKQBWJyINJCbnRnCpvpjvH9jAMCaI7cxICQacfczpQ1GRFTHsTgRaShDfTlm92uF1eN9YGNqgLgHmegXEoU10be44zgRkURYnIg0XDe3Bvh9RpfSwfE5uy5hAgfHiYgkweJEpAXqm5UMjs/p5w4DPRkOcnCciEgSLE5EWkIQBIzzd8HOaf5ws+XgOBGRFFiciLRMCztz7Jjmj3F+jQGUDI4PDOXgOBFRbWBxItJChvpyzOnfCqvHlQyOX7lfMji+9shtDo4TEdUgFiciLdatRcngeNcng+Ozd17Em2tP4VEWB8eJiGoCixORlqtvpsDqZwbH/7qSjFd/iEREHAfHiYjUjcWJSAc8Ozje3NYUj7LyMW71SXy+i4PjRETqxOJEpENa2Jlj57QAjPVtBABYHV0yOH71AQfHiYjUgcWJSMcY6svx+YDW+GmcN6xNngyOL47CuqMcHCciqi4WJyId9XILW/w+szNeal4f+UVK/HfHRUzk4DgRUbWwOBHpsAZmhlg9zgf/7esOA7kMB54Mjh+6+lDqaEREWonFiUjHyWQCJgS4YMc0f7g2KBkcH/vTCXyx6xIHx4mIKonFSQ1CQ0Ph7u4OHx8fqaMQvVBLe3Psmv734PhP0bcwMDQa1zg4TkSkMkHktKjaZGRkwMLCAunp6TA3N5c6DtEL/XXlAd7fch6Pswug0JPhsz4tMbpTIwiCIHU0IqJaV5n3b15xIqqDnh8c/78dF/HWulN4zMFxIqJ/xOJEVEc9HRz/vyeD439eTsariyJxmIPjREQvxOJEVIfJZALeDHBB2NSSwfGHmfkY89MJfLn7EvKLODhORPQ8FicigrtDyeD4G51KBsd/jLqFgaFHODhORPQcFiciAlCy4/iXA1tj1RhvWJkY4HJSBvoujsLPx+5wx3EioidYnIiojB7utvhjRmd0drUpGRwPu4C31p3m4DgREViciKgCDcwNsXZ8h2cGxx/g1UWRiLzGwXEiqttYnIioQhUNjr/x4wl8xcFxIqrDWJyI6B+5O5hj57S/B8dXRd3CoNAjuJ7MwXEiqntYnIjoXxkZlAyOr3wyOH7pyeD4LxwcJ6I6hsWJiFTW85nB8bxCJT4Lu4BJP59GSnaB1NGIiGoFixMRVcrTwfHP+rSEgVyG8EsP8OoPhxF17ZHU0YiIahyLExFVmkwmYGLnJtg+1Q9N65sgOTMfo388jq/3cHCciHQbixMRVVkrBwvsnt4Zr3d0BgCsjHw6OJ4lcTIioprB4kRE1WJkIMfXg9pgxRteqGes/2RwPBLrj3NwnIh0D4sTEalFr1Z22DezS+ng+KfbOThORLqHxYmI1ObZwXF9ucDBcSLSOSxORKRWpYPj7/iXGRyfu/cyCoqUUscjIqoWFiciqhGtG5YMjo96Mji+4vBNDFoSzcFxItJqLE5EVGOMDOSYO6gNlj8ZHL+YWDI4/uvxeA6OE5FWYnEiohr3Sis7/DGzC/ybWSOvUIlPtsdiyi+nkcrBcSLSMixORFQrbM0N8fOEjvikdwvoywXsu/gAry46jOjrHBwnIu3B4kREtUYmEzCpS1Nsf8cfTeqb4EFGyeD4PA6OE5GWYHEiolpXMjgegFEdnSGKwPLDNzF4aTRuPOTgOBFpNhYnIpKEsYFe6eC4pbE+LtzLQN/gKGw4wcFxItJcLE5EJKlXnuw47t/MGrmFxfh4GwfHiUhzsTgRkeSeDo5/HPj34Hjgokgc4eA4EWkYQVTjNfG7d+/i7t27yMrKgkwmg7m5OVxcXGBtba2uU2i0jIwMWFhYID09Hebm5lLHIdJKF+6lI2jDWdx8lA1BACZ1aYL/9HSDgR7/nkdENaMy79/VKk4FBQXYtGkTtmzZgujoaMhkMtSrVw/16tVDUVERUlNTkZqaCjMzM3Tt2hXDhg1D7969IZPp5h+ALE5E6pFTUIQvd1/ChhMJAIA2DS2waIQnmtQ3lTgZEemiGi9Ooihi9erV+PXXX9GtWzd06tQJPj4+LzxZamoqoqOjER0djZiYGIwYMQJjx46t7Gk1HosTkXr9ceE+Ptp2Hmk5hTDSl2N2P3cM93GCIAhSRyMiHVKjxenhw4f46quvMGjQIHTt2rVKAcPDw3Hw4EG8//77qFevXpWeQxOxOBGp3/30PMzaHIMjNx4DAAJb22He4DawNDaQOBkR6YoaK045OTnYu3cvBg8eXO3lNqVSie3bt+O1116r1vNoEhYnopqhVIpYGXkT/9sfh8JiEXbmhlgw3AN+TW2kjkZEOqDWZpxeRC6Xo7i4WN1Pq/FYnIhqVuzddMzY+Pfg+OQuTTGrZ3MOjhNRtVTm/btG/rTh5nVEVBPaOFpgd1AARnZwgigCyw7dwGtLj+AmdxwnolpSI8VJ2wc3Bw0ahHr16mHIkCFSRyGi5xgb6GHe4LZYNro9LIz0EXsvHX2Co7DpJHccJ6KaV+vXt8PDwzFmzBiMHTsW48aNw/79+2s7wr8KCgrCunXrpI5BRP/g1db2+GNmZ/g2Kdlx/MOtsXhn/Rmk5XDHcSKqOVUqTklJSVAqq/ZJ5hs3bsS6deuwdu1arFmzBlu3bq3S89Skbt26wczMTOoYRPQv7C2M8MvEjvjw1RbQkwn4/cJ9BC6KxNEnd+AREalblYpTr169kJGRUfp9amoqTp06pdLvKpVKHDhwAFevXsWBAwdQWFhYqXMfPnwY/fr1g4ODAwRBQFhYWLnHLFmyBC4uLjA0NISXlxciIyMrdQ4i0h5ymYC3uzbFtnf84GJjgqT0PIxadQzf/nEFhcVV+wseEdGLVKk46enpwdLSsvR7CwsLTJ48WaXfDQkJQVJSErZt24akpCQEBwdX6tzZ2dnw8PBASEhIhT/ftGkTZs6ciU8//RRnz55F586dERgYiPj4+NLHeHl5oXXr1uW+EhMTK5WFiDRHW0dL7J4egOHeJYPjSyNKBsdvPcqWOhoR6RC9qvySo6MjoqOj4e/vDwCQyWQoKCg/V5CRkYHNmzfjxo0bsLKygqenJ1566SWMHj26yoEDAwMRGBj4wp8vWLAAb775JiZOnAgA+OGHH7Bv3z4sXboU8+bNAwCcPn26yud/Vn5+PvLz80u/f/YqHBHVPhOFHr4d0hZd3erjo22xOH83HX2CIzGnXysM9XbU+htXiEh6VbriFBISgilTpmDSpElYuXIlgoKC4OzsXO5xgYGBuHnzJpo1a4YFCxbgu+++Q+PGjfHdd9/VyN0vBQUFOH36NHr16lXmeK9evXDkyBG1n2/evHmwsLAo/XJyclL7OYio8gLb/D04nlNQjA+2nsfUX88gPadyowFERM+rUnFq1KgRzpw5gx49eiA+Ph7NmzfHpk2byj0uLS0Nc+fOxZtvvgk7Ozvs378fly5dQmJiIt59991qh3/eo0ePUFxcDFtb2zLHbW1tcf/+fZWf55VXXsHQoUOxd+9eODo64uTJkxU+7uOPP0Z6enrpV0JCQrXyE5H6PD84vjf2Pl5ddBjHbnJwnIiqrkpLdRs2bEBUVBSMjY3Rvn179OnTB6am5T+13N/fHxs2bMDIkSNLL5FbWlpi4cKFcHNzww8//FCt8C/y/OV4URQrdYl+3759Kj1OoVBAoVBUKhsR1Z6ng+N+Ta0xc1MMbj3KxsiVx/D2S03xbs/m0Jdzx3EiqpxK/6nxxRdf4Ndff0Xv3r2RnZ2NHTt24NVXX0V4eHi5x4aGhiI2NhZdunTB/fv3sX79emzbtg0fffRRhUWrumxsbCCXy8tdXUpOTi53FYqI6g4Pp5LB8WHejhBFYEnEDQxZegS3OThORJVU6eJ06NAhrF27Fn369MGlS5ewceNG/P7775g9e3a5x+rr62Pu3LnYvXs3Fi9ejLi4OERERMDOzg5//PGHWl7AswwMDODl5VWuxIWHh8PPz0/t5yMi7WGi0MP8IR5Y8np7mBvq4dzddPQOjsTmUwnccZyIVFbppboJEyaUbn5pYmKCb7/9Fs2bN//H39mzZ0/p0p6vry/69OlT5Q/BzcrKwvXr10u/v3XrFmJiYmBlZQVnZ2fMmjULb7zxBry9veHr64sVK1YgPj4eU6ZMqdL5iEi39G5jD08nS7y7KQbHb6Xgg9/O41DcQ8wd1AYWxvpSxyMiDSeI1firVnZ2NtatW4ekpCS88cYbcHV1BQDI5XIUFxcDKFnaO3nyJKZMmYI9e/YgJSUF8fHx+Pzzz9GzZ89KnzMiIgLdunUrd3zs2LFYs2YNgJINMOfPn4+kpCS0bt0aCxcuRJcuXar6Mv9VaGgoQkNDUVxcjKtXr6r06cpEJK1ipYhlh25gYfhVFClFOFgYYuFwT3RsYi11NCKqZRkZGbCwsFDp/btaxelFni1O3bt3x5YtW2BlZYWuXbsiIiIC6enpCAwMrJEtAqRUmX/wRKQZziWkYcbGs7j9OAeCALzTtSlm9uDgOFFdUpn37xr/k6Gipb2//vqrpk9LRKQSDydL7AnqXDo4HnrwBoYsO8rBcSKqUI1ccZLJZBV+CPCLlvZ0Ba84EWm3PeeT8PG288jIK4KJgRwfBbZA95YN4GBpXO6xSem5MFHowdyQc1FE2q7GluoyMzNx8+ZNeHh4VDskABw4cADdu3dXy3NpAhYnIu2XmJZbOjgOAMYGcmx7xw8t7MzLPGbEimOwNjXA2gkdWJ6ItFyNLdWZmZnh4cOHWLhwIbKysqoc8Pz58/j444/Rrl27Kj8HEVFNcLA0wq9vdcLkl5oAAHIKitEnOAp7zpd8CPjT0hSfkoPHWQXIzi+SMi4R1bIqLdXdvn0bX3zxBYyNjREYGAh/f39YWlr+4+/cunULERER2LFjB15++WW8/fbb0NfXrb+l8YoTkW7589IDTPnlNIqUJX9Mvta+IU7cSkFCai6crYyxcVInOFgaSZySiKqr1u6qu3jxIlatWoUdO3ZAoVDAwcEB9erVg6WlJYqLi5GWlobU1FRcuXIFjo6OGDlyJCZNmgQzM7OqnlIjcTsCIt117UEmBoZGI7uguPSYg4Uhfnvbj6WJSEdIsh3B3bt3cenSJSQkJCA7OxuCIMDCwgKNGzdG27Zt//WKlC7gFSci3XT6TgpeW3q09HtThRwLhnmiVys7CVMRkbpIvo9TXcXiRKR7np1pet54/8b4KLAFFHpyCZIRkbpo1D5ORETa6tnS5GxljK1v+8Kp3t/Lc6ujb+M1flgwUZ3C4kREVIGk9LKlaeOkTvBqZIVNk33hbFWyr5NMAC7cy0DfxVHYeS5R4sREVBtYnIiIKmCi0IO1qUG5u+ccLI2wcVInOFsZo6W9Odo7WyIrvwhBG87i422xyCss/pdnJiJtxhknNeKME5FuycgrRHZ+Eewtyt8993TncGN9ORYduIaQg9chioCbrRlCRrWDq61u3T1MpMtqdcapsLAQ3bp1w9WrV6v7VEREGsXcUL/C0gQA9hZGMDfUh55chv/0csPPEzrCxlSBuAeZ6B8SjS2nEsC/lxLpnmoXJ319fVy4cAGCIKgjj1YKDQ2Fu7s7fHx8pI5CRBIJcLXB7zM6I6CZDXILi/H+b+cxa/M5ZHFncSKdopaluv/85z/Q19fHN998o45MWotLdUSkVIpYeugGFoRfRbFSRBMbEywe1Q6tHCykjkZEL1CZ9289dZywoKAAq1atQnh4OLy9vWFiYlLm5wsWLFDHaYiINJ5MJmBqt2bo4GKFoA1ncfNRNgYtOYL/6+uO0R2d6/TVeSJdoJYrTt26dXvxCQQBf/31V3VPoRV4xYmInpWaXYD3tpzDgSvJAIDA1nb45rW2sDDSrc/pJNJ23DlcIixORPQ8URTxY9QtfPvHFRQWi3CsZ4SQUe3h6WQpdTQieoI7hxMRaQhBEDCxcxP8NsUPTlZGuJuaiyFLj2Dl4ZtQKvn3ViJto9YrTpcuXUJ8fDwKCgrKHO/fv7+6TqHReMWJiP5JRl4hPt4aiz2xSQCAl1s0wP+GesDKxEDiZER1W60v1d28eRODBg1CbGwsBEEo3bvk6RBkcXHd2EmXxYmI/o0oivj1RDw+33UJBUVK2JkbInhkO3RwsZI6GlGdVetLdTNmzICLiwsePHgAY2NjXLx4EYcPH4a3tzciIiLUcQoiIp0gCAJe79gIYe/4o0l9E9zPyMOIFUex+MA1FHPpjkjjqaU4HT16FF988QXq168PmUwGmUyGgIAAzJs3D0FBQeo4hUbjBphEVFnuDubYNS0Ag9s3hFIEvg+/ijE/HUdyZp7U0YjoH6ilOBUXF8PU1BQAYGNjg8TEkk8Jb9SoEeLi4tRxCo02depUXLp0CSdPnpQ6ChFpEROFHhYM88T/hnrASF+O6OuP0XtRJCKvPZQ6GhG9gFqKU+vWrXH+/HkAQMeOHTF//nxER0fjiy++QJMmTdRxCiIinTXEyxG7pvujhZ0ZHmUVYMxPJ/DdvisoKlZKHY2InqOW4vTZZ59BqSz5D/yrr77CnTt30LlzZ+zduxfBwcHqOAURkU5r1sAMYVP9MaqjM0QRCD14AyNXHkNiWq7U0YjoGdW6qy4mJgaenp4V/iwlJQX16tWrUx8vwLvqiEgddp9PxMdbY5GZXwRLY338b4gHerjbSh2LSGfV2l117du3h5eXF5YuXYr09PQyP7OysqpTpYmISF36tnXA7qAAtGlogbScQkxcdwpf7i7ZvoCIpFWt4hQdHY327dvjo48+gr29PUaPHo2DBw+qKxsRUZ3VyNoEv73tiwn+LgCAH6NuYciyI4h/nCNxMqK6TS0bYObm5mLz5s1YvXo1IiMj0bhxY0yYMAFjx46Fo6OjOnJqBS7VEVFNCL/0AO9tOYf03EKYKfTwzWtt0aetvdSxiHSGpB/ye+PGDaxevRrr1q1DUlISevbsib1796rzFBqLxYmIasq9tFwEbTiL03dSAQCvd3TG//V1h6G+XOJkRNpP0uIEAFlZWVi/fj0++eQTpKWl8SNXiIjUoLBYiYXhV7H00A2IItDCzgwho9qjWQNTqaMRabVa/8iVpw4dOoSxY8fCzs4OH3zwAQYPHozo6Gh1noKIqM7Sl8vwwastsHZ8B9iYGuDK/Uz0D4nC1tN3pY5GVGdUuzglJCTgyy+/RNOmTdGtWzfcuHEDixcvRmJiIlauXIlOnTqpI6dG40euEFFt6tK8PvYGdYZfU2vkFBTjP1vO4T+bzyE7v0jqaEQ6r1pLdT179sTBgwdRv359jBkzBhMmTICbm5s682kVLtURUW0qVopYcvA6Fv55FUoRaFrfBCGj2qOlPf/8IaqMyrx/61XnREZGRti6dSv69u0LuZwDikREtUkuEzC9uys6uFghaONZ3HiYjYGh0fhvP3eM6uDMvfSIakCNDIfXVbziRERSeZyVj/e2nMPBuJIPCO7T1h7zBreBuaG+xMmINJ9kw+FERCQNa1MFfhzrg097t4SeTMCe80noGxyF83fTpI5GpFNYnIiIdIRMJuCtLk2wZYovHOsZIT4lB68tPYIfo26BiwtE6sHiRESkY9o518OeoM54tZUdCotFfLn7Et5adxqp2QVSRyPSeixOREQ6yMJIH0tHt8eXA1rBQC7Dn5cfoHdwJE7dTpE6GpFWY3EiItJRgiDgDd/G2D7VDy42JkhKz8PwFccQevA6lEou3RFVBYsTEZGOa+VggV3TAzDQ0wHFShHf7YvD2NUn8DAzX+poRFqHxYmIqA4wVehh4XBPzB/SFob6MkRee4TewZGIvv5I6mhEWoXFiYiojhAEAcO8nbBrWgCa25riYWY+Rv94HAv2x6GoWCl1PCKtwOJERFTHuNqaYcfUAIzs4ARRBIL/uo5Rq47jfnqe1NGINB6LkxrwQ36JSNsYGcgxb3BbLBrhCRMDOU7cSkHgosM4eCVZ6mhEGo0fuaJG/MgVItJGtx9lY9qGM7hwLwMAMKlLE7zXyw0Gevy7NdUN/MgVIiJSWWMbE2x92w/j/BoDAFYcvolhy48iISVH2mBEGojFiYiIoNCTY07/Vlj+hhfMDfUQk5CG3sGR+D02SepoRBqFxYmIiEq90soOe2d0RjtnS2TmFeHt9Wfw3x0XkFdYLHU0Io3A4kRERGU41jPG5sm+mPJSUwDAuqN3MHjJEdx8mCVxMiLpsTgREVE5+nIZPgpsgTXjfWBtYoBLSRnouzgKYWfvSR2NSFIsTkRE9EJd3Rpg74zO6NTECjkFxZi5KQYf/HYOOQVFUkcjkgSLExER/SNbc0Osn9gJM3u4QiYAm0/dxYCQaMTdz5Q6GlGtY3EiIqJ/JZcJmNmjOdZP7IQGZgpcS85C/5AobDwRD24HSHUJixMREanMt6k19s7ojJea10d+kRIfbYtF0MYYZOYVSh2NqFawOBERUaXYmCqwepwPPgpsAblMwK5ziei7OAoX7qVLHY2oxrE4ERFRpclkAqa81BSbJ/uioaUR7jzOweAlR7Am+haX7kinsTgREVGVeTWqhz1BAejlbouCYiXm7LqEyT+fRnoOl+5IN7E4ERFRtVgaG2D5G16Y088dBnIZ9l96gN7BkTh9J1XqaERqx+JERETVJggCxvm7YOvbfmhkbYx7abkYtvwolh26AaWSS3ekO1iciIhIbdo4WmD39AD083BAsVLEN79fwbg1J/EoK1/qaERqweKkBqGhoXB3d4ePj4/UUYiIJGdmqI/gEZ74ZnAbKPRkOHz1IXovisTRG4+ljkZUbYLI2x/UJiMjAxYWFkhPT4e5ubnUcYiIJBd3PxPTfj2Da8lZkAlAUHdXTH/ZFXKZIHU0olKVef/mFSciIqoxbnZm2DHNH8O8HaEUgR/+vIbXVx3Dg4w8qaMRVQmLExER1ShjAz3MH+KBH4Z7wsRAjmM3U9B7USQi4pKljkZUaSxORERUKwa2a4hd0wPgbm+Ox9kFGLf6JL75/QoKi5VSRyNSGYsTERHVmib1TbHtHT+M8W0EAFh26AaGLz+Ku6k5EicjUg2LExER1SpDfTm+GNAaS19vDzNDPZyJT0PvRZHYd/G+1NGI/hWLExERSSKwjT32BnWGh5MlMvKKMPnn05iz8yLyi4qljkb0QixOREQkGScrY2yZ7ItJXZoAANYcuY3Xlh7B7UfZEicjqhiLExERScpAT4ZPerfE6nE+qGesjwv3MtB3cRR2nkuUOhpROSxORESkEbq1aIC9Mzqjg4sVsvKLELThLD7aeh65BVy6I83B4kRERBrD3sIIv07siKDurhAEYOPJBAwIjcK1B5lSRyMCwOJEREQaRk8uw6yezfHLmx1R30yBqw+y0C8kCptPJoCfEkZSY3EiIiKN5N/MBnuDOqOzqw3yCpX4YOt5vLspBln5RVJHozqMxYmIiDRWfTMF1o7vgA9edYNcJiAsJhH9FkfhYmK61NGojmJxIiIijSaTCXinazNsmtQJDhaGuPUoG4OWHMG6o7e5dEe1jsWJiIi0gndjK+wJ6oweLRugoEiJ/+64iLd/OYP03EKpo1EdwuJERERao56JAVaO8cb/9XWHvlzAHxfvo09wJM7Gp0odjeoIFiciItIqgiDgzQAXbH3bD85Wxribmouhy45ixeEbUCq5dEc1i8WJiIi0UltHS+wOCkCftvYoUoqYu/cK3lx7EinZBVJHIx3G4kRERFrL3FAfISPb4etBraHQk+Fg3EMELjqM4zcfSx2NdBSLExERaTVBEPB6x0YIm+qPpvVN8CAjHyNXHsPiA9eQmlOApPTcCn8vKT0XGXkcLKfKEUTey6k2GRkZsLCwQHp6OszNzaWOQ0RU52TnF+G/Oy5i65m7AABThR4sjPSxZYovHCyNSh+XmJaLESuOwdrUAGsndIC5ob5UkUkDVOb9m1eciIhIZ5go9PD9MA98P9QDRvpyZOUX4V5aLgaERiMxreTK09PSFJ+Sg8dZBcjmTuRUCSxORESkc17zcsSu6QFoVt8UAPAwMx+v/nAYJ249Li1NzlbG2DipE+wtjP7l2Yj+xuKkBqGhoXB3d4ePj4/UUYiI6IlmDUyxOygAg9o1BABk5BVh2PKypenZ5TsiVXDGSY0440REpJkW/3UN3++/Wvr9J71bYFKXphImIk3CGSciIqInEtNyseXU3TLH5u69gs/CYlFYrJQoFWkrFiciItJZzw6ClyzPdYSpQg8A8MuxeLy29Ejp0DiRKliciIhIJyWlP1+aOqFTExvsf7cLbEwNAADn76YjcFEkDl5JljgtaQsWJyIi0kkmCj1YmxqUGwR3sDTCzmkBsLcwhJGBHOm5hRi/5iS++f0Kl+7oX3E4XI04HE5EpFky8gqRnV9U4ZYDSem50JfLsPjANaw9egcA4N2oHhaPasctCuoYDocTERGh5LPsXlSC7C2MYGOqwOcDWmPJ6+1hptDDqTup6L0oEgfjuHRHFWNxIiKiOq93G3vsDgpA64bmSM0pxPjVJ/HtH1dQxKU7eg6LExEREYBG1ibY+rYfxvo2AgAsjbiBkSuPvfBDgqluYnEiIiJ6QqEnx+cDWiN0VHuYKvRw8nYq+gRHIYJLd/QEixMREdFz+rS1x+7pAWjlYI6U7AKMW30S87l0R2BxIiIiqlBjm5KluzFPlu6WPFm6u5+eJ3EykhKLExER0QsY6svxxXNLd72DI7l0V4exOBEREf2LipbuvtvHpbu6iMWJiIhIBU+X7t7oVLJ0F3rwBkatOo4HGVy6q0tYnIiIiFRkqC/HlwNbI2RUO5gq9HDiVgp6L4rE4asPpY5GtYTFiYiIqJL6tnXA7ukBcLc3x+PsAoxdfQL/2xfHpbs6gMWJiIioChrbmGDbO34Y3ckZogiEHLzOpbs6gMWJiIioigz15fhqYBssHsmlu7qCxYmIiKia+nk4YNdzS3ff7+fSnS5icSIiIlIDl+eW7hb/dR2vc+lO57A4ERERqcnTpbvgke1gYiDH8SdLd5HXuHSnK1iciIiI1Ky/hwN2B3VGyydLd2N+OoEF++NQrBSljkbVxOJERERUA1xsTLD9HT+83rFk6S74r+t4fdUxJHPpTquxOBEREdUQQ305vh7099LdsZsp6B3MpTttxuJERERUw/o/ueuupb05HmVx6U6bsTgRERHVgib1TbH9HT+M4tKdVmNxIiIiqiWG+nLMHdQGi0Z4llm6i7r2SOpopCIWJyIiolo2wLMhdk0PQAs7MzzKKsAbPx3HgvCrXLrTAixOREREEmhS3xRhU/0xssOTpbsD1zB61XEu3Wk4FiciIiKJGOrLMW/w30t3R28+Ru/gKERf59KdpmJxIiIiktgAz4bYWbp0l4/RPx7HQi7daSQWJyIiIg3Q9Lmlu0VPl+4yuXSnSViciIiINMTTpbsfhnvC+OnS3SIu3WkSFiciIiINM7Dds3fdlSzd/fAnl+40AYvTcxISEtC1a1e4u7ujbdu22LJli9SRiIioDvp76c4Jogj88Oc1vPEjl+6kJoiiyPr6jKSkJDx48ACenp5ITk5G+/btERcXBxMTk3/93YyMDFhYWCA9PR3m5ua1kJaIiOqCsLP38Mn2WOQUFMPGVIHgEZ7wa2YjdSydUZn3b15xeo69vT08PT0BAA0aNICVlRVSUlKkDUVERHXawHYNsXNaANxsS5buXufSnWS0rjgdPnwY/fr1g4ODAwRBQFhYWLnHLFmyBC4uLjA0NISXlxciIyOrdK5Tp05BqVTCycmpmqmJiIiqp1mDkqW7ET5/L92N+ek4HmbmSx2tTtG64pSdnQ0PDw+EhIRU+PNNmzZh5syZ+PTTT3H27Fl07twZgYGBiI+PL32Ml5cXWrduXe4rMTGx9DGPHz/GmDFjsGLFihdmyc/PR0ZGRpkvIiKimmJkIMc3r7XFwuEeMDaQI/r6Y/QOjsSRG7zrrrZo9YyTIAjYvn07Bg4cWHqsY8eOaN++PZYuXVp6rGXLlhg4cCDmzZun0vPm5+ejZ8+eeOutt/DGG2+88HFz5szB559/Xu44Z5yIiKimXU/OwtT1ZxD3IBMyAZjRvTmmvdwMcpkgdTStU2dnnAoKCnD69Gn06tWrzPFevXrhyJEjKj2HKIoYN24cXn755X8sTQDw8ccfIz09vfQrISGhytmJiIgq4+nS3XBvJyhFYOGfV7l0Vwt0qjg9evQIxcXFsLW1LXPc1tYW9+/fV+k5oqOjsWnTJoSFhcHT0xOenp6IjY2t8LEKhQLm5uZlvoiIiGqLkYEc3w5piwXDPGCkz6W72qAndYCaIAhlL1OKolju2IsEBARAqVTWRCwiIqIaMbi9I9o6WmDq+rOIe5CJ0auOY2aP5pjajUt36qZTV5xsbGwgl8vLXV1KTk4udxWKiIhIlzRrYFZm6W5B+FWM/ekEl+7UTKeKk4GBAby8vBAeHl7meHh4OPz8/CRKRUREVDueLt19P7Rk6S7q+iP0Do7E0RuPpY6mM7SuOGVlZSEmJgYxMTEAgFu3biEmJqZ0u4FZs2Zh1apV+Omnn3D58mW8++67iI+Px5QpU2osU2hoKNzd3eHj41Nj5yAiIlLVa16O2DXdH81tTfEwMx+vrzqGxQeuccNMNdC67QgiIiLQrVu3csfHjh2LNWvWACjZAHP+/PlISkpC69atsXDhQnTp0qXGs/EjV4iISJPkFhRj9s4L2HzqLgCgs6sNFg73hI2pQuJkmqUy799aV5w0GYsTERFpoq2n7+KzsAvILSxGAzMFFo1oB9+m1lLH0hh1dh8nIiIiKu81L0fsnOYP1wamSH5m6U7JpbtKY3EiIiKqA1xtzbBjmj+GejlCKQLfh1/F2NUn8CiLd91VBosTERFRHWFsoIfvhnrgf0/uuou89gi9F0Xi2E3edacqFiciIqI6ZshzS3ejVh5DyF9culMFi5MacDsCIiLSNk+X7oY8Wbr7334u3amCd9WpEe+qIyIibbTlVAL+b8cF5BUqYWuuQPCIdujYpO7cdce76oiIiEhlQ72dsHNaAFwbmOJBRj5GrjyG0IPXuXRXARYnIiIiQvMnS3eD2zeEUgS+2xeHsatP4DGX7spgcSIiIiIAJXfdLRjmie+GtIWhvqzkrrvgSBznXXelWJyIiIiojKHeTtgxNQBN65tw6e45LE5ERERUjpudGXZOC8Dgdn8v3Y1bc7LOL92xOBEREVGFTBR6+H6YB+Y/Wbo7fPUhegdH4sStFKmjSYbFSQ24jxMREekqQRAwjEt3pbiPkxpxHyciItJl2flF+CzsArafvQcAeKl5fSwY5gFrU4XEyaqH+zgRERGR2pko9LBgmAfmv9YWCj0ZDl19iD7BUTh5u+4s3bE4ERERkcoEQcAwHyfsmOaPpvVNcD8jDyNWHMOSiLqxdMfiRERERJXWws4cO6cFYFC7hihWipj/RxzGrzmJlOwCqaPVKBYnIiIiqpKnS3ffvtamdOmu96JInV66Y3EiIiKiKhMEAcN9nLFjmj+aPLN0tzTihk4u3bE4ERERUbW1sDPHrmkBGOjpgGKliG//uIIJa3Vv6Y7FiYiIiNTCRKGHhcM9S5fuIuJKlu5O6dDSHYsTERERqc3Tpbuwqf5oYlOydDd8xTEsO6QbS3csTmrAncOJiIjKamlvjp3TAzDgydLdN79fwZs6sHTHncPViDuHExERlSWKIjadTMDsnReRX6SEvYUhFo9sB+/GVlJHK8Wdw4mIiEgjCIKAER3+XrpLStfupTsWJyIiIqpxT5fu+nv8vXQ3cd0ppGrZ0h2LExEREdUKU4UeFo3wxLzBbWCgJ8NfV5LROzgSp+9oz113LE5ERERUawRBwMgOzgh75++lu2HLj2G5lizdsTgRERFRrXN3KLt0N09Llu5YnIiIiEgST5fu5g76e+muT3AkTt9JlTraC7E4ERERkWQEQcCojiVLdy42JkhMz8Pw5Uex4rBmLt2xOBEREZHk3B3MsWt6APp5OKBIKWLu3it4a90pxKdkIyk9t8LfSUrPRUZeYa3mZHFSA+4cTkREVH2mCj0EP7N0d+BKMrp/fwgDQ6KRmFa2PCWm5WL48mMY+9OJWi1P3DlcjbhzOBERkXpcTEzHlJ9PIyG1pDBZGuljT1AAGtYzRmJaLkasOIb4lBw4Wxlj0+ROsLcwqvK5uHM4ERERabVWDhbYO6MzurdsAABIyy1EjwWHERH3oExp2jipeqWpsnjFSY14xYmIiEi9RFHE0ogbmL8vrszxp6XJwbL6pYlXnIiIiEgnCIKAd7o1w/+Gti1zfOFwD7WUpspicSIiIiKNlpiWi+AD18sce3fTuXID47WBxYmIiIg01vOD4Fvf9oWzlTHiU3IwYsWxWi9PLE5ERESkkZLSc8sNgns1ssLGSZ3KlKcX7fNUE1iciIiISCOZKPRgbWpQbhDcwdKotDxZmxrARKFXa5l4V50a8a46IiIi9crIK0R2flGFWw4kpefCRKEHc0P96p2jEu/ftVfRiIiIiCrJ3FD/hcWoNvdveopLdUREREQqYnEiIiIiUhGLExEREZGKWJzUIDQ0FO7u7vDx8ZE6ChEREdUg3lWnRryrjoiISPvws+qIiIiIagCLExEREZGKWJyIiIiIVMTiRERERKQi7hyuRk/n7DMyMiROQkRERKp6+r6tyv1yLE5qlJmZCQBwcnKSOAkRERFVVmZmJiwsLP7xMdyOQI2USiUSExNhZmYGQRDU+twZGRlwcnJCQkKCTm51oOuvD9D918jXp/10/TXq+usDdP811tTrE0URmZmZcHBwgEz2z1NMvOKkRjKZDI6OjjV6DnNzc538j+EpXX99gO6/Rr4+7afrr1HXXx+g+6+xJl7fv11peorD4UREREQqYnEiIiIiUhGLk5ZQKBSYPXs2FAqF1FFqhK6/PkD3XyNfn/bT9deo668P0P3XqAmvj8PhRERERCriFSciIiIiFbE4EREREamIxYmIiIhIRSxORERERCpicdJw8+bNg4+PD8zMzNCgQQMMHDgQcXFxUsdSm6VLl6Jt27alm5n5+vri999/lzpWjZk3bx4EQcDMmTOljqI2c+bMgSAIZb7s7OykjqVW9+7dw+jRo2FtbQ1jY2N4enri9OnTUsdSm8aNG5f7dygIAqZOnSp1NLUoKirCZ599BhcXFxgZGaFJkyb44osvoFQqpY6mNpmZmZg5cyYaNWoEIyMj+Pn54eTJk1LHqrLDhw+jX79+cHBwgCAICAsLK/NzURQxZ84cODg4wMjICF27dsXFixdrJRuLk4Y7dOgQpk6dimPHjiE8PBxFRUXo1asXsrOzpY6mFo6Ojvjmm29w6tQpnDp1Ci+//DIGDBhQa/8B1KaTJ09ixYoVaNu2rdRR1K5Vq1ZISkoq/YqNjZU6ktqkpqbC398f+vr6+P3333Hp0iV8//33sLS0lDqa2pw8ebLMv7/w8HAAwNChQyVOph7ffvstli1bhpCQEFy+fBnz58/Hd999h8WLF0sdTW0mTpyI8PBw/Pzzz4iNjUWvXr3Qo0cP3Lt3T+poVZKdnQ0PDw+EhIRU+PP58+djwYIFCAkJwcmTJ2FnZ4eePXuWfmZsjRJJqyQnJ4sAxEOHDkkdpcbUq1dPXLVqldQx1CozM1N0dXUVw8PDxZdeekmcMWOG1JHUZvbs2aKHh4fUMWrMhx9+KAYEBEgdo1bNmDFDbNq0qahUKqWOohZ9+vQRJ0yYUObY4MGDxdGjR0uUSL1ycnJEuVwu7t69u8xxDw8P8dNPP5UolfoAELdv3176vVKpFO3s7MRvvvmm9FheXp5oYWEhLlu2rMbz8IqTlklPTwcAWFlZSZxE/YqLi7Fx40ZkZ2fD19dX6jhqNXXqVPTp0wc9evSQOkqNuHbtGhwcHODi4oIRI0bg5s2bUkdSm507d8Lb2xtDhw5FgwYN0K5dO6xcuVLqWDWmoKAAv/zyCyZMmKD2DyuXSkBAAA4cOICrV68CAM6dO4eoqCj07t1b4mTqUVRUhOLiYhgaGpY5bmRkhKioKIlS1Zxbt27h/v376NWrV+kxhUKBl156CUeOHKnx8/NDfrWIKIqYNWsWAgIC0Lp1a6njqE1sbCx8fX2Rl5cHU1NTbN++He7u7lLHUpuNGzfizJkzWj1v8E86duyIdevWoXnz5njw4AG++uor+Pn54eLFi7C2tpY6XrXdvHkTS5cuxaxZs/DJJ5/gxIkTCAoKgkKhwJgxY6SOp3ZhYWFIS0vDuHHjpI6iNh9++CHS09PRokULyOVyFBcX4+uvv8bIkSOljqYWZmZm8PX1xZdffomWLVvC1tYWGzZswPHjx+Hq6ip1PLW7f/8+AMDW1rbMcVtbW9y5c6fGz8/ipEWmTZuG8+fP69zfINzc3BATE4O0tDRs3boVY8eOxaFDh3SiPCUkJGDGjBnYv39/ub8N6orAwMDS/92mTRv4+vqiadOmWLt2LWbNmiVhMvVQKpXw9vbG3LlzAQDt2rXDxYsXsXTpUp0sTj/++CMCAwPh4OAgdRS12bRpE3755Rf8+uuvaNWqFWJiYjBz5kw4ODhg7NixUsdTi59//hkTJkxAw4YNIZfL0b59e4waNQpnzpyROlqNef6KqCiKtXKVlMVJS0yfPh07d+7E4cOH4ejoKHUctTIwMECzZs0AAN7e3jh58iQWLVqE5cuXS5ys+k6fPo3k5GR4eXmVHisuLsbhw4cREhKC/Px8yOVyCROqn4mJCdq0aYNr165JHUUt7O3ty5X4li1bYuvWrRIlqjl37tzBn3/+iW3btkkdRa3ef/99fPTRRxgxYgSAkoJ/584dzJs3T2eKU9OmTXHo0CFkZ2cjIyMD9vb2GD58OFxcXKSOpnZP79q9f/8+7O3tS48nJyeXuwpVEzjjpOFEUcS0adOwbds2/PXXXzr5H8HzRFFEfn6+1DHUonv37oiNjUVMTEzpl7e3N15//XXExMToXGkCgPz8fFy+fLnMH2jazN/fv9wWIFevXkWjRo0kSlRzVq9ejQYNGqBPnz5SR1GrnJwcyGRl3+7kcrlObUfwlImJCezt7ZGamop9+/ZhwIABUkdSOxcXF9jZ2ZXe/QmUzOYdOnQIfn5+NX5+XnHScFOnTsWvv/6KHTt2wMzMrHRt18LCAkZGRhKnq75PPvkEgYGBcHJyQmZmJjZu3IiIiAj88ccfUkdTCzMzs3LzaCYmJrC2ttaZObX33nsP/fr1g7OzM5KTk/HVV18hIyNDZ/4m/+6778LPzw9z587FsGHDcOLECaxYsQIrVqyQOppaKZVKrF69GmPHjoWenm69NfTr1w9ff/01nJ2d0apVK5w9exYLFizAhAkTpI6mNvv27YMoinBzc8P169fx/vvvw83NDePHj5c6WpVkZWXh+vXrpd/funULMTExsLKygrOzM2bOnIm5c+fC1dUVrq6umDt3LoyNjTFq1KiaD1fj9+1RtQCo8Gv16tVSR1OLCRMmiI0aNRINDAzE+vXri927dxf3798vdawapWvbEQwfPly0t7cX9fX1RQcHB3Hw4MHixYsXpY6lVrt27RJbt24tKhQKsUWLFuKKFSukjqR2+/btEwGIcXFxUkdRu4yMDHHGjBmis7OzaGhoKDZp0kT89NNPxfz8fKmjqc2mTZvEJk2aiAYGBqKdnZ04depUMS0tTepYVXbw4MEK3/vGjh0rimLJlgSzZ88W7ezsRIVCIXbp0kWMjY2tlWyCKIpizdczIiIiIu3HGSciIiIiFbE4EREREamIxYmIiIhIRSxORERERCpicSIiIiJSEYsTERERkYpYnIiIiIhUxOJEREREpCIWJyKq88LCwtCsWTPI5XLMnDmzxs6zZs0aWFpa/uNj5syZA09PzxrLQETVw+JERHXe5MmTMWTIECQkJODLL7+UNMt7772HAwcOlH4/btw4DBw4ULpARFSGbn2SIxHVWYWFhdDX16/072VlZSE5ORmvvPIKHBwcaiBZ5ZiamsLU1FTqGET0ArziREQap2vXrggKCsIHH3wAKysr2NnZYc6cOWUeIwgCli1bhgEDBsDExARfffVVhc+VmpqKMWPGoF69ejA2NkZgYCCuXbsGAIiIiICZmRkA4OWXX4YgCIiIiKjwea5du4YuXbrA0NAQ7u7uCA8PhyAICAsLK30uQRCQlpZW+jsxMTEQBAG3b98u81xhYWFo3rw5DA0N0bNnTyQkJJT+7Nmlujlz5mDt2rXYsWMHBEEozVdQUIBp06bB3t4ehoaGaNy4MebNm6faP1wiqhYWJyLSSGvXroWJiQmOHz+O+fPn44svvkB4eHiZx8yePRsDBgxAbGwsJkyYUOHzjBs3DqdOncLOnTtx9OhRiKKI3r17o7CwEH5+foiLiwMAbN26FUlJSfDz8yv3HEqlEoMHD4ZcLsexY8ewbNkyfPjhh1V6XTk5Ofj666+xdu1aREdHIyMjAyNGjKjwse+99x6GDRuGV199FUlJSaX5goODsXPnTmzevBlxcXH45Zdf0Lhx4yrlIaLK4VIdEWmktm3bYvbs2QAAV1dXhISE4MCBA+jZs2fpY0aNGvXCwgSUXCXauXMnoqOjSwvR+vXr4eTkhLCwMAwdOhQNGjQAgNIrWxX5888/cfnyZdy+fRuOjo4AgLlz5yIwMLDSr6uwsBAhISHo2LEjgJKC2LJlS5w4cQIdOnQo81hTU1MYGRkhPz+/TLb4+Hi4uroiICAAgiCgUaNGlc5BRFXDK05EpJHatm1b5nt7e3skJyeXOebt7f2Pz3H58mXo6emVlhQAsLa2hpubGy5fvqxylsuXL8PZ2bm0NAGAr6+vyr//LD09vTK5W7RoAUtLy0rlGTduHGJiYuDm5oagoCDs37+/SlmIqPJYnIhIIz0/6C0IApRKZZljJiYm//gcoii+8LggCCpnqeh5nv99mUxW7rGFhYUVPl9F565Mnvbt2+PWrVv48ssvkZubi2HDhmHIkCEq/z4RVR2LExHpLHd3dxQVFeH48eOlxx4/foyrV6+iZcuWlXqe+Ph4JCYmlh47evRomcfUr18fAJCUlFR6LCYmptxzFRUV4dSpU6Xfx8XFIS0tDS1atKjw3AYGBiguLi533NzcHMOHD8fKlSuxadMmbN26FSkpKSq/JiKqGhYnItJZrq6uGDBgAN566y1ERUXh3LlzGD16NBo2bIgBAwao/Dw9evSAm5sbxowZg3PnziEyMhKffvppmcc0a9YMTk5OmDNnDq5evYo9e/bg+++/L/dc+vr6mD59Oo4fP44zZ85g/Pjx6NSpU7n5pqcaN26M8+fPIy4uDo8ePUJhYSEWLlyIjRs34sqVK7h69Sq2bNkCOzu7f91ck4iqj8WJiHTa6tWr4eXlhb59+8LX1xeiKGLv3r2V2vNJJpNh+/btyM/PR4cOHTBx4kR8/fXXZR6jr6+PDRs24MqVK/Dw8MC3335b4RYJxsbG+PDDDzFq1Cj4+vrCyMgIGzdufOG533rrLbi5ucHb2xv169dHdHQ0TE1N8e2338Lb2xs+Pj64ffs29u7dW7pcSEQ1RxBfNARARET/SBAEbN++nTt7E9Uh/OsJERERkYpYnIiIiIhUxA0wiYiqiJMORHUPrzgRERERqYjFiYiIiEhFLE5EREREKmJxIiIiIlIRixMRERGRiliciIiIiFTE4kRERESkIhYnIiIiIhX9PzWzojqbyngJAAAAAElFTkSuQmCC", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "plt.plot(n_qubits, results_list)\n", + "plt.scatter(n_qubits, results_list, marker=\"x\")\n", + "plt.yscale(\"log\")\n", + "plt.xlabel(\"nr of qubits\")\n", + "plt.ylabel(r\"Var $\\left(\\frac{\\partial c}{\\partial \\theta_0}\\right)$\");" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "quantum", + "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.11.4" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/examples/classification.ipynb b/examples/classification.ipynb index 29168f8..80c8a6d 100644 --- a/examples/classification.ipynb +++ b/examples/classification.ipynb @@ -8,7 +8,7 @@ "source": [ "from jax import numpy as jnp, random, vmap, value_and_grad, jit\n", "import qujax\n", - "import matplotlib.pyplot as plt\n" + "import matplotlib.pyplot as plt" ] }, { @@ -31,7 +31,7 @@ "\n", "def classification_function(x, y):\n", " r = jnp.sqrt(x**2 + y**2)\n", - " return jnp.where((r > inner_rad) * (r < outer_rad), 1, 0)\n" + " return jnp.where((r > inner_rad) * (r < outer_rad), 1, 0)" ] }, { @@ -41,7 +41,7 @@ "outputs": [], "source": [ "linsp = jnp.linspace(-1, 1, 1_000)\n", - "Z = vmap(lambda x: vmap(lambda y: classification_function(x, y))(linsp))(linsp)\n" + "Z = vmap(lambda x: vmap(lambda y: classification_function(x, y))(linsp))(linsp)" ] }, { @@ -61,10 +61,10 @@ } ], "source": [ - "plt.figure(figsize=(4,4))\n", - "plt.xlabel(f'$x_0$')\n", - "plt.ylabel(f'$x_1$')\n", - "plt.contourf(linsp, linsp, Z, cmap='Purples');" + "plt.figure(figsize=(4, 4))\n", + "plt.xlabel(f\"$x_0$\")\n", + "plt.ylabel(f\"$x_1$\")\n", + "plt.contourf(linsp, linsp, Z, cmap=\"Purples\");" ] }, { @@ -81,8 +81,10 @@ "outputs": [], "source": [ "n_data = 1_000\n", - "x = random.uniform(random.PRNGKey(0), shape=(n_data, 2), minval=-1, maxval=1) # bivariate (x_0,x_1) data\n", - "y = classification_function(x[:, 0], x[:, 1]) # ground truth labels\n" + "x = random.uniform(\n", + " random.PRNGKey(0), shape=(n_data, 2), minval=-1, maxval=1\n", + ") # bivariate (x_0,x_1) data\n", + "y = classification_function(x[:, 0], x[:, 1]) # ground truth labels" ] }, { @@ -102,10 +104,10 @@ } ], "source": [ - "plt.figure(figsize=(4,4))\n", - "plt.xlabel(f'$x_0$')\n", - "plt.ylabel(f'$x_1$')\n", - "plt.scatter(x[:,0], x[:,1], alpha=jnp.where(y, 1, 0.2), s=10);" + "plt.figure(figsize=(4, 4))\n", + "plt.xlabel(f\"$x_0$\")\n", + "plt.ylabel(f\"$x_1$\")\n", + "plt.scatter(x[:, 0], x[:, 1], alpha=jnp.where(y, 1, 0.2), s=10);" ] }, { @@ -149,7 +151,7 @@ " for qi in range(layer, layer + n_qubits - 1, 2):\n", " gate_seq_seq += [\"CZ\"]\n", " qubit_inds_seq += [[qi % n_qubits, (qi + 1) % n_qubits]]\n", - " param_inds_seq += [[]]\n" + " param_inds_seq += [[]]" ] }, { @@ -194,7 +196,7 @@ "source": [ "angles_to_st = qujax.get_params_to_statetensor_func(\n", " gate_seq_seq, qubit_inds_seq, param_inds_seq\n", - ")\n" + ")" ] }, { @@ -222,7 +224,7 @@ " )\n", "\n", " angles = biases + weights_times_data\n", - " return angles\n" + " return angles" ] }, { @@ -231,7 +233,9 @@ "metadata": {}, "outputs": [], "source": [ - "param_and_x_to_st = lambda param, x_single: angles_to_st(param_and_x_to_angles(param, x_single))" + "param_and_x_to_st = lambda param, x_single: angles_to_st(\n", + " param_and_x_to_angles(param, x_single)\n", + ")" ] }, { @@ -249,9 +253,13 @@ "source": [ "def param_and_x_to_probability(param, x_single):\n", " st = param_and_x_to_st(param, x_single)\n", - " all_probs = jnp.square(jnp.abs(st)) # probabilities of all states given one x datapoint: it is a tensor of shape [2] x n_qubits\n", - " first_qubit_probs = jnp.sum(all_probs, axis=range(1, n_qubits)) # sum over all probabilities related to first qubit\n", - " return first_qubit_probs[1] # probability first qubit is 1\n" + " all_probs = jnp.square(\n", + " jnp.abs(st)\n", + " ) # probabilities of all states given one x datapoint: it is a tensor of shape [2] x n_qubits\n", + " first_qubit_probs = jnp.sum(\n", + " all_probs, axis=range(1, n_qubits)\n", + " ) # sum over all probabilities related to first qubit\n", + " return first_qubit_probs[1] # probability first qubit is 1" ] }, { @@ -301,9 +309,11 @@ "outputs": [], "source": [ "def param_to_cost(param): # use the expected distance as loss function\n", - " donut_probs = vmap(param_and_x_to_probability, in_axes=(None, 0))(param, x) # mapping over x\n", + " donut_probs = vmap(param_and_x_to_probability, in_axes=(None, 0))(\n", + " param, x\n", + " ) # mapping over x\n", " costs = jnp.where(y, 1 - donut_probs, donut_probs)\n", - " return costs.mean()\n" + " return costs.mean()" ] }, { @@ -321,7 +331,7 @@ "metadata": {}, "outputs": [], "source": [ - "param_to_cost_and_grad = jit(value_and_grad(param_to_cost))\n" + "param_to_cost_and_grad = jit(value_and_grad(param_to_cost))" ] }, { @@ -346,7 +356,7 @@ " cost, grad = param_to_cost_and_grad(param)\n", " costs = costs.at[i].set(cost)\n", " param = param - stepsize * grad\n", - " print(i, \"Cost: \", cost, end=\"\\r\")\n" + " print(i, \"Cost: \", cost, end=\"\\r\")" ] }, { @@ -366,10 +376,10 @@ } ], "source": [ - "plt.figure(figsize=(6,4))\n", + "plt.figure(figsize=(6, 4))\n", "plt.plot(costs)\n", - "plt.xlabel('Iteration')\n", - "plt.ylabel('Cost');" + "plt.xlabel(\"Iteration\")\n", + "plt.ylabel(\"Cost\");" ] }, { @@ -390,7 +400,7 @@ " lambda a: vmap(lambda b: param_and_x_to_probability(param, jnp.array([a, b])))(\n", " linsp\n", " )\n", - ")(linsp)\n" + ")(linsp)" ] }, { @@ -415,8 +425,8 @@ "circle_linsp = jnp.linspace(0, 2 * jnp.pi, 100)\n", "plt.plot(inner_rad * jnp.cos(circle_linsp), inner_rad * jnp.sin(circle_linsp), c=\"red\")\n", "plt.plot(outer_rad * jnp.cos(circle_linsp), outer_rad * jnp.sin(circle_linsp), c=\"red\")\n", - "plt.xlabel(f'$x_0$')\n", - "plt.ylabel(f'$x_1$');" + "plt.xlabel(f\"$x_0$\")\n", + "plt.ylabel(f\"$x_1$\");" ] }, { diff --git a/examples/generative_modelling.ipynb b/examples/generative_modelling.ipynb index c699e7d..b0657ee 100644 --- a/examples/generative_modelling.ipynb +++ b/examples/generative_modelling.ipynb @@ -40,7 +40,7 @@ "circuit_depth = 3\n", "init_rad = 0.001 / jnp.pi\n", "beta = 1e3\n", - "get_stepsize = lambda step: (step + 10) ** (-1/3)\n", + "get_stepsize = lambda step: (step + 10) ** (-1 / 3)\n", "n_steps = 500\n", "random_key = random.PRNGKey(0)" ] @@ -70,8 +70,8 @@ "metadata": {}, "outputs": [], "source": [ - "with open('hidalgo_stamp.txt', \"r\") as f:\n", - " all_lines = [line[:-1].rstrip().split(' ') for line in f if line != '\\n']\n", + "with open(\"hidalgo_stamp.txt\", \"r\") as f:\n", + " all_lines = [line[:-1].rstrip().split(\" \") for line in f if line != \"\\n\"]\n", "data = jnp.hstack([[float(a) for a in line] for line in all_lines]) * 1000\n", "data = jnp.array(data, dtype=int)" ] @@ -97,8 +97,8 @@ ], "source": [ "plt.hist(data, bins=50, alpha=0.5)\n", - "plt.ylabel('Frequency')\n", - "plt.xlabel('Stamp thickness $\\\\mu$m');" + "plt.ylabel(\"Frequency\")\n", + "plt.xlabel(\"Stamp thickness $\\\\mu$m\");" ] }, { @@ -137,9 +137,11 @@ "outputs": [], "source": [ "def mmd(kernel, data1, weights1, data2, weights2):\n", - " return expected_kernel(kernel, data1, weights1, data1, weights1) \\\n", - " - 2 * expected_kernel(kernel, data1, weights1, data2, weights2) \\\n", - " + expected_kernel(kernel, data2, weights2, data2, weights2)" + " return (\n", + " expected_kernel(kernel, data1, weights1, data1, weights1)\n", + " - 2 * expected_kernel(kernel, data1, weights1, data2, weights2)\n", + " + expected_kernel(kernel, data2, weights2, data2, weights2)\n", + " )" ] }, { @@ -160,28 +162,26 @@ "outputs": [], "source": [ "def get_circuit(n_qubits, depth):\n", - " \n", " n_params = 2 * n_qubits * (depth + 1)\n", - " \n", - " \n", - " gates = ['H'] * n_qubits + ['Rx'] * n_qubits + ['Ry'] * n_qubits\n", + "\n", + " gates = [\"H\"] * n_qubits + [\"Rx\"] * n_qubits + [\"Ry\"] * n_qubits\n", " qubit_inds = [[i] for i in range(n_qubits)] * 3\n", " param_inds = [[]] * n_qubits + [[i] for i in range(n_qubits * 2)]\n", - " \n", + "\n", " k = 2 * n_qubits\n", "\n", " for _ in range(depth):\n", " for i in range(0, n_qubits - 1):\n", - " gates.append('CZ')\n", + " gates.append(\"CZ\")\n", " qubit_inds.append([i, i + 1])\n", " param_inds.append([])\n", " for i in range(n_qubits):\n", - " gates.append('Rx')\n", + " gates.append(\"Rx\")\n", " qubit_inds.append([i])\n", " param_inds.append([k])\n", " k += 1\n", " for i in range(n_qubits):\n", - " gates.append('Ry')\n", + " gates.append(\"Ry\")\n", " qubit_inds.append([i])\n", " param_inds.append([k])\n", " k += 1\n", @@ -260,7 +260,7 @@ "\n", "\n", "def gaussian_kernel(s1, s2):\n", - " return jnp.exp(- jnp.square(s1 - s2) / bandwidth_sq)" + " return jnp.exp(-jnp.square(s1 - s2) / bandwidth_sq)" ] }, { @@ -280,6 +280,7 @@ "source": [ "data_probs = jnp.ones(len(data)) / len(data)\n", "\n", + "\n", "def param_to_mmd(param):\n", " st = param_to_st(param)\n", " probs = jnp.square(jnp.abs(st.flatten()))\n", @@ -313,7 +314,9 @@ "metadata": {}, "outputs": [], "source": [ - "init_param = random.uniform(init_key, shape=(n_params,), minval=-init_rad, maxval=init_rad)" + "init_param = random.uniform(\n", + " init_key, shape=(n_params,), minval=-init_rad, maxval=init_rad\n", + ")" ] }, { @@ -349,14 +352,18 @@ "for step in range(1, n_steps):\n", " cost_val, cost_grad = param_to_mmd_and_grad(params[step - 1])\n", " cost_vals = cost_vals.at[step - 1].set(cost_val)\n", - " \n", + "\n", " stepsize = get_stepsize(step)\n", "\n", - " new_param = params[step - 1] - stepsize * cost_grad\\\n", - " + jnp.sqrt(2 * stepsize / beta) * random.normal(train_keys[step - 1], shape=(n_params,))\n", + " new_param = (\n", + " params[step - 1]\n", + " - stepsize * cost_grad\n", + " + jnp.sqrt(2 * stepsize / beta)\n", + " * random.normal(train_keys[step - 1], shape=(n_params,))\n", + " )\n", " params = params.at[step].set(new_param)\n", "\n", - " print('Iteration:', step, '\\tCost:', cost_val, end='\\r')" + " print(\"Iteration:\", step, \"\\tCost:\", cost_val, end=\"\\r\")" ] }, { @@ -390,8 +397,8 @@ ], "source": [ "plt.plot(cost_vals)\n", - "plt.xlabel('Iteration')\n", - "plt.ylabel('MMD')\n", + "plt.xlabel(\"Iteration\")\n", + "plt.ylabel(\"MMD\")\n", "plt.ylim(0, 0.1)" ] }, @@ -426,11 +433,11 @@ "final_params = params[-1]\n", "final_st = param_to_st(final_params)\n", "\n", - "plt.hist(data, bins=50, density=True, alpha=0.5, label='Data');\n", - "plt.plot(jnp.square(jnp.abs(final_st.flatten())), label='Final Parameter')\n", + "plt.hist(data, bins=50, density=True, alpha=0.5, label=\"Data\")\n", + "plt.plot(jnp.square(jnp.abs(final_st.flatten())), label=\"Final Parameter\")\n", "plt.xlim(data.min(), data.max())\n", - "plt.ylabel('Probability')\n", - "plt.xlabel('Stamp thickness $\\\\mu$m')\n", + "plt.ylabel(\"Probability\")\n", + "plt.xlabel(\"Stamp thickness $\\\\mu$m\")\n", "plt.legend();" ] }, @@ -450,7 +457,9 @@ "outputs": [], "source": [ "burn_in = 100\n", - "av_probs = vmap(lambda p: jnp.square(jnp.abs(param_to_st(p).flatten())))(params[burn_in:]).mean(0)" + "av_probs = vmap(lambda p: jnp.square(jnp.abs(param_to_st(p).flatten())))(\n", + " params[burn_in:]\n", + ").mean(0)" ] }, { @@ -473,11 +482,11 @@ } ], "source": [ - "plt.hist(data, bins=50, density=True, alpha=0.5, label='Data');\n", - "plt.plot(av_probs, label='Averaged over parameters')\n", + "plt.hist(data, bins=50, density=True, alpha=0.5, label=\"Data\")\n", + "plt.plot(av_probs, label=\"Averaged over parameters\")\n", "plt.xlim(data.min(), data.max())\n", - "plt.ylabel('Probability')\n", - "plt.xlabel('Stamp thickness $\\\\mu$m')\n", + "plt.ylabel(\"Probability\")\n", + "plt.xlabel(\"Stamp thickness $\\\\mu$m\")\n", "plt.legend();" ] }, diff --git a/examples/heisenberg_vqe.ipynb b/examples/heisenberg_vqe.ipynb index c385c2e..0d19da0 100644 --- a/examples/heisenberg_vqe.ipynb +++ b/examples/heisenberg_vqe.ipynb @@ -43,28 +43,26 @@ "outputs": [], "source": [ "def get_circuit(n_qubits, depth):\n", - " \n", " n_params = 2 * n_qubits * (depth + 1)\n", - " \n", - " \n", - " gates = ['H'] * n_qubits + ['Rx'] * n_qubits + ['Ry'] * n_qubits\n", + "\n", + " gates = [\"H\"] * n_qubits + [\"Rx\"] * n_qubits + [\"Ry\"] * n_qubits\n", " qubit_inds = [[i] for i in range(n_qubits)] * 3\n", " param_inds = [[]] * n_qubits + [[i] for i in range(n_qubits * 2)]\n", - " \n", + "\n", " k = 2 * n_qubits\n", "\n", " for _ in range(depth):\n", " for i in range(0, n_qubits - 1):\n", - " gates.append('CZ')\n", + " gates.append(\"CZ\")\n", " qubit_inds.append([i, i + 1])\n", " param_inds.append([])\n", " for i in range(n_qubits):\n", - " gates.append('Rx')\n", + " gates.append(\"Rx\")\n", " qubit_inds.append([i])\n", " param_inds.append([k])\n", " k += 1\n", " for i in range(n_qubits):\n", - " gates.append('Ry')\n", + " gates.append(\"Ry\")\n", " qubit_inds.append([i])\n", " param_inds.append([k])\n", " k += 1\n", @@ -256,7 +254,7 @@ } ], "source": [ - "params = random.uniform(random.PRNGKey(0), shape=(n_params,), minval=0., maxval=2.)\n", + "params = random.uniform(random.PRNGKey(0), shape=(n_params,), minval=0.0, maxval=2.0)\n", "statetensor = param_to_st(params)\n", "statetensor.shape" ] @@ -395,9 +393,11 @@ }, "outputs": [], "source": [ - "hamiltonian_gates = [['X', 'X'], ['Y', 'Y'], ['Z', 'Z']] * (n_qubits - 1)\n", - "hamiltonian_qubit_inds = [[int(i), int(i) + 1] for i in jnp.repeat(jnp.arange(n_qubits), 3)]\n", - "coefficients = [1.] * len(hamiltonian_qubit_inds) " + "hamiltonian_gates = [[\"X\", \"X\"], [\"Y\", \"Y\"], [\"Z\", \"Z\"]] * (n_qubits - 1)\n", + "hamiltonian_qubit_inds = [\n", + " [int(i), int(i) + 1] for i in jnp.repeat(jnp.arange(n_qubits), 3)\n", + "]\n", + "coefficients = [1.0] * len(hamiltonian_qubit_inds)" ] }, { @@ -421,9 +421,9 @@ } ], "source": [ - "print('Gates:\\t', hamiltonian_gates)\n", - "print('Qubits:\\t', hamiltonian_qubit_inds)\n", - "print('Coefficients:\\t', coefficients)" + "print(\"Gates:\\t\", hamiltonian_gates)\n", + "print(\"Qubits:\\t\", hamiltonian_qubit_inds)\n", + "print(\"Coefficients:\\t\", coefficients)" ] }, { @@ -449,9 +449,9 @@ }, "outputs": [], "source": [ - "st_to_expectation = qujax.get_statetensor_to_expectation_func(hamiltonian_gates,\n", - " hamiltonian_qubit_inds,\n", - " coefficients)" + "st_to_expectation = qujax.get_statetensor_to_expectation_func(\n", + " hamiltonian_gates, hamiltonian_qubit_inds, coefficients\n", + ")" ] }, { @@ -577,7 +577,9 @@ } ], "source": [ - "new_params = random.uniform(random.PRNGKey(1), shape=(n_params,), minval=0., maxval=2.)\n", + "new_params = random.uniform(\n", + " random.PRNGKey(1), shape=(n_params,), minval=0.0, maxval=2.0\n", + ")\n", "param_to_expectation(new_params)" ] }, @@ -680,10 +682,10 @@ "def vqe(init_param, n_steps, stepsize):\n", " params = jnp.zeros((n_steps, n_params))\n", " params = params.at[0].set(init_param)\n", - " \n", + "\n", " cost_vals = jnp.zeros(n_steps)\n", " cost_vals = cost_vals.at[0].set(param_to_expectation(init_param))\n", - " \n", + "\n", " for step in range(1, n_steps):\n", " cost_val, cost_grad = cost_and_grad(params[step - 1])\n", " cost_vals = cost_vals.at[step].set(cost_val)\n", @@ -691,8 +693,8 @@ " new_param = params[step - 1] - stepsize * cost_grad\n", " params = params.at[step].set(new_param)\n", "\n", - " print('Iteration:', step, '\\tCost:', cost_val, end='\\r')\n", - " print('\\n')\n", + " print(\"Iteration:\", step, \"\\tCost:\", cost_val, end=\"\\r\")\n", + " print(\"\\n\")\n", " return params, cost_vals" ] }, @@ -768,8 +770,8 @@ ], "source": [ "plt.plot(vqe_cost_vals)\n", - "plt.xlabel('Iteration')\n", - "plt.ylabel('Cost');" + "plt.xlabel(\"Iteration\")\n", + "plt.ylabel(\"Cost\");" ] }, { @@ -886,8 +888,8 @@ ], "source": [ "plt.plot(new_vqe_cost_vals)\n", - "plt.xlabel('Iteration')\n", - "plt.ylabel('Cost');" + "plt.xlabel(\"Iteration\")\n", + "plt.ylabel(\"Cost\");" ] } ], diff --git a/examples/maxcut_vqe.ipynb b/examples/maxcut_vqe.ipynb index 3a73df9..574d814 100644 --- a/examples/maxcut_vqe.ipynb +++ b/examples/maxcut_vqe.ipynb @@ -75,7 +75,9 @@ "metadata": {}, "outputs": [], "source": [ - "nx_graph = networkx.generators.random_graphs.random_regular_graph(degree, n_nodes, seed=int(graph_key[-1]))\n", + "nx_graph = networkx.generators.random_graphs.random_regular_graph(\n", + " degree, n_nodes, seed=int(graph_key[-1])\n", + ")\n", "edge_inds = jnp.array(nx_graph.edges)\n", "weights = random.uniform(weights_key, shape=(len(edge_inds),))" ] @@ -100,11 +102,15 @@ "source": [ "pos = networkx.spring_layout(nx_graph)\n", "networkx.draw_networkx_nodes(nx_graph, pos)\n", - "networkx.draw_networkx_labels(nx_graph, pos, font_color='white')\n", + "networkx.draw_networkx_labels(nx_graph, pos, font_color=\"white\")\n", "networkx.draw_networkx_edges(nx_graph, pos)\n", - "networkx.draw_networkx_edge_labels(nx_graph, pos, label_pos=0.35,\n", - " edge_labels={tuple(e):f'{w:.2f}' for e, w in zip(edge_inds.tolist(), weights)})\n", - "plt.axis('off');" + "networkx.draw_networkx_edge_labels(\n", + " nx_graph,\n", + " pos,\n", + " label_pos=0.35,\n", + " edge_labels={tuple(e): f\"{w:.2f}\" for e, w in zip(edge_inds.tolist(), weights)},\n", + ")\n", + "plt.axis(\"off\");" ] }, { @@ -127,7 +133,7 @@ "\n", " def get_z_edge_weight(ei, w):\n", " i, j = ei\n", - " return jnp.where(z[i] != z[j], w, 0.)\n", + " return jnp.where(z[i] != z[j], w, 0.0)\n", "\n", " return -vmap(get_z_edge_weight)(edge_inds, weights).sum()" ] @@ -158,31 +164,29 @@ "outputs": [], "source": [ "def get_circuit(n_qubits, depth):\n", - " \n", " n_params = 2 * n_qubits * (depth + 1)\n", - " \n", - " \n", - " gates = ['H'] * n_qubits + ['Rx'] * n_qubits + ['Ry'] * n_qubits\n", + "\n", + " gates = [\"H\"] * n_qubits + [\"Rx\"] * n_qubits + [\"Ry\"] * n_qubits\n", " qubit_inds = [[i] for i in range(n_qubits)] * 3\n", " param_inds = [[]] * n_qubits + [[i] for i in range(n_qubits * 2)]\n", - " \n", + "\n", " k = 2 * n_qubits\n", "\n", " for _ in range(depth):\n", " for i in range(0, n_qubits - 1):\n", - " gates.append('CZ')\n", + " gates.append(\"CZ\")\n", " qubit_inds.append([i, i + 1])\n", " param_inds.append([])\n", - " gates.append('CZ')\n", - " qubit_inds.append([n_qubits-1, 0])\n", + " gates.append(\"CZ\")\n", + " qubit_inds.append([n_qubits - 1, 0])\n", " param_inds.append([])\n", " for i in range(n_qubits):\n", - " gates.append('Rx')\n", + " gates.append(\"Rx\")\n", " qubit_inds.append([i])\n", " param_inds.append([k])\n", " k += 1\n", " for i in range(n_qubits):\n", - " gates.append('Ry')\n", + " gates.append(\"Ry\")\n", " qubit_inds.append([i])\n", " param_inds.append([k])\n", " k += 1\n", @@ -284,7 +288,6 @@ "outputs": [], "source": [ "def param_to_sampled_cost_and_grad(param, random_key):\n", - "\n", " random_keys = random.split(random_key, 2 * param.size + 1)\n", " cost_key = random_keys[0]\n", " grad_keys = random_keys[1:]\n", @@ -296,17 +299,24 @@ " def sample_grad_k(k, random_key_plus, random_key_minus):\n", " param_plus = param.at[k].set(param[k] + 0.5)\n", " st_plus = param_to_st(param_plus)\n", - " sample_bitstrings_plus = qujax.sample_bitstrings(random_key_plus, st_plus, n_samps=n_shots)\n", + " sample_bitstrings_plus = qujax.sample_bitstrings(\n", + " random_key_plus, st_plus, n_samps=n_shots\n", + " )\n", "\n", " param_minus = param.at[k].set(param[k] - 0.5)\n", " st_minus = param_to_st(param_minus)\n", - " sample_bitstrings_minus = qujax.sample_bitstrings(random_key_minus, st_minus, n_samps=n_shots)\n", + " sample_bitstrings_minus = qujax.sample_bitstrings(\n", + " random_key_minus, st_minus, n_samps=n_shots\n", + " )\n", "\n", - " return (vmap(cost_of_bitstring)(sample_bitstrings_plus).mean() - vmap(cost_of_bitstring)(sample_bitstrings_minus).mean()) / 2\n", + " return (\n", + " vmap(cost_of_bitstring)(sample_bitstrings_plus).mean()\n", + " - vmap(cost_of_bitstring)(sample_bitstrings_minus).mean()\n", + " ) / 2\n", "\n", - " sampled_grad = vmap(sample_grad_k)(jnp.arange(param.size),\n", - " grad_keys[:param.size],\n", - " grad_keys[param.size:])\n", + " sampled_grad = vmap(sample_grad_k)(\n", + " jnp.arange(param.size), grad_keys[: param.size], grad_keys[param.size :]\n", + " )\n", " return sampled_cost, sampled_grad" ] }, @@ -337,7 +347,9 @@ "metadata": {}, "outputs": [], "source": [ - "init_param = random.uniform(init_key, shape=(n_params,), minval=-init_rad, maxval=init_rad)" + "init_param = random.uniform(\n", + " init_key, shape=(n_params,), minval=-init_rad, maxval=init_rad\n", + ")" ] }, { @@ -375,9 +387,9 @@ "metadata": {}, "outputs": [], "source": [ - "_, (params_path, cost_path) = scan(gd_iteration,\n", - " (init_param, init_opt_state, train_key),\n", - " jnp.arange(n_steps))" + "_, (params_path, cost_path) = scan(\n", + " gd_iteration, (init_param, init_opt_state, train_key), jnp.arange(n_steps)\n", + ")" ] }, { @@ -395,7 +407,9 @@ "metadata": {}, "outputs": [], "source": [ - "opt_cost = vmap(cost_of_bitstring)(qujax.integers_to_bitstrings(jnp.arange(2 ** n_qubits), n_qubits)).min()" + "opt_cost = vmap(cost_of_bitstring)(\n", + " qujax.integers_to_bitstrings(jnp.arange(2**n_qubits), n_qubits)\n", + ").min()" ] }, { @@ -418,10 +432,10 @@ } ], "source": [ - "plt.plot(cost_path, label='Training')\n", - "plt.axhline(opt_cost, color='green', label='Optimal Solution')\n", - "plt.xlabel('Iteration')\n", - "plt.ylabel('(Sampled) Maxcut Cost')\n", + "plt.plot(cost_path, label=\"Training\")\n", + "plt.axhline(opt_cost, color=\"green\", label=\"Optimal Solution\")\n", + "plt.xlabel(\"Iteration\")\n", + "plt.ylabel(\"(Sampled) Maxcut Cost\")\n", "plt.legend();" ] } diff --git a/examples/noise_channel.ipynb b/examples/noise_channel.ipynb index db74420..313df2d 100644 --- a/examples/noise_channel.ipynb +++ b/examples/noise_channel.ipynb @@ -68,7 +68,12 @@ } ], "source": [ - "paulis = {'I': qujax.gates.I, 'X': qujax.gates.X, 'Y': qujax.gates.Y, 'Z': qujax.gates.Z}\n", + "paulis = {\n", + " \"I\": qujax.gates.I,\n", + " \"X\": qujax.gates.X,\n", + " \"Y\": qujax.gates.Y,\n", + " \"Z\": qujax.gates.Z,\n", + "}\n", "paulis" ] }, @@ -87,7 +92,7 @@ } ], "source": [ - "two_paulis_strs = [a+b for a in paulis.keys() for b in paulis.keys()]\n", + "two_paulis_strs = [a + b for a in paulis.keys() for b in paulis.keys()]\n", "print(two_paulis_strs)" ] }, @@ -119,7 +124,9 @@ } ], "source": [ - "two_paulis_matrices = jnp.array([jnp.kron(paulis[a], paulis[b]) for a, b in two_paulis_strs])\n", + "two_paulis_matrices = jnp.array(\n", + " [jnp.kron(paulis[a], paulis[b]) for a, b in two_paulis_strs]\n", + ")\n", "two_paulis_matrices.shape" ] }, @@ -211,7 +218,7 @@ "metadata": {}, "outputs": [], "source": [ - "noiseless_gate_seq = ['H', 'CX', 'CX', 'CX']\n", + "noiseless_gate_seq = [\"H\", \"CX\", \"CX\", \"CX\"]\n", "noiseless_qubit_seq = [[0], [0, 1], [1, 2], [2, 3]]\n", "noiseless_param_ind_seq = [None] * len(noiseless_qubit_seq)" ] @@ -255,9 +262,12 @@ "metadata": {}, "outputs": [], "source": [ - "noisy_gate_seq = [depolarise(getattr(qujax.gates, noiseless_gate_seq[i])) if len(noiseless_qubit_seq[i]) == 2\n", - " else noiseless_gate_seq[i]\n", - " for i in range(len(noiseless_gate_seq))]" + "noisy_gate_seq = [\n", + " depolarise(getattr(qujax.gates, noiseless_gate_seq[i]))\n", + " if len(noiseless_qubit_seq[i]) == 2\n", + " else noiseless_gate_seq[i]\n", + " for i in range(len(noiseless_gate_seq))\n", + "]" ] }, { @@ -287,9 +297,10 @@ "metadata": {}, "outputs": [], "source": [ - "noisy_param_ind_seq = [p_inds if len(noiseless_qubit_seq[i]) == 2\n", - " else None\n", - " for i in range(len(noiseless_gate_seq))]" + "noisy_param_ind_seq = [\n", + " p_inds if len(noiseless_qubit_seq[i]) == 2 else None\n", + " for i in range(len(noiseless_gate_seq))\n", + "]" ] }, { @@ -333,7 +344,9 @@ "metadata": {}, "outputs": [], "source": [ - "pvec_to_densitytensor = qujax.get_params_to_densitytensor_func(noisy_gate_seq, noisy_qubit_seq, noisy_param_ind_seq)" + "pvec_to_densitytensor = qujax.get_params_to_densitytensor_func(\n", + " noisy_gate_seq, noisy_qubit_seq, noisy_param_ind_seq\n", + ")" ] }, { @@ -346,7 +359,7 @@ "def pvec_to_measure_probs(pvec):\n", " dt = pvec_to_densitytensor(pvec)\n", " nqubits = dt.ndim // 2\n", - " dm = dt.reshape(2 ** nqubits, 2 ** nqubits)\n", + " dm = dt.reshape(2**nqubits, 2**nqubits)\n", " return jnp.diag(dm).real" ] }, @@ -367,7 +380,7 @@ "metadata": {}, "outputs": [], "source": [ - "kraus_strs = ['U'] + two_paulis_strs" + "kraus_strs = [\"U\"] + two_paulis_strs" ] }, { @@ -390,8 +403,8 @@ "source": [ "noise_pvec = random.dirichlet(random.PRNGKey(0), jnp.ones(num_pauli_combos + 1))\n", "\n", - "plt.bar(kraus_strs, noise_pvec, color='purple')\n", - "plt.ylabel('p');" + "plt.bar(kraus_strs, noise_pvec, color=\"purple\")\n", + "plt.ylabel(\"p\");" ] }, { @@ -409,7 +422,7 @@ "metadata": {}, "outputs": [], "source": [ - "noiseless_pvec = jnp.zeros(num_pauli_combos+1).at[0].set(1)" + "noiseless_pvec = jnp.zeros(num_pauli_combos + 1).at[0].set(1)" ] }, { @@ -434,13 +447,24 @@ "\n", "noiseless_measure_probs = pvec_to_measure_probs(noiseless_pvec)\n", "noisy_measure_probs = pvec_to_measure_probs(noise_pvec)\n", - "measurements = random.choice(random.PRNGKey(1), len(noiseless_measure_probs), shape=(n_measurements,), p=noisy_measure_probs)\n", + "measurements = random.choice(\n", + " random.PRNGKey(1),\n", + " len(noiseless_measure_probs),\n", + " shape=(n_measurements,),\n", + " p=noisy_measure_probs,\n", + ")\n", "\n", "_x = jnp.arange(noiseless_measure_probs.size)\n", "bar_width = 0.2\n", - "plt.bar(_x - 0.2, noiseless_measure_probs, bar_width, label='Noiseless', color='grey')\n", - "plt.bar(_x, noisy_measure_probs, bar_width, label='Noisy', color='green')\n", - "plt.bar(_x + 0.2, jnp.bincount(measurements) / n_measurements, bar_width, label='Sampled', color='orange')\n", + "plt.bar(_x - 0.2, noiseless_measure_probs, bar_width, label=\"Noiseless\", color=\"grey\")\n", + "plt.bar(_x, noisy_measure_probs, bar_width, label=\"Noisy\", color=\"green\")\n", + "plt.bar(\n", + " _x + 0.2,\n", + " jnp.bincount(measurements) / n_measurements,\n", + " bar_width,\n", + " label=\"Sampled\",\n", + " color=\"orange\",\n", + ")\n", "plt.xticks(_x, [f\"{a:04b}\" for a in _x], rotation=90)\n", "plt.legend();" ] @@ -497,15 +521,18 @@ "def logit(u):\n", " return jnp.log(u) - jnp.log(1 - u)\n", "\n", + "\n", "def inverse_logit(v):\n", " return 1 / (1 + jnp.exp(-v))\n", "\n", + "\n", "def constrain_simplex(y):\n", " k = y.shape[-1] + 1\n", " z = inverse_logit(y - jnp.log(k - jnp.arange(1, k)))\n", - " z1m_cumprod = jnp.append(1., (1 - z[:-1]).cumprod())\n", + " z1m_cumprod = jnp.append(1.0, (1 - z[:-1]).cumprod())\n", " return z * z1m_cumprod\n", "\n", + "\n", "def unconstrain_simplex(x):\n", " k = x.shape[-1] + 1\n", " z = x / (1 - jnp.append(0, x[:-1].cumsum()))\n", @@ -529,7 +556,7 @@ "source": [ "def transformed_log_likelihood(z):\n", " pvec_trunc = constrain_simplex(z)\n", - " pvec = jnp.append(pvec_trunc, 1-pvec_trunc.sum())\n", + " pvec = jnp.append(pvec_trunc, 1 - pvec_trunc.sum())\n", " return log_likelihood(pvec)" ] }, @@ -1081,8 +1108,8 @@ "for i in range(n_steps):\n", " loglik, tloglik_grad = transformed_log_likelihood_and_grad(tp_train)\n", " logliks_train = logliks_train.at[i].set(loglik)\n", - " print('Step:', i+1, '\\t Log-likelihood', loglik, end='\\r')\n", - " \n", + " print(\"Step:\", i + 1, \"\\t Log-likelihood\", loglik, end=\"\\r\")\n", + "\n", " tp_train += stepsize * tloglik_grad" ] }, @@ -1105,8 +1132,8 @@ ], "source": [ "plt.plot(logliks_train)\n", - "plt.xlabel('Iteration')\n", - "plt.ylabel('Log-likelihood');" + "plt.xlabel(\"Iteration\")\n", + "plt.ylabel(\"Log-likelihood\");" ] }, { @@ -1127,7 +1154,7 @@ "outputs": [], "source": [ "mle_pvec = constrain_simplex(tp_train)\n", - "mle_pvec = jnp.append(mle_pvec, 1-mle_pvec.sum())" + "mle_pvec = jnp.append(mle_pvec, 1 - mle_pvec.sum())" ] }, { @@ -1150,10 +1177,10 @@ "source": [ "_x2 = jnp.arange(len(kraus_strs))\n", "bar_width2 = 0.4\n", - "plt.bar(_x2 - 0.2, noise_pvec, bar_width2, label='True', color='purple')\n", - "plt.bar(_x2 + 0.2, mle_pvec, bar_width2, label='Fitted MLE', color='lightblue')\n", + "plt.bar(_x2 - 0.2, noise_pvec, bar_width2, label=\"True\", color=\"purple\")\n", + "plt.bar(_x2 + 0.2, mle_pvec, bar_width2, label=\"Fitted MLE\", color=\"lightblue\")\n", "plt.xticks(_x2, kraus_strs)\n", - "plt.ylabel('p')\n", + "plt.ylabel(\"p\")\n", "plt.legend();" ] }, @@ -1187,11 +1214,17 @@ "source": [ "mle_measure_probs = pvec_to_measure_probs(mle_pvec)\n", "\n", - "plt.bar(_x - 0.2, noisy_measure_probs, bar_width, label='Noisy', color='green')\n", - "plt.bar(_x, jnp.bincount(measurements) / n_measurements, bar_width, label='Sampled', color='orange')\n", - "plt.bar(_x + 0.2, mle_measure_probs, bar_width, label='MLE', color='steelblue')\n", + "plt.bar(_x - 0.2, noisy_measure_probs, bar_width, label=\"Noisy\", color=\"green\")\n", + "plt.bar(\n", + " _x,\n", + " jnp.bincount(measurements) / n_measurements,\n", + " bar_width,\n", + " label=\"Sampled\",\n", + " color=\"orange\",\n", + ")\n", + "plt.bar(_x + 0.2, mle_measure_probs, bar_width, label=\"MLE\", color=\"steelblue\")\n", "plt.xticks(_x, [f\"{a:04b}\" for a in _x], rotation=90)\n", - "plt.ylabel(r'$P(m \\mid p)$')\n", + "plt.ylabel(r\"$P(m \\mid p)$\")\n", "plt.legend();" ] }, diff --git a/examples/qaoa.ipynb b/examples/qaoa.ipynb index 5976dd4..7a775ca 100644 --- a/examples/qaoa.ipynb +++ b/examples/qaoa.ipynb @@ -52,13 +52,13 @@ "source": [ "n_qubits = 4\n", "hamiltonian_qubit_inds = [(0, 1), (1, 2), (0, 2), (1, 3)]\n", - "hamiltonian_gates = [['Z', 'Z']] * (len(hamiltonian_qubit_inds))\n", + "hamiltonian_gates = [[\"Z\", \"Z\"]] * (len(hamiltonian_qubit_inds))\n", "\n", "\n", "# Notice that in order to use the random package from jax we first need to define a seeded key\n", "seed = 13\n", "key = random.PRNGKey(seed)\n", - "coefficients = random.uniform(key, shape=(len(hamiltonian_qubit_inds), ))" + "coefficients = random.uniform(key, shape=(len(hamiltonian_qubit_inds),))" ] }, { @@ -78,9 +78,9 @@ } ], "source": [ - "print('Gates:\\t', hamiltonian_gates)\n", - "print('Qubits:\\t', hamiltonian_qubit_inds)\n", - "print('Coefficients:\\t', coefficients)" + "print(\"Gates:\\t\", hamiltonian_gates)\n", + "print(\"Qubits:\\t\", hamiltonian_qubit_inds)\n", + "print(\"Coefficients:\\t\", coefficients)" ] }, { @@ -119,21 +119,20 @@ "source": [ "circuit_gates = []\n", "circuit_qubit_inds = []\n", - "circuit_param_inds =[]\n", + "circuit_param_inds = []\n", "\n", "param_ind = 0\n", "\n", "# Initial State\n", "for i in range(n_qubits):\n", - " circuit_gates.append('H')\n", + " circuit_gates.append(\"H\")\n", " circuit_qubit_inds.append([i])\n", " circuit_param_inds.append([])\n", - " \n", + "\n", "for d in range(depth):\n", - " \n", " # Mixing Unitary\n", " for i in range(n_qubits):\n", - " circuit_gates.append('Rx')\n", + " circuit_gates.append(\"Rx\")\n", " circuit_qubit_inds.append([i])\n", " circuit_param_inds.append([param_ind])\n", " param_ind += 1\n", @@ -142,16 +141,16 @@ " for index in range(len(hamiltonian_qubit_inds)):\n", " pair = hamiltonian_qubit_inds[index]\n", " coef = coefficients[index]\n", - " \n", - " circuit_gates.append('CX')\n", + "\n", + " circuit_gates.append(\"CX\")\n", " circuit_qubit_inds.append([pair[0], pair[1]])\n", " circuit_param_inds.append([])\n", - " \n", + "\n", " circuit_gates.append(lambda p: qujax.gates.Rz(p * coef))\n", " circuit_qubit_inds.append([pair[1]])\n", " circuit_param_inds.append([param_ind])\n", - " \n", - " circuit_gates.append('CX')\n", + "\n", + " circuit_gates.append(\"CX\")\n", " circuit_qubit_inds.append([pair[0], pair[1]])\n", " circuit_param_inds.append([])\n", " param_ind += 1" @@ -188,8 +187,9 @@ } ], "source": [ - "qujax.print_circuit(circuit_gates, circuit_qubit_inds, circuit_param_inds,\n", - " gate_ind_max=20);" + "qujax.print_circuit(\n", + " circuit_gates, circuit_qubit_inds, circuit_param_inds, gate_ind_max=20\n", + ");" ] }, { @@ -207,7 +207,9 @@ "metadata": {}, "outputs": [], "source": [ - "param_to_st = qujax.get_params_to_statetensor_func(circuit_gates, circuit_qubit_inds, circuit_param_inds)" + "param_to_st = qujax.get_params_to_statetensor_func(\n", + " circuit_gates, circuit_qubit_inds, circuit_param_inds\n", + ")" ] }, { @@ -225,9 +227,9 @@ "metadata": {}, "outputs": [], "source": [ - "st_to_expectation = qujax.get_statetensor_to_expectation_func(hamiltonian_gates,\n", - " hamiltonian_qubit_inds,\n", - " coefficients)\n", + "st_to_expectation = qujax.get_statetensor_to_expectation_func(\n", + " hamiltonian_gates, hamiltonian_qubit_inds, coefficients\n", + ")\n", "\n", "param_to_expectation = lambda param: st_to_expectation(param_to_st(param))" ] @@ -268,7 +270,7 @@ "source": [ "seed = 123\n", "key = random.PRNGKey(seed)\n", - "init_param = random.uniform(key, shape=(param_ind, ))\n", + "init_param = random.uniform(key, shape=(param_ind,))\n", "\n", "n_steps = 150\n", "stepsize = 0.01" @@ -446,7 +448,7 @@ " cost_val, cost_grad = cost_and_grad(param)\n", " cost_vals = cost_vals.at[step].set(cost_val)\n", " param = param - stepsize * cost_grad\n", - " print('Iteration:', step, '\\tCost:', cost_val, end='\\r')" + " print(\"Iteration:\", step, \"\\tCost:\", cost_val, end=\"\\r\")" ] }, { @@ -478,8 +480,8 @@ ], "source": [ "plt.plot(cost_vals)\n", - "plt.xlabel('Iteration')\n", - "plt.ylabel('Cost');" + "plt.xlabel(\"Iteration\")\n", + "plt.ylabel(\"Cost\");" ] } ], diff --git a/examples/reducing_jit_compilation_time.ipynb b/examples/reducing_jit_compilation_time.ipynb new file mode 100644 index 0000000..c605f78 --- /dev/null +++ b/examples/reducing_jit_compilation_time.ipynb @@ -0,0 +1,310 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Reducing JIT Compilation time" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## The Problem" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "One characteristic of `jax.jit` is that it asks to XLA to compile every operation it traces. This can lead to very long compilation times when working with large quantum circuits.\n", + "\n", + "If the circuit has a repeating structure, as is common in e.g. parameterized quantum circuits, we can simply ask JAX to compile a single repetition of the circuit and reuse it multiple times.\n", + "\n", + "`qujax` makes this easy through a convenience function called `repeat_circuit`. This notebook illustrates how to use this function and shows what happens when we do not." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Imports" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "No GPU/TPU found, falling back to CPU. (Set TF_CPP_MIN_LOG_LEVEL=0 and rerun for more info.)\n" + ] + } + ], + "source": [ + "from typing import List, Tuple\n", + "import time\n", + "\n", + "import jax\n", + "from jax.random import PRNGKey\n", + "\n", + "import qujax\n", + "from qujax import print_circuit, repeat_circuit, all_zeros_statetensor" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Circuit definition" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "def circuit(\n", + " n_qubits: int, depth: int\n", + ") -> Tuple[List[str], List[List[int]], List[List[int]], int]:\n", + " gates: List[str] = []\n", + " qubit_inds: List[List[int]] = []\n", + " param_inds: List[List[int]] = []\n", + "\n", + " parameter_index = 0\n", + "\n", + " for _ in range(depth):\n", + " # Rx layer\n", + " for i in range(n_qubits):\n", + " gates.append(\"Rx\")\n", + " qubit_inds.append([i])\n", + " param_inds.append([parameter_index])\n", + " parameter_index += 1\n", + "\n", + " # Rz layer\n", + " for i in range(0, n_qubits):\n", + " gates.append(\"Rz\")\n", + " qubit_inds.append([i])\n", + " param_inds.append([parameter_index])\n", + " parameter_index += 1\n", + "\n", + " # CRz layer\n", + " for i in range(n_qubits - 1):\n", + " gates.append(\"CRz\")\n", + " qubit_inds.append([i, i + 1])\n", + " param_inds.append([parameter_index])\n", + " parameter_index += 1\n", + "\n", + " return gates, qubit_inds, param_inds, parameter_index" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Print one repetition of the circuit to visually check for correctness" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [], + "source": [ + "gates, qubit_inds, param_inds, nr_of_parameters = circuit(4, 1)" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "q0: ---Rx[0]---Rz[4]-----◯-------------------\n", + " | \n", + "q1: ---Rx[1]---Rz[5]---CRz[8]----◯-----------\n", + " | \n", + "q2: ---Rx[2]---Rz[6]-----------CRz[9]----◯---\n", + " | \n", + "q3: ---Rx[3]---Rz[7]------------------CRz[10]\n" + ] + } + ], + "source": [ + "print_circuit(gates, qubit_inds, param_inds);" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Compute compilation times" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [], + "source": [ + "n_qubits = 12\n", + "circuit_depths = list(range(4, 25, 4))\n", + "rng_seed = 0" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "**Using** the `repeat_circuit` function" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [], + "source": [ + "optimised_times_list = []\n", + "rng_key = PRNGKey(rng_seed)\n", + "\n", + "for d in circuit_depths:\n", + " # qujax specification of circuit\n", + " gates, qubit_inds, param_inds, nr_of_parameters = circuit(n_qubits, 1)\n", + "\n", + " # Get function that returns one application of the circuit\n", + " params_to_statetensor = qujax.get_params_to_statetensor_func(\n", + " gates, qubit_inds, param_inds\n", + " )\n", + "\n", + " rng_key, parameters_rng = jax.random.split(rng_key)\n", + "\n", + " random_angles = jax.random.uniform(parameters_rng, (d * nr_of_parameters,))\n", + "\n", + " # Allow for an arbitrary number of circuit repetitions while avoiding compilation overhead\n", + " repeated_circuit = jax.jit(repeat_circuit(params_to_statetensor, nr_of_parameters))\n", + "\n", + " initial_state = all_zeros_statetensor(n_qubits)\n", + "\n", + " start = time.time()\n", + " repeated_circuit.lower(random_angles, initial_state).compile()\n", + " elapsed = time.time() - start\n", + "\n", + " optimised_times_list.append(elapsed)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "**Not using** the `repeat_circuit` function" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [], + "source": [ + "unoptimised_times_list = []\n", + "rng_key = PRNGKey(rng_seed)\n", + "\n", + "for d in circuit_depths:\n", + " # qujax specification of circuit\n", + " gates, qubit_inds, param_inds, nr_of_parameters = circuit(n_qubits, d)\n", + "\n", + " params_to_statetensor = jax.jit(\n", + " qujax.get_params_to_statetensor_func(gates, qubit_inds, param_inds)\n", + " )\n", + "\n", + " rng_key, parameters_rng = jax.random.split(rng_key)\n", + "\n", + " random_angles = jax.random.uniform(parameters_rng, (d * nr_of_parameters,))\n", + "\n", + " start = time.time()\n", + " params_to_statetensor.lower(random_angles).compile()\n", + " elapsed = time.time() - start\n", + "\n", + " unoptimised_times_list.append(elapsed)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Plot results" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [], + "source": [ + "import matplotlib.pyplot as plt" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkIAAAGwCAYAAABFFQqPAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABZpklEQVR4nO3dd1xTV/8H8E9YYQhBVDYIKlVxKyqoddTRR637eepspdq6a9W2tnaJHa6nVdvaulpHh9W2jvr8bB21ilqkooJ7oOJiiAMSZhg5vz8igTAUAmR+3q8XLTl35HvvTcyXc773RCKEECAiIiKyQFaGDoCIiIjIUJgIERERkcViIkREREQWi4kQERERWSwmQkRERGSxmAgRERGRxWIiRERERBbLxtAB1DaVSoWkpCQ4OztDIpEYOhwiIiKqBCEEMjIy4O3tDSur2uu3MftEKCkpCX5+foYOg4iIiHRw+/Zt+Pr61tr+zT4RcnZ2BqA+kS4uLgaOhoiIiCpDoVDAz89P8zleW8w+ESoaDnNxcWEiREREZGJqu6yFxdJERERksZgIERERkcViIkREREQWy+xrhCqrsLAQ+fn5hg6DjIStrS2sra0NHQYREdUyi0+EhBBISUlBenq6oUMhI+Pq6gpPT0/OP0VEZMYsPhEqSoLc3d3h6OjIDz2CEALZ2dlITU0FAHh5eRk4IiIiqi0WnQgVFhZqkqB69eoZOhwyIg4ODgCA1NRUuLu7c5iMiMhMWXSxdFFNkKOjo4EjIWNU9Lpg7RgRkfmy6ESoCIfDqDx8XRARmT8mQkRERFRpitx8JMtzyl2WLM+BIte0etGZCBEREVGlKHLzMX79cYxcE42kdO1kKCk9ByPXRGP8+uMmlQwxEaIat3HjRri6ulZ7PxKJBDt37qz2fioSERGBtm3b1tr+iYjMTZayAA8y83DrYTZGrS1OhpLSczBqbTRuPczGg8w8ZCkLDBxp5Rk0EYqIiIBEItH68fT01CwXQiAiIgLe3t5wcHBAz549cf78eQNGrM3cugd1ERAQgBUrVmi1jRw5EleuXKn2vpOTk9G/f/9q74eIiGqGl8wBWyaFwt/NEbceZmPkmmNYE3lNkwT5uzliy6RQeMkcDB1qpRm8R6hFixZITk7W/Jw9e1azbOnSpVi2bBlWrlyJmJgYeHp6om/fvsjIyDBgxGrm2D1YUxwcHODu7l7t/Xh6ekIqldZAREREVFO8XYuTodtpOVj0xyWtJMjb1XSSIMAIEiEbGxt4enpqfho0aABA3Ru0YsUKvPvuuxg+fDhatmyJTZs2ITs7G5s3bzZw1IbvHiyvJ6Zt27aIiIgAoB5W+uabbzBs2DA4OjoiKCgIu3bt0lo/MjISnTp1glQqhZeXF95++20UFBTH27NnT8yYMQMzZsyAq6sr6tWrh/feew9CCM3ymzdvYvbs2ZoePaDs0FjRENT69evh7++POnXqYOrUqSgsLMTSpUvh6ekJd3d3fPLJJ1rxlRway8vLw4wZM+Dl5QV7e3sEBARg0aJFmnXlcjkmTZoEd3d3uLi44JlnnsHp06e19rd48WJ4eHjA2dkZEydORG5ubpXPOxERqZOhkR19tdqWj2xjckkQYASJUHx8PLy9vREYGIhRo0bh+vXrAICEhASkpKSgX79+mnWlUil69OiBqKioCvenVCqhUCi0fmpD6e7BUWujcfLmQ6PqHlywYAGef/55nDlzBgMGDMDYsWPx8OFDAEBiYiIGDBiAjh074vTp01i1ahW+/fZbfPzxx1r72LRpE2xsbPDPP//giy++wPLly/HNN98AALZv3w5fX198+OGHmh69ily7dg1//PEH9uzZg59++gnr16/HwIEDcefOHURGRmLJkiV47733EB0dXe72X3zxBXbt2oWff/4Zly9fxg8//ICAgAAA6qR54MCBSElJwe+//46TJ0+iffv26N27t+Z4f/75Z8yfPx+ffPIJTpw4AS8vL3z99dfVPcVERBbpaPw9fLpXuwRi9tbTZUZITIFBE6HOnTvju+++w969e7Fu3TqkpKSgS5cuePDgAVJSUgAAHh4eWtt4eHholpVn0aJFkMlkmh8/P79ai79k9+Cth9kYseqYUXUPhoeHY/To0WjSpAkWLlyIrKwsHD9+HADw9ddfw8/PDytXrkSzZs0wdOhQLFiwAJ999hlUKpVmH35+fli+fDmaNm2KsWPH4tVXX8Xy5csBAG5ubrC2toazs7OmR68iKpUK69evR3BwMAYNGoRevXrh8uXLWLFiBZo2bYqXXnoJTZs2xaFDh8rd/tatWwgKCkK3bt3QsGFDdOvWDaNHjwYAHDx4EGfPnsUvv/yCkJAQBAUF4dNPP4Wrqyt+/fVXAMCKFSswYcIEvPzyy2jatCk+/vhjBAcH18RpJiKyKJdTFAjfEAMBQGpjhZ8na3cKmFoyZNBEqH///hgxYgRatWqFPn36YPfu3QDUvRBFSk9qJ4R47ER38+bNg1wu1/zcvn27doJ/xNvVActHttFqM5buwdatW2t+d3JygrOzs+b7sy5evIiwsDCtc9m1a1dkZmbizp07mrbQ0FCtdcLCwhAfH4/CwsIqxRIQEABnZ2fNYw8PDwQHB8PKykqrrSi+0sLDwxEXF4emTZti5syZ2Ldvn2bZyZMnkZmZiXr16qFOnTqan4SEBFy7dk3reEsq/ZiIiB7vTlo2hnz1NwpUAtZWEuyY1gWdAuuVGSGp6EYiY2RU3zXm5OSEVq1aIT4+HkOHDgWg/lLUkl96mZqaWqaXqCSpVKrXAtuk9BzM3qpdizJ76+la7xGysrLS1OoUKf1VELa2tlqPJRKJprenvISyaH+1MaNyebE8Lr7S2rdvj4SEBPzxxx/4888/8fzzz6NPnz749ddfoVKp4OXlVW5vUk3cxk9ERGqbom4gN18FCYBvXgxBsLcMQPEIyai10ahXxw5OUqNKLx7L4DVCJSmVSly8eBFeXl4IDAyEp6cn9u/fr1mel5eHyMhIdOnSxYBRFitZGO3v5ohtU8P01j3YoEEDrZochUKBhISESm8fHByMqKgorWQqKioKzs7O8PHx0bSVrtmJjo5GUFCQ5ktI7ezsqtw7pCsXFxeMHDkS69atw9atW7Ft2zY8fPgQ7du3R0pKCmxsbNCkSROtn/r16wMAmjdvXu6xEBFR5ew5l4x1R9SfM/MHBaNXM+27g71dHbB1cig2TegEF3vb8nZhlAyaCL3xxhuIjIxEQkIC/vnnH/z73/+GQqHA+PHjIZFIMGvWLCxcuBA7duzAuXPnEB4eDkdHR4wZM8aQYQNQzxNUujC6Q0M3vXUPPvPMM/j+++9x5MgRnDt3DuPHj6/SN6RPmzYNt2/fxquvvopLly7ht99+w/z58zFnzhyt4arbt29jzpw5uHz5Mn766Sd8+eWXeO211zTLAwICcPjwYSQmJuL+/fs1eowlLV++HFu2bMGlS5dw5coV/PLLL/D09ISrqyv69OmDsLAwDB06FHv37sWNGzcQFRWF9957DydOnAAAvPbaa1i/fj3Wr1+PK1euYP78+UY1JxURkTGLv5uB139Wj35M7BaI8K6B5a7nJXMwqSQIMPDQ2J07dzB69Gjcv38fDRo0QGhoKKKjo9GwYUMAwNy5c5GTk4Np06YhLS0NnTt3xr59+7RqTQzFSWqDenXsAEBrGExf3YPz5s3D9evX8dxzz0Emk+Gjjz6qUo+Qj48Pfv/9d7z55pto06YN3NzcMHHiRLz33nta67344ovIyclBp06dYG1tjVdffRWTJk3SLP/www8xefJkNG7cGEqlssxwXU2pU6cOlixZgvj4eFhbW6Njx474/fffNUnb77//jnfffRcTJkzAvXv34Onpie7du2uGUUeOHIlr167hrbfeQm5uLkaMGIGpU6di7969tRIvEZG5kOfkY9L3J5GVV4gujethXv9mhg6pRklEbX1yGQmFQgGZTAa5XA4XFxetZbm5uUhISEBgYCDs7e2rvu/cfGQpC8q9RT5ZngMnqY3JZcYl9ezZE23bti0zX5GlqO7rg4jI1KlUAhM3xeDg5XvwcXXArhldUa+OfupwH/f5XZNMp5rJCLnY21aY6JjS9OJERETlWf7nFRy8fA9SGyuseaGD3pIgfTKqYmkiIiIyDnvOJePLv64CABaPaIWWPjIDR1Q72CNEFapockMiIjJvJYujJ3QNxLB2vk/YwnSxR4iIiIg0ShZHhzWqh3cGmFdxdGlMhIiIiAiAujh69tY4JNzPgo+rA1aOaQcba/NOFcz76IiIiKjSVvx5BX9dSjXr4ujSmAgRERER9pxLwRcWUBxdGhMhIiIiC3c1NQOv/xwHwPyLo0tjIkRaAgICqj2BYkREBNq2bVsj8VREIpFg586dtfocRESWQJGbj0nfqYujQxu5YZ6ZF0eXxkTIQm3cuLHcb2aPiYnR+goNXbzxxhs4cOBAtfZBRES1T6USmL0lDtfvZ8FbZo+vxrSHrZkXR5fGeYSqI1cOKDMBmU/ZZfJEQFoHsDetMdYGDRpUex916tRBnTp1aiAaIiKqTSsOxOOApjg6xCKKo0uzrLSvJuXKgR9GABsHAPI72svkd9TtP4xQr1cLlEolZs6cCXd3d9jb26Nbt26IiYkBoJ4IUSKRYPfu3WjTpg3s7e3RuXNnnD17VrP8pZdeglwuh0QigUQiQUREBICyQ2MSiQRr1qzBc889B0dHRzRv3hzHjh3D1atX0bNnTzg5OSEsLAzXrl3TbFN6aOzQoUPo1KkTnJyc4Orqiq5du+LmzZua5f/73//QoUMH2Nvbo1GjRliwYAEKCgo0y+Pj49G9e3fY29sjODgY+/fvr4UzSkRkWfaeT8EXB+IBAIuGt0IrX9P6w72mMBHSlTITyLoHpN0ANg4sTobkd9SP026olysza+Xp586di23btmHTpk04deoUmjRpgmeffRYPHz7UrPPmm2/i008/RUxMDNzd3TF48GDk5+ejS5cuWLFiBVxcXJCcnIzk5GS88cYbFT7XRx99hBdffBFxcXFo1qwZxowZg8mTJ2PevHk4ceIEAGDGjBnlbltQUIChQ4eiR48eOHPmDI4dO4ZJkyZBIpEAAPbu3Ytx48Zh5syZuHDhAtasWYONGzfik08+AQCoVCoMHz4c1tbWiI6OxurVq/HWW2/V1GkkIrJIV1MzMGdrHADgpa4BGN7ecoqjyxBmTi6XCwBCLpeXWZaTkyMuXLggcnJydNt5+m0hVrQWYr6L+v83o7Ufp9+uZvTly8zMFLa2tuLHH3/UtOXl5Qlvb2+xdOlScfDgQQFAbNmyRbP8wYMHwsHBQWzdulUIIcSGDRuETCYrs++GDRuK5cuXax4DEO+9957m8bFjxwQA8e2332rafvrpJ2Fvb695PH/+fNGmTRvN8wIQhw4dKvdYnn76abFw4UKttu+//154eXkJIYTYu3evsLa2FrdvF5/LP/74QwAQO3bsqOAM1Yxqvz6IiIyQPCdP9PrvQdHwrf8Tz6+OEnkFhYYOqVyP+/yuSewRqg6ZLxC+G6gboO4BWt9P/f+6Aep2We1k2NeuXUN+fj66du2qabO1tUWnTp1w8eJFTVtYWJjmdzc3NzRt2lRreWW1bt1a87uHhwcAoFWrVlptubm5UCgUZbZ1c3NDeHg4nn32WQwaNAiff/45kpOTNctPnjyJDz/8UFNXVKdOHbzyyitITk5GdnY2Ll68CH9/f/j6Fp/LksdFRESVV6Y4eqzlFUeXZtlHXxNkvsCwtdptw9bWWhIEAEIIANAML5VsL91W2pOWl8fW1rbM9uW1qVSqcrffsGEDjh07hi5dumDr1q146qmnEB0drdlmwYIFiIuL0/ycPXsW8fHxsLe31xxrdY+BiIiKi6PtHhVH17fA4ujSmAhVl/wOsKPU7eY7JpUtoK5BTZo0gZ2dHY4ePappy8/Px4kTJ9C8eXNNW1GyAQBpaWm4cuUKmjVTzw9hZ2eHwsLCWouxtHbt2mHevHmIiopCy5YtsXnzZgBA+/btcfnyZTRp0qTMj5WVFYKDg3Hr1i0kJSVp9nXs2DG9xU1EZC72lSyOHma5xdGlMRGqjpKF0XUDgAn7iofJShZQ1zAnJydMnToVb775Jvbs2YMLFy7glVdeQXZ2NiZOnKhZ78MPP8SBAwdw7tw5hIeHo379+hg6dCgA9d1hmZmZOHDgAO7fv4/s7OxaiTUhIQHz5s3DsWPHcPPmTezbtw9XrlzRJGwffPABvvvuO0REROD8+fO4ePEitm7divfeew8A0KdPHzRt2hQvvvgiTp8+jSNHjuDdd9+tlViJiMzV1dRMzPn5NAAgvEsARnSw4OLoUpgI6UqeqJ0Ehe8G/Dtr1wxtHKherxYsXrwYI0aMwAsvvID27dvj6tWr2Lt3L+rWrau1zmuvvYYOHTogOTkZu3btgp2dHQCgS5cumDJlCkaOHIkGDRpg6dKltRKno6MjLl26hBEjRuCpp57CpEmTMGPGDEyePBkA8Oyzz+L//u//sH//fnTs2BGhoaFYtmwZGjZsCACwsrLCjh07oFQq0alTJ7z88suaO8qIiOjJ1DNHn0CmsgCdA93w7sDmT97IgkhEeUUYZkShUEAmk0Eul8PFxUVrWW5uLhISEhAYGAh7e/uq7bhoHqGse2ULo4t6ipwaAOO26X1SxUOHDqFXr15IS0srd/ZoqpxqvT6IiIyASiUw6fsT+PNiKrxl9tj1ajeTqQt63Od3TeLM0rqyl6mTnPJmlpb5AuG/m+TM0kREZD4+PxCPPy+qi6NXv9DBZJIgfWIiVB32sooTnfK+doOIiEhP9p1PwecliqNb+7oaNiAjxUTIDPXs2bPc286JiMgysDi68lgsTUREZEYUufmY9D2LoyuLiRDA3hMqF18XRGRqVCqBOVvjcP1eFrw4c3SlWPTZKZodubbm0CHTVvS6KDmLNhGRMfvir+Li6DUsjq4Ui64Rsra2hqurK1JTUwGo57zh1zeQEALZ2dlITU2Fq6srrK2tDR0SEdET7b9wFyv+VBdHL2RxdKVZdCIEAJ6engCgSYaIiri6umpeH0RExuxqaiZmb40DoC6O/jeLoyvN4hMhiUQCLy8vuLu7Iz8/39DhkJGwtbVlTxARmYSMEsXRnVgcXWUWnwgVsba25gcfERGZFJVKYPbW05ri6K9ZHF1lPFtEREQmSl0cfVc9c/Q4FkfrgokQERGRCSpZHP3J0JZo4+dq2IBMFBMhIiIiE3PtXibmPCqOHh/WEP8J8TNsQCaMiRAREZEJycjNx6TvTiBDWYBOAW5477lgQ4dk0pgIERERmQiVSmDOz6dxjTNH1xiePSIiIhPx5V9Xsf9CcXF0A2cWR1cXEyEiIiIT8OeFu1j+5xUAwMcsjq4xTISIiIiM3LV7xTNHvxjWEM+zOLrGMBEiIiIyYqWLo99ncXSNYiJERERkpFQqgdcfFUd7urA4ujbwbBIRERmplQevYl9RcfQLLI6uDUyEiIiIjNCBi9rF0W1ZHF0rmAgREREZmev3MjFrSxyEAF4IZXF0bWIiREREZEQycvMx6fuTyFAWoGNAXRZH1zImQkREREaiqDj6amomPF3s8fXYDrCz4Ud1beLZJSIiMhJfFRVHW1th1bj2LI7WAyZCRERERuCvS3exrERxdDv/ugaOyDIwESIiIjKw6/cy8dpPJYqjO7I4Wl+YCBERERlQprKAxdEGxESIiIjIQNTF0XGa4uivxrZncbSe8WwTEREZyFcHr2Lv+eLiaHdne0OHZHGYCBERERlAyeLoj4a2YHG0gTARIiIi0rOE+1l47dHM0eNC/TGyo7+hQ7JYTISIiIj0KFNZgEnfnUBGbgFCGtbFB8+1MHRIFo2JEBERkZ4UFUfHp2bCw0WKr8exONrQePaJiIj05OtDxcXRq8d1YHG0EWAiREREpAd/XbqLz/azONrYGE0itGjRIkgkEsyaNUvTJoRAREQEvL294eDggJ49e+L8+fOGC5KIiEgHJYujx3ZmcbQxMYpEKCYmBmvXrkXr1q212pcuXYply5Zh5cqViImJgaenJ/r27YuMjAwDRUpERFQ1pYuj5w9icbQxMXgilJmZibFjx2LdunWoW7e4m1AIgRUrVuDdd9/F8OHD0bJlS2zatAnZ2dnYvHmzASMmIiKqHCEE3vj5NIujjZjBr8b06dMxcOBA9OnTR6s9ISEBKSkp6Nevn6ZNKpWiR48eiIqKqnB/SqUSCoVC64eIiMgQvj50DXvOpzyaOZrF0cbIxpBPvmXLFpw6dQoxMTFllqWkpAAAPDw8tNo9PDxw8+bNCve5aNEiLFiwoGYDJSIiqqKDl1Lx6b7LAIAPh7RAexZHGyWD9Qjdvn0br732Gn744QfY21ecIUskEq3HQogybSXNmzcPcrlc83P79u0ai5mIiKgyEu5nYeaWWAgBjOnsj1GdWBxtrAzWI3Ty5EmkpqaiQ4cOmrbCwkIcPnwYK1euxOXL6iw6JSUFXl5emnVSU1PL9BKVJJVKIZVKay9wIiKixyhZHN2hYV1EsDjaqBmsR6h37944e/Ys4uLiND8hISEYO3Ys4uLi0KhRI3h6emL//v2abfLy8hAZGYkuXboYKmwiIqIKCSHw5i/q4mh3ZylWjWVxtLEzWI+Qs7MzWrZsqdXm5OSEevXqadpnzZqFhQsXIigoCEFBQVi4cCEcHR0xZswYQ4RMRET0WF8fuoY/zqXA1lqiLo52YXG0sTNosfSTzJ07Fzk5OZg2bRrS0tLQuXNn7Nu3D87OzoYOjYiISMvByyWLo1uiQ0MWR5sCiRBCGDqI2qRQKCCTySCXy+Hi4mLocIiIyAzduJ+FwSuPQpFbgDGd/bFwWCtDh2Ty9PX5zYFLIiKiashUFmDS9yegyC1Ae39XzB8UbOiQqAqYCBEREemoqDj6yl11cfTqcR0gtbE2dFhUBUyEiIiIdMTiaNPHRIiIiEgHJYujFwxmcbSpYiJERERURTfuZ+G1n9QzR4/u5I8xnTlztKliIkRERFQFWaWKoyMGszjalDERIiIiqiQhBN78VV0c3cBZilUsjjZ5TISIiIgqaVXkNfx+Vl0cvXpce3iwONrkMREiIiKqhEOXU/HfvSWLo90MHBHVBCZCRERET3DzQRZmaoqj/VgcbUaYCBERET1GlrIAk747WaI4uoWhQ6IaxESIiIioAkXF0ZfvZrA42kwxESIiIqrA6sjrLI42c0yEiIiIyhF55R6W7r0EAIgY3ILF0WaKiRAREVEpNx9k4dXNpzTF0WM7NzR0SFRLmAgRERGVULI4uh2Lo80eEyEiIqJHhBCYu+2Mpjh6NYujzR4TISIiokdWR17H7jPJsLWWYNVYFkdbAiZCRERE0C6Onj+oBUICWBxtCZgIERGRxSs5c/Sojn4Yy5mjLQYTISIismjZeQWY/P1JyHPy0dbPFQuGtIBEIjF0WKQnTISIiMhiqWeOPoNLKSyOtlRMhIiIyGKtOawujraxUhdHe8pYHG1pmAgREZFFOnzlHpbueVQcPZjF0ZaKiRAREVmcWw+y8epPsVAJYGSIH8axONpiMREiIiKLkp1XgEnfn9AUR384lMXRloyJEBERWYySxdH167A4mpgIERGRBVlbsjh6HIujiYkQERFZiCPx97CkRHF0RxZHE5gIERGRBbj1IBszNquLo58P8WVxNGkwESIiIrNWsji6jZ8rPhzSksXRpMFEiIiIzIIiNx/J8hytNiEE5j4qjq7nZIc14zrA3pbF0VSMiRAREZk8RW4+xq8/jpFropGUXpwMrT18Hf93JhkA4OZkB0cpkyDSxkSIiIhMXpayAA8y83DrYTZGrVUnQyWLowFAWaBClrLAgFGSMZIIIYShg6hNCoUCMpkMcrkcLi4uhg6HiIhqSVJ6Dkatjcath9nwltlDkVuAzEeJj19dB2ydHAZvVwcDR0mVpa/Pb/YIERGRWfB2dcCWSaHwktkjSZ7LJIgqhYkQERGZjRRFcQJUZMWotkyCqEJMhIiIyCzsv3AXo9dGIyNXOxGavfW0VgE1UUlMhIiIyOT9EH0Tk78/AWWBCgDgW9cB26aGwd/NUauAmqg0JkJERGSyhBD4dO9lvLfzHFSPbv3xq+uAnyeHoUNDN2yZFKqVDJWeZ4iIiRAREZmk/EIV3vjlDFYevAoA8JTZlymMLiqg9ndzRL06dnCS2hgyZDJCfEUQEZHJyVQWYNqPp3D4yj1YW0mwcFhL9G/lhSxlAbxk2oXR3q4O2Do5FE5SG7jY2xooYjJWTISIiMikpGbkYsLGGJxLVMDB1hpfj22PXs3cAaDCRKd0ckRUhIkQERGZjGv3MjF+/XHcSctBPSc7rA/viDZ+roYOi0wYEyEiIjIJJ2+m4eVNMUjLzkdAPUdsmtAJDes5GTosMnFMhIiIyOjtO5+CV3+KhbJAhTZ+rlg/PgT16kgNHRaZASZCRERk1L6Pvon5v6lvj3+mmTtWjmkHRzt+fFHN4CuJiIiMkhACn+67jK8OXgMAjO7kh4+GtISNNWd+oZrDRIiIiIxOfqEKb207g+2nEgEAs/s8hZm9m0AikRg4MjI3TISIiMioZCoLMPWHkzgSfx/WVhIsGtYKz3f0M3RYZKaYCBERkdFIzcjFSxticD7p0RxB49qjV1N3Q4dFZoyJEBERGYWScwTVr6OeI6i1r6uhwyIzx0SIiIgM7uTNNEzcFIN0zhFEesZEiIiIDGrv+RTM5BxBZCBMhIiIyGBKzhHUu5k7vuQcQaRnfLUREZHeCSHw372X8fUhzhFEhsVEiIiI9CqvQIW3txfPETSn71N49RnOEUSGYdDUe9WqVWjdujVcXFzg4uKCsLAw/PHHH5rlQghERETA29sbDg4O6NmzJ86fP2/AiImIqDoylQWYuCkG208lwtpKgqUjWmNm7yAmQWQwBk2EfH19sXjxYpw4cQInTpzAM888gyFDhmiSnaVLl2LZsmVYuXIlYmJi4Onpib59+yIjI8OQYRMRkQ5SFbkYueYYjsTfh6OdNb4ZH8KJEsngJEIIYeggSnJzc8N///tfTJgwAd7e3pg1axbeeustAIBSqYSHhweWLFmCyZMnl7u9UqmEUqnUPFYoFPDz84NcLoeLi4tejoGIiLRdTVXPEZSYzjmCqHIUCgVkMlmtf34bTVVaYWEhtmzZgqysLISFhSEhIQEpKSno16+fZh2pVIoePXogKiqqwv0sWrQIMplM8+Pnx782iIgM6eTNh/j36igkpucgsL4Ttk/tyiSIjIbBE6GzZ8+iTp06kEqlmDJlCnbs2IHg4GCkpKQAADw8PLTW9/Dw0Cwrz7x58yCXyzU/t2/frtX4iYioYnvPp2DMun+Qnp2Ptn6u+HVKGPzrORo6LCINg9811rRpU8TFxSE9PR3btm3D+PHjERkZqVleuoBOCPHYojqpVAqplBNxEREZ2vfHbmD+rvNQCaBPc3d8Obo9HOysDR0WkRaDJ0J2dnZo0qQJACAkJAQxMTH4/PPPNXVBKSkp8PLy0qyfmppappeIiIiMR9k5gvzx0ZAWnCOIjJLRvSqFEFAqlQgMDISnpyf279+vWZaXl4fIyEh06dLFgBESEVFF8gpUeP3n05ok6PW+T2HhME6USMbLoD1C77zzDvr37w8/Pz9kZGRgy5YtOHToEPbs2QOJRIJZs2Zh4cKFCAoKQlBQEBYuXAhHR0eMGTPGkGETEVE5MnLzMfWHUzh69T6srSRYNLwVng/hDStk3HRKhPbs2YM6deqgW7duAICvvvoK69atQ3BwML766ivUrVu3Uvu5e/cuXnjhBSQnJ0Mmk6F169bYs2cP+vbtCwCYO3cucnJyMG3aNKSlpaFz587Yt28fnJ2ddQmbiIhqSaoiF+EbYnAhWQFHO2t8NbY9ejV1N3RYRE+k0zxCrVq1wpIlSzBgwACcPXsWHTt2xJw5c/DXX3+hefPm2LBhQ23EqhN9zUNARGSpOEcQ1QZ9fX7r1COUkJCA4OBgAMC2bdvw3HPPYeHChTh16hQGDBhQowESEZHxOnHjIV7+7gTSs/MRWN8Jm17qxNvjyaToVL1mZ2eH7OxsAMCff/6pmfTQzc0NCoWi5qIjIiKjtedcCsZ+UzxH0LapXZgEkcnRqUeoW7dumDNnDrp27Yrjx49j69atAIArV67A19e3RgMkIiLj892jOYKEAPo098CXo9txjiAySTr1CK1cuRI2Njb49ddfsWrVKvj4+AAA/vjjD/zrX/+q0QCJiMh4CCGwZM8lfPCbOgka09kfq8dxokQyXUb3pas1jcXSREQ1I69Ahbe3ncH22EQAwBv9nsL0Xk0eO9s/ka6Mrli6KrU/TDiIiMxL6TmCFg9vhf9wjiAyA5VOhFxdXSud9RcWFuocEBERGZe7j+YIuvhojqCvx7ZHT84RRGai0onQwYMHNb/fuHEDb7/9NsLDwxEWFgYAOHbsGDZt2oRFixbVfJRERGQQV1MzMH59jGaOoA3hndDKV2bosIhqjE41Qr1798bLL7+M0aNHa7Vv3rwZa9euxaFDh2oqvmpjjRARkW5O3HiIiZtOQJ7DOYJI//T1+a3TXWPHjh1DSEhImfaQkBAcP3682kEREZFhFc0RJM/hHEFk3nRKhPz8/LB69eoy7WvWrIGfH4vniIhM2aaoG5j640koC1To09wDP70SCjcnO0OHRVQrdJpQcfny5RgxYgT27t2L0NBQAEB0dDSuXbuGbdu21WiARESkHyqVwNK9l7E68hoA9RxBHw5uARtrnf5mJjIJOr26BwwYgPj4eAwePBgPHz7EgwcPMGTIEFy5coXfNUZEZILyClR4/ZfTmiTozWeb4pOhLZkEkdnjhIpERBau5BxBNlYSLOIcQWQEjG5CxdLS09Nx/PhxpKamQqVSaS178cUXqx0YERHVvtJzBK0a1wE9nmpg6LCI9EanROh///sfxo4di6ysLDg7O2tNtCiRSJgIERGZAO05gqTYEN6RcwSRxdFp8Pf111/HhAkTkJGRgfT0dKSlpWl+Hj58WNMxEhFRDYu58RAjVh1DYnoOGtV3wo5pXZgEkUXSqUcoMTERM2fOhKMj55QgIjI1e84lY+aWOOQVqNDO3xXfju/I2+PJYunUI/Tss8/ixIkTNR0LERHVMvUcQaeQ92iOoM0vc44gsmw69QgNHDgQb775Ji5cuIBWrVrB1tZWa/ngwYNrJDgiIqoZKpXAkr2XsCbyOgBgbGd/LOAcQUS63T5vZVXxG0cikRjVt8/z9nkisnR5BSrM/fU0dsYlAVDPETStZ2OtG12IjI1R3z5f+nZ5IiIyThm5+Zjyw0n8ffUBbKwkWDyiNf7dwdfQYREZDZ3nESIiIuNWco4gp0dzBHXnHEFEWnQeHI6MjMSgQYPQpEkTBAUFYfDgwThy5EhNxkZERDq6mpqB4V9H4WKyAvXrSLF1chiTIKJy6JQI/fDDD+jTpw8cHR0xc+ZMzJgxAw4ODujduzc2b95c0zESEVEVlDdHUEsfzhFEVB6diqWbN2+OSZMmYfbs2Vrty5Ytw7p163Dx4sUaC7C6WCxNRJbkj7PJeG2reo6g9v6u+IZzBJGJ0tfnt049QtevX8egQYPKtA8ePBgJCQnVDoqIiKpu498JmLZZPUdQ32AP/Mg5goieSKdEyM/PDwcOHCjTfuDAAfj58RuLiYj0SaUSWPT7RUT87wKEAMaF+mP1uA5wsLM2dGhERk+nu8Zef/11zJw5E3FxcejSpQskEgmOHj2KjRs34vPPP6/pGImIqAJ5BSq8+etp/MY5goh0olMiNHXqVHh6euKzzz7Dzz//DEBdN7R161YMGTKkRgMkIqLyKXLzMZVzBBFVi07F0qaExdJEZI7uKnIxfv1xXErJ4BxBZJaMembpmJgYqFQqdO7cWav9n3/+gbW1NUJCQmokOCIiKiv+bgbCN8QgMT0HDZyl2BDekbfHE+lIp2Lp6dOn4/bt22XaExMTMX369GoHRURE5Tue8BAjVkWp5whq4ITtUzlHEFF16NQjdOHCBbRv375Me7t27XDhwoVqB0VERGWVniPo2/EdUZe3xxNVi049QlKpFHfv3i3TnpycDBsbfn0ZEVFN21BijqB+wR7Y/EookyCiGqBTItS3b1/MmzcPcrlc05aeno533nkHffv2rbHgiIgsXdEcQQsezRH0QmhDrBrXAfa2nCOIqCbo1H3z2WefoXv37mjYsCHatWsHAIiLi4OHhwe+//77Gg2QiMhSKQsK8eYvZ7DrNOcIIqotOiVCPj4+OHPmDH788UecPn0aDg4OeOmllzB69GjY2trWdIxERBZHkZuPKd+fRNQ19RxBS0a0xgjOEURU43Qu6HFycsKkSZNqMhYiIgKQIs9F+AbOEUSkDzrVCAHA999/j27dusHb2xs3b94EACxfvhy//fZbjQVHRGRp4u9mYPjXf+NSSgYaOEuxdXIYkyCiWqRTIrRq1SrMmTMH/fv3R1paGgoLCwEAdevWxYoVK2oyPiIii1E0R1CSPJdzBBHpiU6J0Jdffol169bh3Xff1bpdPiQkBGfPnq2x4IiILMXvZ5Mx7tt/oMgtQHt/V2yb0gV+bo6GDovI7OlUI5SQkKC5W6wkqVSKrKysagdFRGRJ1h9NwEe71bfH9wv2wBej2/H2eCI90alHKDAwEHFxcWXa//jjDwQHB1c3JiIis6LIzUeyPKdMu0ol8M6Os/jw/zhHEJGh6NQj9Oabb2L69OnIzc2FEALHjx/HTz/9hEWLFuGbb76p6RiJiEyWIjcf49cfx4PMPGyZFApvVwcA6jmCpv94Cn9eTAUAvNY7CLP6BHGOICI90ykReumll1BQUIC5c+ciOzsbY8aMga+vLz7//HOMGjWqpmMkIjJZWcoCPMjMw62H2Ri1NhpbJoWijr0Nwtcfx6lb6QCAek52GNXJj0kQkQFIhBCiqhvl5ORACAFHR0fcv38f169fx99//43g4GA8++yztRGnzhQKBWQyGeRyOVxcXAwdDhFZoKT0HIxaG41bD7Ph7WoPGysr3HqYDQBo4CzFb9O7anqKiEhNX5/fOtUIDRkyBN999x0AwMbGBoMHD8ayZcswdOhQrFq1qkYDJCIydd6uDtj8Sme4OdkhKT1XkwR5utgzCSIyMJ0SoVOnTuHpp58GAPz666/w8PDAzZs38d133+GLL76o0QCJiEzduUQ5ZmyOxcOsPK32r8a2YxJEZGA6JULZ2dlwdnYGAOzbtw/Dhw+HlZUVQkNDNbNMExFZuvTsPLy38ywGrTyKuNvpKF0BNHvraSSll72bjIj0R6dEqEmTJti5cydu376NvXv3ol+/fgCA1NRU1uEQkcVTqQR+jrmNZz6LxA/RtyAE4GhnDQHA380R26aGwd/NUVNAzWSIyHB0SoQ++OADvPHGGwgICEDnzp0RFhYGQN07VN5Ei0REluJcohwjVkdh7rYzeJiVh8D6TnB3liI7rxD+bo7YMikUHRq6YcukUK1kqLx5hoio9ul01xgApKSkIDk5GW3atIGVlTqfOn78OFxcXNCsWbMaDbI6eNcYEemDPDsfn+67jB//uQmVAJzsrDG771MY1t4HL286UWYeIaD4brJ6deywaUInuNjbGvAIiIyLvj6/dU6ETAUTISKqTSqVwK+n7mDxH5c0xdBD2nrjnQHN4eFiD0A9qWKWsgBesrKF0cnyHDhJbZgEEZWir89vnSZUJCIi9TDYB7+d00yMGOReBx8OaYmwxvW01nOxt60w0SkvOSIi/WEiRERURfLsfHy2/zJ+iC4eBpvV5ymEdw2ArbVOpZdEZCBMhIiIKkmlEtj2aBjswaNhsEFtvPHugObwlNkbODoi0oVB/3RZtGgROnbsCGdnZ7i7u2Po0KG4fPmy1jpCCERERMDb2xsODg7o2bMnzp8/b6CIichSnU+S4z9rjuHNX8/gQVYemrjXweZXOuPL0e2YBBGZMIMmQpGRkZg+fTqio6Oxf/9+FBQUoF+/fsjKytKss3TpUixbtgwrV65ETEwMPD090bdvX2RkZBgwciKyFPKcfMz/7RwGfXkUJ2+mwdHOGu8MaIbfZz6NLo3rGzo8Iqomo7pr7N69e3B3d0dkZCS6d+8OIQS8vb0xa9YsvPXWWwAApVIJDw8PLFmyBJMnTy6zD6VSCaVSqXmsUCjg5+fHu8aIqEpUKoHtsYlY/MdF3M9UD4M919oL7w5szgJnIj0w6i9drS1yuRwA4ObmBgBISEhASkqKZuZqAJBKpejRoweioqLK3ceiRYsgk8k0P35+frUfOBGZlQtJCjy/5hje+OU07mfmoXEDJ/z4cmesHNOeSRCRmTGaYmkhBObMmYNu3bqhZcuWANSTNgKAh4eH1rpFX/Jannnz5mHOnDmax0U9QkRETyLPycfy/Vfw3bEbUD36WozXegfhpa6BsLMxqr8biaiGGE0iNGPGDJw5cwZHjx4ts0wi0f6qQiFEmbYiUqkUUqm0VmIkIvMkhMD2U4lYxGEwIotjFInQq6++il27duHw4cPw9fXVtHt6egJQ9wx5eXlp2lNTU8v0EhER6eJCkgLzd51DzI00AEDjBk74cEhLdG3CQmgiS2DQREgIgVdffRU7duzAoUOHEBgYqLU8MDAQnp6e2L9/v+bLXPPy8hAZGYklS5YYImQiMhOK3Hws26c9DDazdxAmcBiMyKIYNBGaPn06Nm/ejN9++w3Ozs6amiCZTAYHBwdIJBLMmjULCxcuRFBQEIKCgrBw4UI4OjpizJgxhgydiEyUEAI7YhOx8PdLuJ+pvsN0YCv1MFjJL0QlIstg0ERo1apVAICePXtqtW/YsAHh4eEAgLlz5yInJwfTpk1DWloaOnfujH379sHZ2VnP0RKRqbuYrMAHvxUPgzVq4IQPB7dEtyAOgxFZKqOaR6g28NvniUiRW3Q32E0UqgQcbNXDYBO7cRiMyFjx2+eJiKpJCIGdcYn4ZHfxMNiAVp54b2Awh8GICAATISIyU5dSFPhg53kcv/EQANCovhMWDGmBp4MaGDgyIjImTISIyKwocvOxYn88Nh27oRkGe7V3E0zsFgipjbWhwyMiI8NEiIjMghACv8Ul4ZPfL+JeRvEw2LsDg+HDYTAiqgATISIyeZdTMvD+b+dwPKF4GCxicAt0f4rDYET0eEyEiMhkZeTmY8Wf8dgYpR4Gs7e1wqvPBOHlpzkMRkSVw0SIiEyOEAK7Tifh493Fw2D/auGJ9wdxGIyIqoaJEBGZlCt3M/D+znP459EwWOCjYbAeHAYjIh0wESIik5CRm4/P/4zHBg6DEVENYiJEREataBjsk90XkfpoGOzZFh54/7lg+NZ1NHB0RGTqmAgRkdEqPQwWUM8REYNboGdTdwNHRkTmgokQERmdTGUBPv/zCjb8fQMFj4bBZvRqgle6N+IwGBHVKCZCRGQ0hBD435lkfLL7Au4qOAxGRLWPiRARGYX4uxn44LfzOHb9AQCg4aNhsF4cBiOiWsREiIgMKlNZgC8OxGP90QQUqASkNsXDYPa2HAYjotrFRIiIDEIIgf87k4yPSwyD9Q32wAfPBcPPjcNgRKQfTISISO/i72Zg/q7ziLpWYhhsUAv0asZhMCLSLyZCRKQ3mcoCfHkgHt+WGAab3qsJJnEYjIgMhIkQEdW6omGwT3ZfRIoiFwDQp7kH5g/iMBgRGRYTISKqVVdT1XeDFQ2D+bs5ImJwMJ5p5mHgyIiImAgRUS3JUhbgi7/i8e2R4mGwaT2bYHIPDoMRkfFgIkRENUoIgd/PpuCj/7vAYTAiMnpMhIioxlxNzUTErvM4evU+AMDPzQERg1qgd3MOgxGRcWIiRETVlqUswJd/XcW3R68jv1DAzsYK03o2xpQejTkMRkRGjYkQEemsaBjs490XkCxXD4P1buaO+YNawL8eh8GIyPgxESIinZQ3DDb/uRboE8xhMCIyHUyEiKhKsvPUw2DfHCkeBpvaozGm9uQwGBGZHiZCRFQpQgj8cS4FH//fBSQ9GgZ7ppk75g8KRsN6TgaOjohIN0yEiOiJrt1TD4MdiVcPg/nWVd8NxmEwIjJ1TISILJgiNx9ZygJ4yRzKLEuW50AiAb6Luol1JYbBpvRojGkcBiMiM8FEiMhCKXLzMX79cTzIzMOWSaHwdi1OhhLTsjHkq78hz8lHfqEAAPRq2gARg1twGIyIzIqVoQMgIsPIUhbgQWYebj3Mxqi10UhKzwEARF9/gN7LInE/Mw/5hQKeMnusezEE68M7MgkiIrMjEUIIQwdRmxQKBWQyGeRyOVxcXAwdDpFRSUrPwai10bj1MBu+dR3QKdAN208lapa/1CUAc//VDA52HAYjIv3S1+c3h8aILJi3qwO+HR+C/6w+hjtpObiTpk6C7G2t8N2ETugUWM/AERIR1S4mQkQW6kGmEpuibmDTsZuQ5+RrLfthYieEBDAJIiLzx0SIyMLcfpiNb45cx9YTt5GbrwIA2FhJUKAqHiWf8/OZMgXURETmiMXSRBbiUooCs7fGoeenh7Dp2E3k5qvQzNMZ9evYoUAl4O/miG1Tw+Dv5limgJqIyFwxESIyczE3HmLCxhj8a8UR7IhNRKFK4Omg+vhidFtkKQtwPzMP/m6O2DIpFB0aumHLpFCtZChZzmSIiMwXh8aIzJBKJXDwcipWHbqGEzfTAAASCTCgpRem9GiMVr4yKHLzseHvG5BIJFrDYN6uDtgyKRSj1kajXh07OEn5zwQRmS/ePk9kRvILVfjf6SSsjryGK3czAQB21lYY0cEXk7o3QmB97XmAnjSztJPUBi72tnqJnYioJN4+T0SVlp1XgK0xt/HNkQQkPqrrqSO1wdhQf0zsGgh3F/tyt3Oxt60w0SkvOSIiMjdMhIhMWFpWHr47dhMboxKQlq2+Bb5+HTtM6BaIsZ0bQubA3hwiosdhIkRkgpLSc/Dt0QT8dPwWsvMKAQD+bo6Y1L0R/t3Bl1+ISkRUSUyEiEzI1dQMrI68jp2xiZp5f4K9XDC1Z2P0b+kJG2veCEpEVBVMhIhMwKlbaVh16Br2X7iraQtt5IapPZuge1B9SCQSA0ZHRGS6mAgRGSkhBCKv3MOqQ9fwT8JDTfuzLTwwpUdjtPOva8DoiIjMAxMhIiNTUKjC7rPJWB15HReTFQAAW2sJhrb1weQejdDE3dnAERIRmQ8mQkRGIje/EL+cvIO1h6/h9kP1LfCOdtYY08kfE58O5O3sRES1gIkQkYHJc/LxQ/RNbPg7Afcz8wAAbk52CO8SgBfDGsLV0c7AERIRmS8mQkQGcleRi/VHE/DjP7eQqSwAAPi4OmBS90Z4PsQPDna8BZ6IqLYxESLSs+v3MrH28HVsP5WIvEIVAKCphzOm9myMga29YMtb4ImI9IaJEJGenLmTjtWR1/DHuRQUfcNfx4C6mNqzMXo1dect8EREBsBEiKgWCSHw99UHWBV5FX9ffaBp79PcHVN6NEZIgJsBoyMiIiZCRLWgUCWw93wKVh26hrOJcgCAtZUEQ9p4Y3KPxmjqyVvgiYiMARMhohqkLCjE9lOJWHv4OhLuZwEA7G2tMKqjP15+OhC+dR0NHCEREZXERIioBmTk5uPHf25h/dEEpGYoAQAyB1uM7xKA8C4BcHPiLfBERMbIoLenHD58GIMGDYK3tzckEgl27typtVwIgYiICHh7e8PBwQE9e/bE+fPnDRMsUTnuZSixdM8ldFn8Fxb/cQmpGUp4yezx/nPBiHr7Gczp+xSTICIiI2bQHqGsrCy0adMGL730EkaMGFFm+dKlS7Fs2TJs3LgRTz31FD7++GP07dsXly9fhrMzayzIcG4+yMLaw9fxy8k7yCtQ3wLfxL0OJndvhCFtfWBnw1vgiYhMgUETof79+6N///7lLhNCYMWKFXj33XcxfPhwAMCmTZvg4eGBzZs3Y/LkyfoMlQgAcD5JjtWR17H7TBJUj26Bb+vnimk9G6NPcw9YWfEWeCIiU2K0NUIJCQlISUlBv379NG1SqRQ9evRAVFRUhYmQUqmEUqnUPFYoFLUeK5k3IQSirz/EqshrOHzlnqa9Z9MGmNKjMToHunEOICIiE2W0iVBKSgoAwMPDQ6vdw8MDN2/erHC7RYsWYcGCBbUaG1kGlUpg/8W7WHXoGuJupwMArCTAc629MaVHYwR7uxg2QCIiqjajTYSKlP5LWwjx2L++582bhzlz5mgeKxQK+Pn51Vp8ZH7yClTYGZeINZHXcO2e+hZ4qY0Vng/xwytPN4J/Pd4CT0RkLow2EfL09ASg7hny8vLStKemppbpJSpJKpVCKpXWenxkfrKUBfjp+C18cyQBKYpcAICzvQ1eDGuI8C6BaODM1xURkbkx2kQoMDAQnp6e2L9/P9q1awcAyMvLQ2RkJJYsWWLg6MicPMhUYlPUDWw6dhPynHwAgLuzFC8/HYjRnfzhbG9r4AiJiKi2GDQRyszMxNWrVzWPExISEBcXBzc3N/j7+2PWrFlYuHAhgoKCEBQUhIULF8LR0RFjxowxYNRkLu6kZeObIwnYEnMLufnqW+AD6zthcvdGGNbeB1IbawNHSEREtc2gidCJEyfQq1cvzeOi2p7x48dj48aNmDt3LnJycjBt2jSkpaWhc+fO2LdvH+cQomq5lKLAmsjr2HU6CYWP7oFv5SPDtJ6N0a+FJ6x5CzwRkcWQCCGEoYOoTQqFAjKZDHK5HC4uvMvHksXceIjVh67hwKVUTVu3JvUxtWdjdGlcj7fAExEZEX19fhttjRBRTVCpBA5eTsWqQ9dw4mYaAEAiAQa09MKUHo3Ryldm4AiJiMiQmAiRWcovVOF/p5OwJvI6Lt/NAADYWVthRAcfTOreGIH1nQwcIRERGQMmQmRWcvIKsTXmFtYdSUBieg4AoI7UBmND/TGxayDcXewNHCERERkTJkJkFtKz87Ap6iY2HbuBh1l5AID6dewwoVsgxnZuCJkDb4EnIqKymAiRUVLk5iNLWQAvmUOZZcnyHDhJbeBib4uk9Bx8ezQBPx2/hey8QgCAv5sjJnVvhH938IW9LW+BJyKiijERIqOjyM3H+PXH8SAzD1smhcLbtTgZSkrPwai10XCSWuMpD2fsPpOMgke3wAd7uWBKz8YY0NITNtZWhgqfiIhMCBMhMjpZygI8yMzDrYfZGLU2WpMMJaXnYNjXf+OuQgkAuJisLoIObeSGqT2boHtQfd4CT0REVcJ5hMgoFfX83HqYDb+6DhjfpSH+u/cKlAUqzTrPtvDAlB6N0c6/rgEjJSKi2qCvz28mQmS0ziXKMe7bf5Cena/VPrCVF2b3DUITd84wTkRkrjihIlmknLxC7LuQgh2xiTgSf1/zFRhF1r7QHv1aeBkoOiIiMjdMhMjgVCqB6OsPsD02EXvOpSBTWaBZZmdthbzC4uGwj3dfQksfV60CaiIiIl0xESKDuXI3A9tPJeK3uEQky3M17b51HdCnuQf2nk9BsjwX/m6OWD6yDWZvPV2mgJqIiKg6WCNEepWakYtdcUnYEZuI80kKTbuzvQ2ea+2FYe184VPXHqPX/oNbD7Ph7+aodddYUQG1v5sjtk4OLXeeISIiMn2sESKzUVT3s/1UIo5eLa77sbGSoGdTdwxv74NnmrlrJj9U5OajXh07ANDq+fF2dcCWSaEYtTYa9erYwUnKly8REVUPe4SoVhQW1f2cSsSec8nIejTrMwC09XPF8PY+eK61N9yc7MrdvrIzSxMRkXlijxCZpMspGdgeewe/xSYhRVFc9+Pn5oBhbX0wtJ0PGjWo88T9uNjbVpjocDiMiIhqChMhqraiup/tpxJxIbm47sfF3gYDW3tjeHsfhDSsy1mfiYjI6DARIp1k5xVg3/m72B6biKPx91A03Y+t9aO6n3Y+6FWi7oeIiMgYMRGiSitUCRy79gDbY+9g77kUrbqfdv6uGN5OXfdTt4K6HyIiImPDRIie6FKKAjtOJeK3OO26H383Rwxt54Nh7XwQWN/JgBESERHphokQlStVkYvf4pKwPTYRF0vU/cgcbDGwtReGt/NBB9b9EBGRiWMiRBrZeQXYe14938/fV+9r1f30ejTfT69m7pDasO6HiIjMAxMhC1eoEoi6dh87TiViz/kUZJeo+2nv74ph7X3xXCsv1v0QEZFZYiJkoS4mK7AjVv09X3cVSk17w3qOGNpWXfcTwLofIiIyc0yELMhdRS5+i0vE9lOJuJSSoWmXOdjiudZeGN7eB+39WfdDRESWg4mQmctSqut+dsSWrft5ppk7hrXzRa9mDVj3Q0REFomJkBkqVAn8ffU+dsQmYs+5FOTkF9f9hDSsi2HtfTCwlRdcHVn3Q0RElo2JkBm5kKTAjtg7+C0uCakZxXU/AfUcMaydL4a184F/PUcDRkhERGRcmAiZuBS5uu5nR6x23Y+roy0GtfbGsPY+aOfnyrofIiKicjARMkFZygLsOfeo7ufafYhHdT921lbo3dwdQ9v5oFdTd9jZWBk2UCIiIiPHRMhEFBSq8Pe1B9hx6g72nr+rVffTMaAuhrXzxcBWXpA52howSiIiItPCRMiICSFwIfnR93ydTsK9EnU/gfWdMKydD4a2Zd0PERGRrpgIGaFkeQ5+i0vCjlOJuHy3uO6nrqMtBrXxxrB2PmjLuh8iIqJqYyJkJDI1dT93EHXtQXHdj40V+jRXz/fT46kGrPshIiKqQUyEDKigUIWjj+b72Xs+Bbn5Ks2yTgFuGNbeBwNaeUHmwLofIiKi2sBESM+EEDifVPQ9X0m4n1lc99OoqO6nnQ/83Fj3Q0REVNuYCOlJsjwHO2OTsCP2Dq7czdS0uznZYVBrLwxr74s2vjLW/RAREekRE6EqUuTmI0tZAC+ZQ5llyfIcOElt4GKvHsrKyM3XzPdz7Lp23U/f5h4Y1s4HPZo2gK01636IiIgMgYlQFShy8zF+/XE8yMzDlkmh8HYtToaS0nMwam003Jxs8Ur3Rth77i72XShV9xPohuHtfNCfdT9ERERGgYlQFWQpC/AgMw+3HmZj1NpoTTKUmJaNEauOIUWRiztpwPQfYzXbNGrghOHtfDCkLet+iIiIjI1EiKIBG/OkUCggk8kgl8vh4uJS7f0V9fzcepgNb5k9ejVzx88nbiO/sPg01nOy08z305p1P0RERFVW05/fFWEipIOk9Bz8a8VhKHILtNp7N3PHmM7+6P4U636IiIiqQ1+JEIfGdODt6oAJXQOw4sBVTdt3Ezqi+1PuBoyKiIiIqordFjpISs/Br6cStdre23keSek5BoqIiIiIdMFEqIqKaoTupOXA380R26aGwd/NUVNAzWSIiIjIdDARqoJkeXGhtL+bI7ZMCkWHhm7YMilUKxlKljMZIiIiMgVMhKrASWqDenXsNElQ0TxC3q4OmmSoXh07OElZekVERGQKeNdYVfdXhZmliYiISDe8a8xIudjbVpjolJccERERkfHi0BgRERFZLCZCVZUrB+SJ5S+TJ6qXE5kKvp6JyMIxEaqKXDnwwwhg4wBAfkd7mfyOuv2HEfzwqAn8gK59fD3rD1/P+sHzrB9mdp6ZCFWFMhPIugek3QA2Diz+8JDfUT9Ou6Fersw0ZJSmjx/Q+sHXs37w9awfPM/6YYbnmYlQVch8gPDdQN2A4g+PW/8Uf2jUDVAvl/kYNk5Txw9o/eDrWT/4etYPnmf9MMPzzNvndSG/A6ztBWSlFrfZOQGNnlH/XyIBICnxf+DRf8pZVvL/eMyykttXtJ/Kbl9ObI9dp7LPgWpuX+L5cx4Ch/8LZN0HnBoAHV8BYtY+elwf6P4m4OCmXl/zEhZV/12zfUW/4zHtNf18j1unJp+v1DEpM4ALO9X/LyJ1AVoOB6TO0FwfiRW0r5mV9vXTLMcTlle0fenXRGW3RxXjq85zPTq+qjwXJEDWXWDbK4AiEXDxAZ5dBOydV/z43+sBZ0/t64PS16oclV6/omW10K7TvqoTa4llmanAH28CGSnq8/n0m+p/RzJTgDqeQP+lQJ0GeCydPhJ12MaUnyfrPrD3HSDzLlDHA+g2B4j6Qv161vwB5avD82rjt8/XkFo7kf+bBZzcUHP7IyIiMmU1mAQBnEfIuMnvAFf2aLc51AVCJgL2LsV/aZf3l3/RX/3lLivv/3jCssdsX+VlQOWft6JlpXoeqhJTefvOlQPJp4uP1astYC9T/15eT1m5PU9P+F1r+/L2W9Hv1Xm+0s9dW89XiWNSZgJnf9Ye05e6AK3+re7hLHkthQra10xVwe/lXd/Sy1Xl/w48ft0nPdcT1y1ajirGVdXnUmm/p4QKEAVAgbL4PFtLtXu2il8cpa57ZdtRQXtN7b+221FBexX3U6AEFCXqV2S+gE0V5nkr/f6p3EY6bGLiz5OfDaQlFD8etrbGkiB9MolE6Ouvv8Z///tfJCcno0WLFlixYgWefvppwwRTNA6akazOfoetBXZMUo+Lnvu1RrNhi1d0rkvKTQdG/chzXFOKznGuvOzr+dpffD3XpJI1FEVcvHiOa1p5/25Y2QAv7uR5rknlnecdk0zy9Wz0xdJbt27FrFmz8O677yI2NhZPP/00+vfvj1u3buk/GHli2UJS/85lC04ruq2QKq/kh0bdAGDCvlLn+M7jt6cn4+tZf/h61g+eZ/0ws/Ns9InQsmXLMHHiRLz88sto3rw5VqxYAT8/P6xatUr/wUjrqAt3S4+DynyLPzycGqjXI93xA1o/+HrWD76e9YPnWT/M8Dwb9dBYXl4eTp48ibffflurvV+/foiKiip3G6VSCaWyeAxeoVDUXED2MmDcNnVNRelbimW+QPjv6g+NohoW0k3RBzRQ/gf0xoH8gK4JfD3rB1/P+sHzrB9meJ6N+q6xpKQk+Pj44O+//0aXLl007QsXLsSmTZtw+fLlMttERERgwYIFZdpru+qcaliuvPwPaED9lwY/oMmU8PWsHzzP+qGn86yvu8aMfmgMACSlKt6FEGXaisybNw9yuVzzc/v2bX2ESDXNXlbxRH4yH/5jRqaFr2f94HnWDzM7z0Y9NFa/fn1YW1sjJSVFqz01NRUeHh7lbiOVSiGVSvURHhEREZk4o+4RsrOzQ4cOHbB//36t9v3792sNlRERERHpwqh7hABgzpw5eOGFFxASEoKwsDCsXbsWt27dwpQpUwwdGhEREZk4o0+ERo4ciQcPHuDDDz9EcnIyWrZsid9//x0NGzY0dGhERERk4oz6rrGaoK+qcyIiIqo5vGuMiIiIqJYxESIiIiKLxUSIiIiILBYTISIiIrJYRn/XWHUV1YLX6HeOERERUa0q+tyu7Xu6zD4RysjIAAD4+fkZOBIiIiKqqgcPHkAmq72v7TD72+dVKhWSkpLg7Oxc4feTGRuFQgE/Pz/cvn3bLG/5N/fjA8z/GM39+ADzP0Yen+kz92OUy+Xw9/dHWloaXF1da+15zL5HyMrKCr6+voYOQycuLi5m+eIuYu7HB5j/MZr78QHmf4w8PtNn7sdoZVW75cwsliYiIiKLxUSIiIiILBYTISMklUoxf/58SKVSQ4dSK8z9+ADzP0ZzPz7A/I+Rx2f6zP0Y9XV8Zl8sTURERFQR9ggRERGRxWIiRERERBaLiRARERFZLCZCREREZLGYCOlZREQEJBKJ1o+np+djt4mMjESHDh1gb2+PRo0aYfXq1XqKVjcBAQFljlEikWD69Onlrn/o0KFy17906ZKeIy/f4cOHMWjQIHh7e0MikWDnzp1ay4UQiIiIgLe3NxwcHNCzZ0+cP3/+ifvdtm0bgoODIZVKERwcjB07dtTSETze444vPz8fb731Flq1agUnJyd4e3vjxRdfRFJS0mP3uXHjxnKvaW5ubi0fTfmedA3Dw8PLxBoaGvrE/ZrCNQRQ7rWQSCT473//W+E+jekaLlq0CB07doSzszPc3d0xdOhQXL58WWsdU34fPun4zOF9WJlraKj3IRMhA2jRogWSk5M1P2fPnq1w3YSEBAwYMABPP/00YmNj8c4772DmzJnYtm2bHiOumpiYGK3j279/PwDgP//5z2O3u3z5stZ2QUFB+gj3ibKystCmTRusXLmy3OVLly7FsmXLsHLlSsTExMDT0xN9+/bVfM9deY4dO4aRI0fihRdewOnTp/HCCy/g+eefxz///FNbh1Ghxx1fdnY2Tp06hffffx+nTp3C9u3bceXKFQwePPiJ+3VxcdG6nsnJybC3t6+NQ3iiJ11DAPjXv/6lFevvv//+2H2ayjUEUOY6rF+/HhKJBCNGjHjsfo3lGkZGRmL69OmIjo7G/v37UVBQgH79+iErK0uzjim/D590fObwPqzMNQQM9D4UpFfz588Xbdq0qfT6c+fOFc2aNdNqmzx5sggNDa3hyGrPa6+9Jho3bixUKlW5yw8ePCgAiLS0NP0GpgMAYseOHZrHKpVKeHp6isWLF2vacnNzhUwmE6tXr65wP88//7z417/+pdX27LPPilGjRtV4zFVR+vjKc/z4cQFA3Lx5s8J1NmzYIGQyWc0GV0PKO8bx48eLIUOGVGk/pnwNhwwZIp555pnHrmPM1zA1NVUAEJGRkUII83sflj6+8pj6+7C8YzTU+5A9QgYQHx8Pb29vBAYGYtSoUbh+/XqF6x47dgz9+vXTanv22Wdx4sQJ5Ofn13ao1ZaXl4cffvgBEyZMeOKX3rZr1w5eXl7o3bs3Dh48qKcIqychIQEpKSla10gqlaJHjx6IioqqcLuKruvjtjEWcrkcEonkiV+CmJmZiYYNG8LX1xfPPfccYmNj9ROgjg4dOgR3d3c89dRTeOWVV5CamvrY9U31Gt69exe7d+/GxIkTn7iusV5DuVwOAHBzcwNgfu/D0sdX0Tqm/D6s6BgN8T5kIqRnnTt3xnfffYe9e/di3bp1SElJQZcuXfDgwYNy109JSYGHh4dWm4eHBwoKCnD//n19hFwtO3fuRHp6OsLDwytcx8vLC2vXrsW2bduwfft2NG3aFL1798bhw4f1F6iOUlJSAKDca1S0rKLtqrqNMcjNzcXbb7+NMWPGPPZLHps1a4aNGzdi165d+Omnn2Bvb4+uXbsiPj5ej9FWXv/+/fHjjz/ir7/+wmeffYaYmBg888wzUCqVFW5jqtdw06ZNcHZ2xvDhwx+7nrFeQyEE5syZg27duqFly5YAzOt9WN7xlWbq78OKjtFQ70Oz//Z5Y9O/f3/N761atUJYWBgaN26MTZs2Yc6cOeVuU7onRTyaDPxJPSzG4Ntvv0X//v3h7e1d4TpNmzZF06ZNNY/DwsJw+/ZtfPrpp+jevbs+wqy28q7Rk66PLtsYUn5+PkaNGgWVSoWvv/76seuGhoZqFTl27doV7du3x5dffokvvviitkOtspEjR2p+b9myJUJCQtCwYUPs3r37sQmDqV1DAFi/fj3Gjh37xDoRY72GM2bMwJkzZ3D06NEyy8zhffi44wPM431Y0TEa6n3IHiEDc3JyQqtWrSrM0D09PctktqmpqbCxsUG9evX0EaLObt68iT///BMvv/xylbcNDQ01+F8tlVF0x19516j0Xymlt6vqNoaUn5+P559/HgkJCdi/f/9j/wotj5WVFTp27GgS1xRQ91I2bNjwsfGa2jUEgCNHjuDy5cs6vSeN4Rq++uqr2LVrFw4ePAhfX19Nu7m8Dys6viLm8D580jGWpK/3IRMhA1Mqlbh48SK8vLzKXR4WFqa566rIvn37EBISAltbW32EqLMNGzbA3d0dAwcOrPK2sbGxFZ4TYxIYGAhPT0+ta5SXl4fIyEh06dKlwu0quq6P28ZQiv7xjY+Px59//qlTAi6EQFxcnElcUwB48OABbt++/dh4TekaFvn222/RoUMHtGnTpsrbGvIaCiEwY8YMbN++HX/99RcCAwO1lpv6+/BJxweY/vuwMsdYmt7eh1Uqz6Zqe/3118WhQ4fE9evXRXR0tHjuueeEs7OzuHHjhhBCiLffflu88MILmvWvX78uHB0dxezZs8WFCxfEt99+K2xtbcWvv/5qqEOolMLCQuHv7y/eeuutMstKH+Py5cvFjh07xJUrV8S5c+fE22+/LQCIbdu26TPkCmVkZIjY2FgRGxsrAIhly5aJ2NhYzd0aixcvFjKZTGzfvl2cPXtWjB49Wnh5eQmFQqHZxwsvvCDefvttzeO///5bWFtbi8WLF4uLFy+KxYsXCxsbGxEdHW1Ux5efny8GDx4sfH19RVxcnEhOTtb8KJXKCo8vIiJC7NmzR1y7dk3ExsaKl156SdjY2Ih//vlH78cnxOOPMSMjQ7z++usiKipKJCQkiIMHD4qwsDDh4+NjFtewiFwuF46OjmLVqlXl7sOYr+HUqVOFTCYThw4d0noNZmdna9Yx5ffhk47PHN6HTzpGQ74PmQjp2ciRI4WXl5ewtbUV3t7eYvjw4eL8+fOa5ePHjxc9evTQ2ubQoUOiXbt2ws7OTgQEBFT4D5kx2bt3rwAgLl++XGZZ6WNcsmSJaNy4sbC3txd169YV3bp1E7t379ZjtI9XdHt/6Z/x48cLIdS37s6fP194enoKqVQqunfvLs6ePau1jx49emjWL/LLL7+Ipk2bCltbW9GsWTODJX6PO76EhIRylwEQBw8e1Oyj9PHNmjVL+Pv7Czs7O9GgQQPRr18/ERUVpf+De+Rxx5idnS369esnGjRoIGxtbYW/v78YP368uHXrltY+TPUaFlmzZo1wcHAQ6enp5e7DmK9hRa/BDRs2aNYx5ffhk47PHN6HTzpGQ74PJY8CJCIiIrI4rBEiIiIii8VEiIiIiCwWEyEiIiKyWEyEiIiIyGIxESIiIiKLxUSIiIiILBYTISIiIrJYTISIiIjIYjERIqIquXHjBiQSCeLi4mr1eTZu3AhXV9ca2Vd4eDiGDh1aI/t6koCAAKxYsUIvz0VE1cdEiIiqxM/PD8nJyWjZsmWtPs/IkSNx5coVzeOIiAi0bdu2Vp+zKmoyUSMiw7ExdABEZFqsra3h6elZ4XIhBAoLC2FjU71/XhwcHODg4FCtfRARPQl7hIioDJVKhSVLlqBJkyaQSqXw9/fHJ598AqDs0NihQ4cgkUiwd+9ehISEQCqV4siRI4/dR9E26enpmueMi4uDRCLBjRs3AGj3uGzcuBELFizA6dOnIZFIIJFIsHHjxnJjLywsxJw5c+Dq6op69eph7ty5KP2VikIILF26FI0aNYKDgwPatGmDX3/9VbO8KL7du3ejTZs2sLe3R+fOnXH27FnN8pdeeglyuVwTT0REhGb77OxsTJgwAc7OzvD398fatWt1vBJEVNuYCBFRGfPmzcOSJUvw/vvv48KFC9i8eTM8PDweu83cuXOxaNEiXLx4Ea1bt9ZpHxUZOXIkXn/9dbRo0QLJyclITk7GyJEjy133s88+w/r16/Htt9/i6NGjePjwIXbs2KG1znvvvYcNGzZg1apVOH/+PGbPno1x48YhMjJSa70333wTn376KWJiYuDu7o7BgwcjPz8fXbp0wYoVK+Di4qKJ54033tCKISQkBLGxsZg2bRqmTp2KS5cu6XTsRFTLqvRd9URk9hQKhZBKpWLdunXlLk9ISBAARGxsrBBCiIMHDwoAYufOnZXeR9E2aWlpmrbY2FgBQCQkJAghhNiwYYOQyWSa5fPnzxdt2rR5YvxeXl5i8eLFmsf5+fnC19dXDBkyRAghRGZmprC3txdRUVFa202cOFGMHj1aK74tW7Zolj948EA4ODiIrVu3lhtfkYYNG4px48ZpHqtUKuHu7i5WrVr1xNiJSP9YI0REWi5evAilUonevXtXabuQkJBq76O65HI5kpOTERYWpmmzsbFBSEiIZnjswoULyM3NRd++fbW2zcvLQ7t27bTaSu7Hzc0NTZs2xcWLF58YR+vWrTW/SyQSeHp6IjU1VadjIqLaxUSIiLToWqDs5ORU6X1YWalH5UWJ2p38/HydnreqVCoVAGD37t3w8fHRWiaVSp+4vUQieeI6tra2ZbYpel4iMi6sESIiLUFBQXBwcMCBAwdqbR8NGjQAACQnJ2vanjQvkZ2dHQoLCx+7jkwmg5eXF6KjozVtBQUFOHnypOZxcHAwpFIpbt26hSZNmmj9+Pn5ae2v5H7S0tJw5coVNGvWrNLxEJHxY48QEWmxt7fHW2+9hblz58LOzg5du3bFvXv3cP78eUycOLFG9lGUdERERODjjz9GfHw8Pvvss8fuMyAgAAkJCYiLi4Ovry+cnZ3L7cF57bXXsHjxYgQFBaF58+ZYtmyZ1t1pzs7OeOONNzB79myoVCp069YNCoUCUVFRqFOnDsaPH69Z98MPP0S9evXg4eGBd999F/Xr19dMzBgQEIDMzEwcOHAAbdq0gaOjIxwdHSt1fojIeDARIqIy3n//fdjY2OCDDz5AUlISvLy8MGXKlBrbh62tLX766SdMnToVbdq0QceOHfHxxx/jP//5T4X7GzFiBLZv345evXohPT0dGzZsQHh4eJn1Xn/9dSQnJyM8PBxWVlaYMGEChg0bBrlcrlnno48+gru7OxYtWoTr16/D1dUV7du3xzvvvKO1r8WLF+O1115DfHw82rRpg127dsHOzg4A0KVLF0yZMgUjR47EgwcPMH/+fK1b6InINEiEKDXBBhGRhTt06BB69eqFtLQ0zh5NZOZYI0REREQWi4kQERERWSwOjREREZHFYo8QERERWSwmQkRERGSxmAgRERGRxWIiRERERBaLiRARERFZLCZCREREZLGYCBEREZHFYiJEREREFuv/AaltnqLhIdwNAAAAAElFTkSuQmCC", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "plt.plot(circuit_depths, unoptimised_times_list)\n", + "plt.scatter(circuit_depths, unoptimised_times_list, marker=\"x\", label=\"unoptimised\")\n", + "plt.plot(circuit_depths, optimised_times_list)\n", + "plt.scatter(circuit_depths, optimised_times_list, marker=\"x\", label=\"optimised\")\n", + "plt.xlabel(\"circuit depth\")\n", + "plt.ylabel(r\"seconds\")\n", + "plt.legend();" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Awesome!" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "qujax", + "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.11.4" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/examples/variational_inference.ipynb b/examples/variational_inference.ipynb index 88a5375..358af0a 100644 --- a/examples/variational_inference.ipynb +++ b/examples/variational_inference.ipynb @@ -74,11 +74,13 @@ } ], "source": [ - "target_prob_func = lambda x: 0.3 + x + 0.3 * jnp.sin(2*jnp.pi*x) + 0.3 * jnp.sin(4*jnp.pi*x)\n", - "target_pmf = target_prob_func(jnp.arange(2 ** n_qubits) / (2**n_qubits))\n", + "target_prob_func = (\n", + " lambda x: 0.3 + x + 0.3 * jnp.sin(2 * jnp.pi * x) + 0.3 * jnp.sin(4 * jnp.pi * x)\n", + ")\n", + "target_pmf = target_prob_func(jnp.arange(2**n_qubits) / (2**n_qubits))\n", "target_pmf /= target_pmf.sum()\n", "\n", - "plt.bar(jnp.arange(2 ** n_qubits), target_pmf);" + "plt.bar(jnp.arange(2**n_qubits), target_pmf);" ] }, { @@ -124,28 +126,26 @@ "outputs": [], "source": [ "def get_circuit(n_qubits, depth):\n", - " \n", " n_params = 2 * n_qubits * (depth + 1)\n", - " \n", - " \n", - " gates = ['H'] * n_qubits + ['Rx'] * n_qubits + ['Ry'] * n_qubits\n", + "\n", + " gates = [\"H\"] * n_qubits + [\"Rx\"] * n_qubits + [\"Ry\"] * n_qubits\n", " qubit_inds = [[i] for i in range(n_qubits)] * 3\n", " param_inds = [[]] * n_qubits + [[i] for i in range(n_qubits * 2)]\n", - " \n", + "\n", " k = 2 * n_qubits\n", "\n", " for _ in range(depth):\n", " for i in range(0, n_qubits - 1):\n", - " gates.append('CZ')\n", + " gates.append(\"CZ\")\n", " qubit_inds.append([i, i + 1])\n", " param_inds.append([])\n", " for i in range(n_qubits):\n", - " gates.append('Rx')\n", + " gates.append(\"Rx\")\n", " qubit_inds.append([i])\n", " param_inds.append([k])\n", " k += 1\n", " for i in range(n_qubits):\n", - " gates.append('Ry')\n", + " gates.append(\"Ry\")\n", " qubit_inds.append([i])\n", " param_inds.append([k])\n", " k += 1\n", @@ -270,7 +270,9 @@ "metadata": {}, "outputs": [], "source": [ - "init_param = random.uniform(random_key, shape=(n_params,), minval=-init_rad, maxval=init_rad)" + "init_param = random.uniform(\n", + " random_key, shape=(n_params,), minval=-init_rad, maxval=init_rad\n", + ")" ] }, { @@ -307,9 +309,9 @@ "metadata": {}, "outputs": [], "source": [ - "_, (params_path, cost_path) = scan(gd_iteration,\n", - " (init_param, init_opt_state),\n", - " jnp.arange(n_steps))" + "_, (params_path, cost_path) = scan(\n", + " gd_iteration, (init_param, init_opt_state), jnp.arange(n_steps)\n", + ")" ] }, { @@ -333,8 +335,8 @@ ], "source": [ "plt.plot(cost_path)\n", - "plt.xlabel('Iteration')\n", - "plt.ylabel('KL');" + "plt.xlabel(\"Iteration\")\n", + "plt.ylabel(\"KL\");" ] }, { @@ -360,8 +362,8 @@ "final_params = params_path[-1]\n", "final_st = param_to_st(final_params)\n", "\n", - "plt.plot(target_pmf, label='Target PMF')\n", - "plt.plot(jnp.square(jnp.abs(final_st.flatten())), label='Trained Circuit')\n", + "plt.plot(target_pmf, label=\"Target PMF\")\n", + "plt.plot(jnp.square(jnp.abs(final_st.flatten())), label=\"Trained Circuit\")\n", "plt.legend();" ] }, diff --git a/qujax/__init__.py b/qujax/__init__.py index c273b63..1b565b4 100644 --- a/qujax/__init__.py +++ b/qujax/__init__.py @@ -7,6 +7,7 @@ from qujax import gates +from qujax.statetensor import all_zeros_statetensor from qujax.statetensor import apply_gate from qujax.statetensor import get_params_to_statetensor_func from qujax.statetensor import get_params_to_unitarytensor_func @@ -15,6 +16,7 @@ from qujax.statetensor_observable import get_statetensor_to_expectation_func from qujax.statetensor_observable import get_statetensor_to_sampled_expectation_func +from qujax.densitytensor import all_zeros_densitytensor from qujax.densitytensor import _kraus_single from qujax.densitytensor import kraus from qujax.densitytensor import get_params_to_densitytensor_func @@ -33,6 +35,7 @@ from qujax.utils import print_circuit from qujax.utils import integers_to_bitstrings from qujax.utils import bitstrings_to_integers +from qujax.utils import repeat_circuit from qujax.utils import sample_integers from qujax.utils import sample_bitstrings from qujax.utils import statetensor_to_densitytensor diff --git a/qujax/densitytensor.py b/qujax/densitytensor.py index b355ae9..9f0e6ee 100644 --- a/qujax/densitytensor.py +++ b/qujax/densitytensor.py @@ -2,9 +2,11 @@ from typing import Callable, Iterable, Sequence, Tuple, Union +import jax from jax import numpy as jnp from jax.lax import scan - +from jax._src.dtypes import canonicalize_dtype +from jax._src.typing import DTypeLike from qujax.statetensor import ( UnionCallableOptionalArray, _arrayify_inds, @@ -128,6 +130,15 @@ def partial_trace( return densitytensor +def all_zeros_densitytensor(n_qubits: int, dtype: DTypeLike = complex) -> jax.Array: + """ + Returns a densitytensor representation of the all-zeros state |00...0> on `n_qubits` qubits + """ + densitytensor = jnp.zeros((2,) * 2 * n_qubits, canonicalize_dtype(dtype)) + densitytensor = densitytensor.at[(0,) * 2 * n_qubits].set(1.0) + return densitytensor + + def get_params_to_densitytensor_func( kraus_ops_seq: Sequence[KrausOp], qubit_inds_seq: Sequence[Sequence[int]], @@ -194,8 +205,7 @@ def params_to_densitytensor_func( """ if densitytensor_in is None: - densitytensor = jnp.zeros((2,) * 2 * n_qubits) - densitytensor = densitytensor.at[(0,) * 2 * n_qubits].set(1.0) + densitytensor = all_zeros_densitytensor(n_qubits) else: densitytensor = densitytensor_in params = jnp.atleast_1d(params) diff --git a/qujax/statetensor.py b/qujax/statetensor.py index a72c2ed..8e83b3f 100644 --- a/qujax/statetensor.py +++ b/qujax/statetensor.py @@ -3,7 +3,10 @@ from functools import partial from typing import Callable, Sequence, Union +import jax from jax import numpy as jnp +from jax._src.dtypes import canonicalize_dtype +from jax._src.typing import DTypeLike from qujax import gates from qujax.utils import Gate, UnionCallableOptionalArray, _arrayify_inds, check_circuit @@ -92,6 +95,15 @@ def _gate_func_to_unitary( return gate_unitary +def all_zeros_statetensor(n_qubits: int, dtype: DTypeLike = complex) -> jax.Array: + """ + Returns a statetensor representation of the all-zeros state |00...0> on `n_qubits` qubits + """ + statetensor = jnp.zeros((2,) * n_qubits, dtype=canonicalize_dtype(dtype)) + statetensor = statetensor.at[(0,) * n_qubits].set(1.0) + return statetensor + + def get_params_to_statetensor_func( gate_seq: Sequence[Gate], qubit_inds_seq: Sequence[Sequence[int]], @@ -147,8 +159,7 @@ def params_to_statetensor_func( """ if statetensor_in is None: - statetensor = jnp.zeros((2,) * n_qubits) - statetensor = statetensor.at[(0,) * n_qubits].set(1.0) + statetensor = all_zeros_statetensor(n_qubits) else: statetensor = statetensor_in params = jnp.atleast_1d(params) diff --git a/qujax/utils.py b/qujax/utils.py index 459a88a..0078c47 100644 --- a/qujax/utils.py +++ b/qujax/utils.py @@ -5,6 +5,7 @@ from typing import Callable, Iterable, List, Optional, Protocol, Sequence, Tuple, Union from warnings import warn +import jax from jax import numpy as jnp from jax import random @@ -467,7 +468,6 @@ def sample_bitstrings( Returns: Array with sampled bitstrings, shape=(n_samps, statetensor.ndim). - """ return integers_to_bitstrings( sample_integers(random_key, statetensor, n_samps), statetensor.ndim @@ -491,3 +491,32 @@ def statetensor_to_densitytensor(statetensor: jnp.ndarray) -> jnp.ndarray: 2 for _ in range(2 * n_qubits) ) return dt + + +def repeat_circuit( + circuit: Callable[[jax.Array, jax.Array], jax.Array], nr_of_parameters: int +) -> Callable[[jax.Array, jax.Array], jax.Array]: + """ + Repeats circuit encoded by `circuit` an arbitrary number of times. + Avoids compilation overhead with increasing circuit depth. + + Args: + circuit: The function encoding the circuit. + nr_of_parameters: The number of parameters that `circuit` takes. + + Returns: + A function taking an arbitrary number of parameters and returning as + many applications of `circuit` as the number of parameters allows. + An exception is thrown if this function is supplied with a parameter array + of size not divisible by `nr_of_parameters`. + """ + + def repeated_circuit(params: jax.Array, statetensor_in: jax.Array) -> jax.Array: + def f(state, p): + return circuit(p, state), None + + reshaped_parameters = params.reshape(-1, nr_of_parameters) + result, _ = jax.lax.scan(f, statetensor_in, reshaped_parameters) + return result + + return repeated_circuit diff --git a/tests/test_expectations.py b/tests/test_expectations.py index b77cdd9..5128785 100644 --- a/tests/test_expectations.py +++ b/tests/test_expectations.py @@ -129,7 +129,7 @@ def big_hermitian_matrix(hermitian_str_seq, qubit_inds): big_h = jnp.kron(big_h, hermitian_arrs[k]) return big_h - sum_big_hs = jnp.zeros((2**n_qubits, 2**n_qubits), dtype="complex") + sum_big_hs = jnp.zeros((2**n_qubits, 2**n_qubits), dtype=complex) for i in range(len(hermitian_str_seq_seq)): sum_big_hs += coefs[i] * big_hermitian_matrix( hermitian_str_seq_seq[i], qubit_inds_seq[i]