Demonstrator-3-Local_cost_function.ipynb

{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<br>\n",
    "\n",
    "#  Local cost function\n",
    "\n",
    "<br>"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As we mentioned, there are two types of cost functions (i) global and (ii) local. We have seen the global cost function and because of it's shortcomings like barrean plateau as the size of the system under investigation increases. There comes local cost function to our rescue, as it gives us the best possible chance to tackle the ubiquitous barren plateau problem in the variational algorithms. \n",
    "\n",
    "<br>\n",
    "\n",
    "The theoretical foundation for the local cost function is given below and the way in which the expectation value take the form. To our advantage we go from the need to calculate three expectation values to just calculating only one expectation value."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We have seen the global cost function in the second demonstrator, but now we will look at the local cost function theory and implementation. \n",
    "\n",
    "To this end, we map from the <b> global </b> to the <b> local </b> cost function, using the following transformation (same minimum)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{align*}\n",
    "\\sum \\limits_{m,n} c_m^\\ast c_n \\langle 0 \\vert V^\\dagger (\\alpha) A_m^\\dagger U {\\vert 0 \\rangle  \\langle 0 \\vert} U^\\dagger A_n V (\\alpha) \\vert 0 \\rangle \\quad \\rightarrow \\quad  \\sum \\limits_{m,n} c_m^\\ast c_n \\langle 0 \\vert V^\\dagger (\\alpha) A_m^\\dagger U {P} U^\\dagger A_n V (\\alpha) \\vert 0 \\rangle\n",
    "\\end{align*}"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "with \n",
    "\n",
    "\\begin{align*}\n",
    "{P}=\\dfrac{1}{2} - \\dfrac{1}{2n} \\sum \\limits_{j=0}^{n-1} Z_j\n",
    "\\end{align*}"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In <a href=\"https://arxiv.org/pdf/1909.05820.pdf\">Bravo-Prieto et al.</a> this mapping relies on Eq.(D7). This local cost function is although easier to implement, has a worse convergence and success rate."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The big advatage of the mapping is, that after the transformation we only need to determine the following expectation values"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{align*}\n",
    "& (i) \\quad \\  \\langle 0 \\vert V^\\dagger (\\alpha) A_m^\\dagger  A_n V (\\alpha) \\vert 0 \\rangle = \\langle \\psi (\\alpha) \\vert A_m^\\dagger  A_n \\vert \\psi(\\alpha) \\rangle \\\\\n",
    "& (ii) \\quad \\langle 0 \\vert V^\\dagger (\\alpha) A_m^\\dagger U Z_j U^\\dagger A_n V (\\alpha) \\vert 0 \\rangle = \\langle \\psi (\\alpha) \\vert A_m^\\dagger U Z_j U^\\dagger A_n \\vert \\psi(\\alpha) \\rangle\n",
    "\\end{align*}"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since unitary matrices satisfy\n",
    "\n",
    "\\begin{align}\n",
    "U^\\dagger = U^{-1} \\quad \\Rightarrow U U^\\dagger = U U^{-1} = \\mathbb{I}\n",
    "\\end{align}"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "let us define $ Z_{j=-1} = \\mathbb{I}$. With this nomenclatur, we only need to evaluate one single expectation value for the <b> local </b> cost function, namely:\n",
    "\n",
    "\\begin{align*}\n",
    "\\quad \\langle 0 \\vert V^\\dagger (\\alpha) A_m^\\dagger U Z_j U^\\dagger A_n V (\\alpha) \\vert 0 \\rangle\n",
    "\\end{align*}\n",
    "\n",
    "for \n",
    "\n",
    "\\begin{align*}\n",
    "m,n \\in \\lbrace 0,1,2, \\dots, N \\rbrace  \\quad \\wedge \\quad j \\in \\lbrace-1, 0, 1, \\dots, n-1\\rbrace\n",
    "\\end{align*}"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Moreover, we <b> do not </b> need to control the variational gate $V(\\alpha)$, since"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{align*}\n",
    "\\quad \\langle 0 \\vert V^\\dagger (\\alpha) A_m^\\dagger U Z_j U^\\dagger A_n V (\\alpha) \\vert 0 \\rangle =  \\langle \\psi (\\alpha) \\vert A_m^\\dagger U Z_j U^\\dagger A_n  \\vert \\psi(\\alpha) \\rangle\n",
    "\\end{align*}"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Even better, the unitariness of $U$ also allows us to <b> omit </b> the fact that we need to control $U$ and $U^\\dagger$ it in the Hadamard test, as well (will be proven below). <br>\n",
    "\n",
    "In what follows, we in particular prove, that the following circuit enables us the evaluate the expectation value above:\n",
    "\n",
    "<br>\n",
    "\n",
    "Figure 4: Hadamard test of the VQLS algorithm using the local cost function. The figure is taken from the <a href=\"https://pennylane.ai/qml/demos/tutorial_vqls.html\">Pennylane documentation</a> and for internal use only.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from IPython.display import Image\n",
    "Image(r\"Images/Local-Hadamard.png\")"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Analogously to the calculation above, we apply the ${\\text{gates}}$ sucessively in the steps $(i)-(iii)$ and measure eventually.\n",
    "\n",
    "<br>\n",
    "\n",
    "Let us start with the Hadamard test for evaluating the real part of our expectation value of interest:"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{align*}\n",
    "\\vert 0 \\rangle \\otimes  \\vert  \\psi \\rangle & \\overset{(i)}{\\rightarrow}  {H} \\vert 0 \\rangle \\otimes  \\vert  \\psi \\rangle = \\dfrac{1}{\\sqrt{2}} \\left(  \\vert 0 \\rangle + \\vert 1 \\rangle \\right) \\otimes \\vert \\psi \\rangle  \\quad \\\\ & \\overset{(ii)}{\\rightarrow} \\dfrac{1}{\\sqrt{2}} \\left(  \\vert 0 \\rangle + \\vert 1 \\rangle \\right) \\otimes   {CA^\\dagger_{m,0,1} U C Z_{j,0,1} U^\\dagger CA_{n,0,1}} \\vert \\psi \\rangle  \\\\\n",
    "& =  \\dfrac{1}{\\sqrt{2}} \\left(  \\vert 0  \\rangle \\otimes {U U^\\dagger} \\vert  \\psi \\rangle + \\vert 1 \\rangle \\otimes A^\\dagger_{m} U Z_{j} U^\\dagger A_{n} \\vert  \\psi \\rangle \\right) \\\\\n",
    "& =  \\dfrac{1}{\\sqrt{2}} \\left(  \\vert 0  \\rangle \\otimes  \\vert  \\psi \\rangle + \\vert 1 \\rangle \\otimes A^\\dagger_{m} U Z_{j} U^\\dagger A_{n} \\vert  \\psi \\rangle \\right)   \n",
    "\\\\ & \\overset{(iii)}{\\rightarrow} \\dfrac{1}{\\sqrt{2}} \\left(  {H}  \\vert 0  \\rangle \\otimes  \\vert  \\psi \\rangle +  {H} \\vert 1 \\rangle \\otimes A^\\dagger_{m} U Z_{j} U^\\dagger A_{n} \\vert  \\psi \\rangle \\right) \\\\\n",
    "& =  \\dfrac{1}{2} \\left( \\left[   \\vert 0  \\rangle + \\vert 1  \\rangle \\right] \\otimes  \\vert  \\psi \\rangle + \\left[   \\vert 0  \\rangle - \\vert 1  \\rangle \\right] \\otimes A^\\dagger_{m} U Z_{j} U^\\dagger A_{n} \\vert  \\psi \\rangle \\right) \\\\\n",
    "&  =\\dfrac{1}{2} \\left( \\vert 0  \\rangle \\otimes \\left[ \\mathbb{I} + A^\\dagger_{m} U Z_{j} U^\\dagger A_{n} \\right]   \\vert  \\psi \\rangle + \\vert 1 \\rangle \\otimes \\left[  \\mathbb{I}  - A^\\dagger_{m} U Z_{j} U^\\dagger A_{n} \\right]  \\vert  \\psi \\rangle \\right) \n",
    "\\end{align*}"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "By taking the difference of the probablity of measuring $0$ and $1$, we obtain the real value of the expectation value. Similar to the calculation of global cost function expectation values.\n",
    "\n",
    "\\begin{align*}\n",
    "{P(0)-P(1)} & = \\Big \\vert \\Big \\vert \\, [\\mathbb{I} +  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}] \\vert \\psi \\rangle \\, \\Big \\vert \\Big \\vert^2 - \\Big \\vert \\Big \\vert \\, [\\mathbb{I} -  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}] \\vert \\psi \\rangle \\, \\Big \\vert \\Big \\vert^2 \\\\ & = \\dfrac{1}{4} \\langle \\psi \\vert [\\mathbb{I} +  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}]^\\dagger [\\mathbb{I} +  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}] \\vert \\psi \\rangle \\\\ & \\quad - \\dfrac{1}{4}\\langle \\psi \\vert [\\mathbb{I} -  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}]^\\dagger [\\mathbb{I} -  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}] \\vert \\psi \\rangle \\\\\n",
    "&= \\dots =  {\\mathrm{Re} \\langle \\psi \\vert  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}  \\vert \\psi \\rangle}\n",
    "\\end{align*}"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Analogously, we can evaluate the imaginary part by introducing an additional phase gate for preperation:"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{align*}\n",
    "\\vert 0 \\rangle \\otimes  \\vert  \\psi \\rangle & \\overset{(i)}{\\rightarrow} {P\\left(-\\frac{\\pi}{2}\\right) H}   \\vert 0 \\rangle \\otimes  \\vert  \\psi \\rangle = \\dfrac{1}{\\sqrt{2}} \\left(  \\vert 0 \\rangle - i \\vert 1 \\rangle \\right) \\otimes \\vert \\psi \\rangle  \\quad \\\\ & \\overset{(ii)}{\\rightarrow} \\dfrac{1}{\\sqrt{2}} \\left(  \\vert 0 \\rangle - i \\vert 1 \\rangle \\right) \\otimes   {CA^\\dagger_{m,0,1} U C Z_{j,0,1} U^\\dagger CA_{n,0,1}} \\vert \\psi \\rangle  \\\\\n",
    "& =  \\dfrac{1}{\\sqrt{2}} \\left(  \\vert 0  \\rangle \\otimes {U U^\\dagger} \\vert  \\psi \\rangle - i \\vert 1 \\rangle \\otimes A^\\dagger_{m} U Z_{j} U^\\dagger A_{n} \\vert  \\psi \\rangle \\right) \\\\\n",
    "& =  \\dfrac{1}{\\sqrt{2}} \\left(  \\vert 0  \\rangle \\otimes  \\vert  \\psi \\rangle -i \\vert 1 \\rangle \\otimes A^\\dagger_{m} U Z_{j} U^\\dagger A_{n} \\vert  \\psi \\rangle \\right)   \n",
    "\\\\ & \\overset{(iii)}{\\rightarrow} \\dfrac{1}{\\sqrt{2}} \\left(  {H}  \\vert 0  \\rangle \\otimes  \\vert  \\psi \\rangle -  {H} i \\vert 1 \\rangle \\otimes A^\\dagger_{m} U Z_{j} U^\\dagger A_{n} \\vert  \\psi \\rangle \\right) \\\\\n",
    "& =  \\dfrac{1}{2} \\left( \\left[   \\vert 0  \\rangle + \\vert 1  \\rangle \\right] \\otimes  \\vert  \\psi \\rangle -i \\left[   \\vert 0  \\rangle - \\vert 1  \\rangle \\right] \\otimes A^\\dagger_{m} U Z_{j} U^\\dagger A_{n} \\vert  \\psi \\rangle \\right) \\\\\n",
    "&  =\\dfrac{1}{2} \\left( \\vert 0  \\rangle \\otimes \\left[ \\mathbb{I} -i A^\\dagger_{m} U Z_{j} U^\\dagger A_{n} \\right]   \\vert  \\psi \\rangle + \\vert 1 \\rangle \\otimes \\left[  \\mathbb{I}  +i A^\\dagger_{m} U Z_{j} U^\\dagger A_{n} \\right]  \\vert  \\psi \\rangle \\right) \n",
    "\\end{align*}"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{align*}\n",
    "{P(0)-P(1)} & = \\Big \\vert \\Big \\vert \\, [\\mathbb{I} -i A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}] \\vert \\psi \\rangle \\, \\Big \\vert \\Big \\vert^2 - \\Big \\vert \\Big \\vert \\, [\\mathbb{I} + i  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}] \\vert \\psi \\rangle \\, \\Big \\vert \\Big \\vert^2 \\\\ & = \\dfrac{1}{4} \\langle \\psi \\vert [\\mathbb{I} - i  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}]^\\dagger [\\mathbb{I} -i  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}] \\vert \\psi \\rangle \\\\ & \\quad - \\dfrac{1}{4}\\langle \\psi \\vert [\\mathbb{I} + i  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}]^\\dagger [\\mathbb{I} + i  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}] \\vert \\psi \\rangle \\\\\n",
    "&= \\dots =  {\\mathrm{Im} \\langle \\psi \\vert  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}  \\vert \\psi \\rangle}\n",
    "\\end{align*}"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$\\dots$ and evaulte the overall expectation value by adding real and imaginary part\n",
    "\n",
    "<br>\n",
    "\n",
    "\\begin{align*}\n",
    "\\langle \\psi \\vert  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}  \\vert \\psi \\rangle = \\mathrm{Re} \\langle \\psi \\vert  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}  \\vert \\psi \\rangle + i \\, \\mathrm{Im} \\langle \\psi \\vert  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}  \\vert \\psi \\rangle\n",
    "\\end{align*}"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Notice how $ {U U^\\dagger = \\mathbb{I}}$ allows us to drop the control line of the Hadamard test."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<b> Last (optimization) step </b>:\n",
    "\n",
    "- Adjust the parateters $\\alpha$ classically and optimize the costfunction (<a href=\"https://link.springer.com/content/pdf/10.1007/s11128-021-03140-x.pdf\">ongoing research topic</a>)\n",
    "\n",
    "    - Use SPSA if their is noise \n",
    "    - Use COBYLA if you are noise-free (no shot-noise, as well)\n",
    "\n",
    "    <br>"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- Return the optimal result \n",
    "\\begin{align*}\n",
    "\\vert \\psi (\\alpha_\\mathrm{opt})\\rangle = A \\vert x (\\alpha_\\mathrm{opt}) \\rangle \\approx \\vert b \\rangle\n",
    "\\end{align*}\n",
    "\n"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the case of local cost function, we have all the functions (matrices and vectors involved in calculating the expectation values) are all defined in the main program. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Grab functions and modules from dependencies\n",
    "import math\n",
    "import numpy as np\n",
    "from scipy.optimize import minimize, fmin_cg, OptimizeResult\n",
    "from scipy.sparse import csr_matrix\n",
    "from qiskit.algorithms.optimizers import SPSA\n",
    "import mthree\n",
    "from scipy.sparse import random\n",
    "from scipy import stats\n",
    "from numpy.random import default_rng\n",
    "from numpy import linalg as LA\n",
    "\n",
    "#from collections.abc import Iterable\n",
    "import functools\n",
    "import inspect\n",
    "import itertools\n",
    "#import numbers\n",
    "from operator import matmul\n",
    "\n",
    "# Grab functions and modules from Qiskit needed\n",
    "import qiskit\n",
    "from qiskit import QuantumCircuit, QuantumRegister, ClassicalRegister, execute, transpile, Aer, assemble\n",
    "from qiskit.circuit import ParameterVector\n",
    "import qiskit.circuit.library.n_local as lib_local\n",
    "from qiskit.circuit import Gate\n",
    "\n",
    "#from qiskit.providers.aer import UnitarySimulator\n",
    "from qiskit.quantum_info import Statevector\n",
    "from qiskit.visualization import plot_bloch_multivector\n",
    "\n",
    "#from qiskit.opflow import X, Y, Z, Zero, One, I\n",
    "from qiskit.quantum_info.operators.symplectic.pauli_list import PauliList\n",
    "from qiskit.quantum_info.operators.symplectic.pauli import Pauli\n",
    "\n",
    "from qiskit.quantum_info import SparsePauliOp\n",
    "\n",
    "from qiskit.circuit.library import EfficientSU2, TwoLocal, NLocal, PauliTwoDesign\n",
    "#from qiskit_nature.circuit.library import UCCSD, PUCCD, SUCCD\n",
    "\n",
    "from qiskit.providers.ibmq.runtime import UserMessenger\n",
    "\n",
    "# Plotting\n",
    "import matplotlib.pyplot as plt\n",
    "from matplotlib.ticker import LinearLocator, FormatStrFormatter\n",
    "from qiskit import IBMQ\n",
    "\n",
    "from typing import List, Union"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the Local cost function, we only need to calculate one expectation value $\\langle \\psi \\vert  A^\\dagger_{m} U Z_{j} U^\\dagger A_{n}  \\vert \\psi \\rangle$ unlike in the global cost funcion. The logic of the program is as follows: \n",
    "\n",
    "- First, we determine the number of the qubits depending on the vector\n",
    "- We decompose the matrix into string og Pauli operatos along with the coefficients\n",
    "- We then take the ansatz input which can be custom made (Hardware efficient) or import ansatzes from the qiskit library\n",
    "- we then calculate the expectation value and the cost function. We do that by seeing the optimization routine, the user has requested, eg. SPSA verses SciPy ones, and dispatch to the correct optimizer. The selected optimizer starts at x0 and calls 'cost_loc' everytime the optimizer needs to evaluate the cost function. The result is returned as a SciPy OptimizerResult object.\n",
    "- Finally, we end up with the result, cost function and the quantum solution."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Once the matrix and vector are given as an input, we can start building the quantum circuit. The number of qubits required to build the circuit is log of the length of the vector."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "vector=np.ones(4)\n",
    "num_qubits = int(np.log2(len(vector))) # number of system qubits\n",
    "tot_qubits = num_qubits + 1  # Addition of an ancillary qubit\n",
    "ancilla_idx = num_qubits  # Index of the ancillary qubit for the Hadamard test(last position)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Decomposing the given matrix into pauli operators and coefficients"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Convert the matrix to Hamiltonian\n",
    "matrix = np.array([[ 1,0,  0.  ,  0. ],\n",
    "       [-0.81,  2.62, -0.81,  0. ,],\n",
    "       [ 0.  , -0.81,  2.62, -0.81],\n",
    "       [ 0.  ,  0.  , 0,1]])\n",
    "pauli_op = SparsePauliOp.from_operator(matrix) \n",
    "coeffs = pauli_op.coeffs\n",
    "obs_list = pauli_op.paulis"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since $A_m, U, A_m^{\\dagger}$ are the terms in the above expextation value. We can look at them one by one and build the quantum circuit for each of them.  \n",
    "\n",
    "Controlled $CA$ takes in the input of number of qubits, PauliList of pauli operators and idx to select the pauli operator string. \n",
    "\n",
    "We get the quantum circuit with the number of qubits needed and with the Pauli gates applied on each of the qubit and we then convert the quantum circuit to gate."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "def CA(qc, obs_list, idx):\n",
    "\n",
    "    operators = obs_list[idx]   \n",
    "\n",
    "    for j in reversed(range(len(operators))):\n",
    "        if operators[j] == Pauli('I'):\n",
    "            None\n",
    "        elif operators[j] == Pauli('X'):\n",
    "            qc.cx(ancilla_idx, j)\n",
    "        elif operators[j] == Pauli('Y'):\n",
    "            qc.cy(ancilla_idx, j)\n",
    "        elif operators[j] == Pauli('Z'):\n",
    "            qc.cz(ancilla_idx, j)\n",
    "    return qc"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "the transpiled circuit for the backend is:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAKoAAADWCAYAAABBlhk4AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAMOklEQVR4nO3df2yUhR3H8XevxRYoRaSTassPkQNKLQwKiDMWdEBYNNOIIAwwQiIdNs2mMWiG6JIq/lhjdJNpnFGIyHQgihOJqRNKWfmpjIIdFoVQwA6xVKGlhfbu9sfZQqFyd9je83zL55U0hoc77qu+eZ67Qp5vTCAQCCDich6nBxAJh0IVExSqmKBQxQSFKiYoVDFBoYoJClVMUKhigkIVExSqmKBQxQSFKiYoVDFBoYoJClVMUKhigkIVExSqmKBQxQSFKiYoVDFBoYoJClVMUKhigkIVE+KcHkDax6rtcLg6+q+b2gPuHNn2v65C7aAOV8NX3zg9RdvRpV9MUKhigkIVExSqmKBQxQSFKiYoVDFBoYoJClVMcHWofr+fgoICvF4vCQkJDBs2jKKiIgYNGsTcuXOdHq9DKXxlDu8suoWA3998LOD3syI/m3+99lsHJwtydahz5swhPz+fnJwc1q5dy9SpU5k+fTr79u0jKyvL6fE6lLH3/JkTVQf5bO1zzce2f/AsJ78/QvaM5y7wzOhwbajLly9n6dKlvP/++zz00EPcfPPNLFiwgBtuuIHGxsbmUI8cOcLEiRPp0qULw4YNY8eOHQ5PbtNlCYlMmreMLav+yNGKUo4e2Mm21U8wKXc5neK7OD2ee/9SylNPPcWkSZMYO3Zsi+MDBgygU6dOZGZmAjBv3jwGDx7M6tWreeONN7jrrrsoLy8nNjbWibFNSxlwPVm3PcxHf50BBBh1+6P0usYdVy5XnlEPHTrE7t27mTJlynk/V1FRQUZGBvHx8Zw4cYI1a9bw2GOP0blzZ+bOnYvP52Pz5s0hXyMmJqZDfxUVrb+o//ajbv8DsZ0S6BSfyMjb5kf8/KKi9RHNGS7XhgqQkpLS4nhdXR1FRUXNl/29e/fSs2dPkpOTmx+TmZlJWVlZ9IbtYDyeWHqmZdAz7TpiPO7Jwz2TnKUpvPLy8hbHn332WSorKxkxYgQAtbW1JCUltXhMUlISNTU1IV8jEAh06K+xY8e1zf+MCI0dOy6iOcPlyveo/fv3Z+jQoSxatIgrrriC1NRUVq5cyYcffgjQfEbt2rUrJ06caPHc48ePk5iYGPWZpX258ozq8XhYsWIFGRkZzJs3j9mzZ5OcnExubi5xcXEMHToUAK/Xy7fffktVVVXzc3fv3s2QIUOcGl3aiSvPqAADBw5k3bp1LY7NmjWL9PR0OnfuDEC3bt249dZbyc/P5+mnn2bZsmXExMQwZswYJ0buMCbmLHF6hPO4NtTWbN++/bwIX3rpJWbOnEmPHj3wer288847+tZUB2Qm1JqaGsrLy7n//vtbHO/VqxeFhYUOTSXRYibUxMREfD6f02OIQ1z5YUrkXApVTFCoYoJCFRMUqpigUMUEhSommPk+qkQmtcfFPa/pDoDXXhnd1w0lJhDJ37WSDu/3bwb/+fwMZ+c4ly79YoJCFRMUqpigUMUEhSomKFQxQaGKCQpVTFCoYoJCFRMUqpigUMUEhSomKFQxQaGKCQpVTFCoYoJCFRMUqpigUMUEhSomKFQxQaGKCa4O1erSXr8fyg7DB/+Bf+6A0oPg84d8mlyAq++UMmfOHN59910WLlxIVlYWJSUlTJ8+naNHj/Lggw86PV6rKqpgSTEcq215PKkz3HMjDOjlzFzWufaMGu7S3scff5whQ4bg8XhYuXKlozMf+R4WfwzVtef/3Ik6ePkTOFh1/s9JaK4NNdylvV6vlxdeeIHRo0c7MWYLhbvhdCO0do+kAMHL/9rSaE/VMbgy1HCX9gLMnDmTCRMmkJCQEO0xW6hvgB0HWo+0SQAo+xq+r4vWVB2Ha0OF0Et7f4q23ubcK7U/vjBvNzcwY6Tj26dDbXuO9uuF4spQw13a6yan646H/9iT37fjJB2TKz/1h7u096doj7ttvvgxfHXkxy//McBVl8Oxyr1EcDKJqqbbTrrtbqSuPKOGu7TXbcZnhH6POuE6XBupm7nyjArhLe0FaGhowOfz4ff7aWhooL6+nvj4+Ije/7SVwVfB1NGwYmvrwf56OAzvG/WxOgTXhtqa1pb23nfffSxduhSA4uJiAPbv30+/fv2iPR4Av/DCwBT4915Y99/gsZsGwY1eSOnuyEgdgisv/a1pWtp77gepJUuWEAgEWnw5FWmT5G5w+1ljTh6pSH8qM2dULe29tJk5o8qlTaGKCQpVTFCoYoJCFRMUqpigUMUEhSomKFQxQaGKCQpVTFCoYoJCFRMUqpigUMUEhSomKFQxQaGKCQpVTFCoYoJCFRMUqpigUMUEhSomKFQxQaGKCQpVTFCoYoJCFRMUqpigUMUEhSomuDpUq0t7Ac5eKuJ314KRsDS47J7Jrr7jtMWlvX4/lHwJRXvOHHtidfA+/jcNhLhY52a7kEPH4JOyMz9+5O3gYoxfZgRXDjktJuC2hUI/WL58OTNmzGD9+vUt9qFOnjyZVatWsW3bNjIzM8nJyaGwsJDa2lpGjBjB4sWLSU9Pd2Rmnx+WbgyuPY/h/M0oA1Ng7jj3xbqnEv62Pvib7OyZY4DYWMi5GbwOb8V27aU/nKW9jY2NDBgwgK1bt1JVVcX48eO5++67HZoYir8IRgqtr+8p/x8Ufh7VkUKqb4DXN5wfKfywaNgHr20ILiN2kitDDXdpb9euXXn00UdJTU0lNjaWvLw8du3aRX19fdRn9gdgwxehH7exHBpd9P5v+3449SMbsSF4vO50cCGxk1wbKkS+tLekpIR+/fqFtWm6rZfPdk/uzbHa0P9utafgqmuGOr6ct+lr0Yv/wO+78OnS72tk4Z+WamnvuS5maW91dTW5ubk8+eSTUZnxXDGe8N94RvLY9uaJCW8Wj8Mzu/LDlN/vZ/jw4VRWVlJQUNBiaW9FRQWbN2/m+uuvb358XV0dEydOJDs727FQfX54bFXwjHkhl8VB/p0Q3yk6c4VSuBvW7Az9uDtGwDhnPqMCLj2jRrK0t7GxkalTp+L1eh2LFCDWE1wjGcro/u6JFGDMteAJcQWO88Co/tGZ58e48oz6Y2bNmsXOnTspLS1tPnbvvfdy9OhRVq9eTVycs98WPtUAf/k4+D3Jc8UAP0uC302ErvFRH+2CNn0Jb285/1tqTT/+zQ3B32BOMhVqeno6Y8aM4fXXXwfgwIEDzR+eYmPPvIcqKyujT58+jsxY3xC8lG756sy3dOJiYdQ1cNvP3Rdpk10HYW0pfP3dmWO9r4BJQyEj1bGxmpkJtaamhu7du/P888+Tl5fn9Dgh1TdA5XfBP0pNuRy6XOb0RKEFAnDkONTUQ7cE6OWiRcNmQpVLmys/TImcS6GKCQpVTFCoYoJCFRMUqpigUMUEhSomKFQxQaGKCQpVTFCoYoJCFRMUqpigUMUEhSomKFQxQaGKCQpVTFCoYoJCFRMUqpigUMUEhSomKNRW5OXlkZaW5vi9rOQM3SmlFRs3bsTr9ZKamkpjo8P3BBdAoV5QXFycQnUJXfrFBIUqJihUMUGhigkKtRU5OTmkpaXh8/lIS0sjNzfX6ZEuefrUH6FAIBDRfiRpGzqjRqhoy07efK+QRp+L1u9dAsyHumvXLiZPnkxycjIJCQl4vV4WLFjQLq916tRpNmzZyemGRuJi3bPU7FJg+s8IP/30U7Kzs+nbty8FBQX06dOH/fv3U1JS0i6vV/LZ55ysP8X4G1tfcSntx/R71HHjxrFnzx7Ky8tJSkqK6LmPPPNKO00lkXj64blhPc7spf/kyZMUFxczbdq0iCMVe8xe+qurq/H7/aSmXty2rnB/J0PwvekzL/+d3ldfyewpv7qo15OfxmyoPXr0wOPxcPjw4Yt6/sVc+r/Yd1BvGdpYh7/0d+nShezsbN566y2OHz/u9DjSzkx/mDr7U//8+fPp27cvFRUVFBcX8+qrr7bJa6zbtIOPNmwjd9Yd9L76yjb5NSVyZi/9AFlZWWzatImFCxfywAMPUF9fT+/evZk2bVqb/PqnTp2meGspg/r3VqQOM31GbW+V31Sx7L1Cpt12i0J1mEINwe/34/GYfSvfYShUMUGnCjFBoYoJClVMUKhigkIVExSqmKBQxQSFKiYoVDFBoYoJClVMUKhigkIVExSqmKBQxQSFKiYoVDFBoYoJClVMUKhigkIVExSqmKBQxQSFKiYoVDFBoYoJClVM+D/1MY36FoOPXAAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 206.997x264.88 with 1 Axes>"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "backend = Aer.get_backend('qasm_simulator') # using qasm simulator as backend. \n",
    "qc = QuantumCircuit(num_qubits+1, 1)\n",
    "mycirc = CA(qc,obs_list,4)\n",
    "mycirc = transpile(mycirc, backend)\n",
    "mycirc.draw('mpl')"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Similyrly for controlled $CA_{dag}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "def CA_dag(qc, obs_list, idx):\n",
    "\n",
    "    operators = obs_list[idx]   \n",
    "\n",
    "    for j in range(len(operators)):\n",
    "        if operators[j] == Pauli('I'):\n",
    "            None\n",
    "        elif operators[j] == Pauli('X'):\n",
    "            qc.cx(ancilla_idx, j)\n",
    "        elif operators[j] == Pauli('Y'):\n",
    "            qc.cy(ancilla_idx, j)\n",
    "        elif operators[j] == Pauli('Z'):\n",
    "            qc.cz(ancilla_idx, j)\n",
    "    return qc   "
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "transpiled circuit takes the form for the given backend:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAKoAAADWCAYAAABBlhk4AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAMR0lEQVR4nO3df2zU9R3H8efdFVsKFoFO0LaA2ANK1zKoKMZYwCFh00QnijBgURKo2DSbxqEZoktQ/LHG4DaHcU4gIOr4obgpMZ2TUlLKD2X8kGBRGAXWIWIVWij0fuyPW1sKpb0Wet/vG16P5NL0y/d67zbPfr93V/L9eMLhcBgRl/M6PYBINBSqmKBQxQSFKiYoVDFBoYoJClVMUKhigkIVExSqmKBQxQSFKiYoVDFBoYoJClVMUKhigkIVExSqmKBQxQSFKiYoVDFBoYoJClVMUKhigkIVExSqmBDn9ABut2oLHKpy5rFTusM9Nzjz2G6jUFtxqAq++trpKUSnfjFBoYoJClVMUKhigkIVExSqmKBQxQSFKiYoVDHB1aGGQiEKCwvx+/0kJCQwZMgQiouLGThwIDNmzHB6vGYVvTaNlfNuIxwKNWwLh0Isn5vLx2885OBktrk61GnTpjF37lzy8vJYs2YNEyZMYNKkSezdu5ecnBynx2vWyF/8nuNHD/DZmpcatm35+4uc+P4wuZNfauGe0hLXhrps2TIWL17M+++/z2OPPcbo0aOZPXs2N998M4FAoCHUw4cPM3bsWBITExkyZAhbt251dO4rEroybuZSNq76LUcqtnNk/zY2r36GcfnL6BSf6Ohslrn2P6U899xzjBs3jpEjRzbZnp6eTqdOncjKygJg5syZDBo0iNWrV7NkyRLuvfdeysvL8fl8TowNQO/0m8i583E++tNkIMzwu56k13XuPANY4coj6sGDB9m5cyf33XffOf9WUVFBZmYm8fHxHD9+nA8++ICnnnqKzp07M2PGDILBIGVlZa0+hsfjiepWXLy2Xd/D8Lt+g69TAp3iu3LDnbPa9TWKi9dGPafVW7RcGypA7969m2w/efIkxcXFDaf9PXv20LNnT5KTkxv2ycrKYteuXbEb9jy8Xh89UzPpmfpDPF5X/phNceVPsD688vLyJttffPFFKisrGTZsGAA1NTUkJSU12ScpKYnq6upWHyMcDkd1Gzly1MX5ptph5MhRUc9p9RYtVz5H7d+/P9nZ2cybN48ePXqQkpLCihUr+PDDDwEajqhdunTh+PHjTe577NgxunbtGvOZpWO58ojq9XpZvnw5mZmZzJw5kwcffJDk5GTy8/OJi4sjOzsbAL/fzzfffMPRo0cb7rtz504GDx7s1OjSQVx5RAUYMGAAn3zySZNtU6dOJSMjg86dOwNw5ZVXcscddzB37lyef/55li5disfjYcSIEU6MfI6xeYucHuGS4dpQm7Nly5ZzIlywYAFTpkyhe/fu+P1+Vq5c6ehbU9IxzIRaXV1NeXk5Dz/8cJPtvXr1oqioyKGpJFbMhNq1a1eCwaDTY4hDXPliSuRsClVMUKhigkIVExSqmKBQxQSFKiaYeR/VKSnd23/f+qsAXn917B/7UuMJt+X/Wkmb/OrNyMf5k52d41KgU7+YoFDFBIUqJihUMUGhigkKVUxQqGKCQhUTFKqYoFDFBIUqJihUMUGhigkKVUxQqGKCQhUTFKqYoFDFBIUqJihUMUGhigkKVUxQqGKCq0O1uGivZVU18M9d8P5nkY/fnXB6okauvlLKtGnTePfdd5kzZw45OTmUlpYyadIkjhw5wqOPPur0eJeMQBBWboGyLyEMeIh8/Nu/4JZ0+NkN4HP4kObaI2q0i/Y+/fTTDB48GK/Xy4oVKxye2qa3ymDD/yOFMz6GYf0eWL7JqckauTbUaBft9fv9vPzyy9x4441OjGnewW/h03+3vE/ZV3D4+5iMc16uDDXaRXsBpkyZwu23305CQkKsx7wkbNwbOdVHs5+TXBsqtL5o74WI5crJTq/g3NJt4ZvvEgwGWvxZhYIBFvzlLa0ufbZoF+2VC3f6RBTndI8nuv06kCtf9Ue7aO+FiMXVNusvO+nmK3vuOgSvrW15H6/Xx6vzHsL/xkMxmanZGRx75BZEu2ivXLhB18A1V53/eaoHSOsB6b1iOFQzXHlEhegW7QWoq6sjGAwSCoWoq6ujtraW+Pj4Nj3/uZx5vZA3GhZ8DIePNW6vfy+191UwfRQ4/eN0bajNaW7R3unTp7N48WIASkpKANi3bx/9+vWL9XhmXZUIv/4pbKuAJaWRbYOuheHXQXYaxLlgDWRXnvqbU79o79kvpBYtWkQ4HG5yU6RtF+eDnOsaP88bDcP6uSNSMHRE1aK9lzczR1S5vClUMUGhigkKVUxQqGKCQhUTFKqYoFDFBIUqJihUMUGhigkKVUxQqGKCQhUTFKqYoFDFBIUqJihUMUGhigkKVUxQqGKCQhUTFKqYoFDFBIUqJihUMUGhigkKVUxQqGKCQhUTFKqYoFDFBFeHannR3jMXQgm5d1GU86pz2TWTXX3FaYuL9oZCUPolFO9u3PbMarh1INw6wD2XGj/bwW8jK0rXe+IdGNoXfpwZWTXFaZ6wSxdBWrZsGZMnT2bt2rVN1kMdP348q1atYvPmzWRlZZGXl0dRURE1NTUMGzaMV155hYyMDEdmDoZg8XrYfqBxVZEzDegNM0a5L9bdlfDntZFfsjNn9gA+X+R6/n6Hl+9x7ak/mkV7A4EA6enpbNq0iaNHjzJmzBjuv/9+hyaGki8ikcK5kQKU/xeKPo/pSK2qrYOF686NFCKfB4Pwxjo43fIqlB3OlaFGu2hvly5dePLJJ0lJScHn81FQUMCOHTuora2N+cyhMKz7ovX91pdDwEXP/7bsg1OB5n+xILL95GnYuj+WU53LtaFC2xftLS0tpV+/flGtNH2xF5/tlpzGtzWtf281p+Ca67IdX6y3/jbvj38lFMWivXN+t1iL9p6tPYv2VlVVkZ+fz7PPPhuTGc/m8Ub/xLMt+3Y0rye6WbwOz+zKF1OhUIihQ4dSWVlJYWFhk0V7KyoqKCsr46abbmrY/+TJk4wdO5bc3FzHQg2G4KlVkSNmS66Ig7n3QHyn2MzVmqKd8MG21ve7exiMcuY1KuDSI2pbFu0NBAJMmDABv9/vWKQAPi/c4m99vxv7uydSgBHXg7eVM3CcF4b3j8085+PKI+r5TJ06lW3btrF9+/aGbQ888ABHjhxh9erVxMU5+7bwqTr4wz8i70mezQP8IAl+ORa6xMd8tBZt+BLe2XjuW2r1n//85sgvmJNMhZqRkcGIESNYuHAhAPv372948eTzNT6H2rVrF3369HFkxtq6yKl041eNb+nE+SIL4N75I/dFWm/HAVizHf7zXeO2tB4wLhsyUxwbq4GZUKurq+nWrRvz58+noKDA6XFaVVsHld9F/pTa+ypIvMLpiVoXDkeWQq+uhSsToFc3pydqZCZUuby58sWUyNkUqpigUMUEhSomKFQxQaGKCQpVTFCoYoJCFRMUqpigUMUEhSomKFQxQaGKCQpVTFCoYoJCFRMUqpigUMUEhSomKFQxQaGKCQpVTFCoYoJCbUZBQQGpqamOX8tKGulKKc1Yv349fr+flJQUAgGHrwkugEJtUVxcnEJ1CZ36xQSFKiYoVDFBoYoJCrUZeXl5pKamEgwGSU1NJT8/3+mRLnt61d9G4XC4TesjycWhI2obFW/cxpvvFREIumj5vcuA+VB37NjB+PHjSU5OJiEhAb/fz+zZszvksU6dOs26jds4XRcgzueeRc0uB6b/Rvjpp5+Sm5tL3759KSwspE+fPuzbt4/S0tIOebzSzz7nRO0pxtzS/BKX0nFMP0cdNWoUu3fvpry8nKSkpDbd94kXXuugqaQtnn98RlT7mT31nzhxgpKSEiZOnNjmSMUes6f+qqoqQqEQKSntW60r2t9kiDw3feHVt0i79moevO8n7Xo8uTBmQ+3evTter5dDhw616/7tOfV/sfeAnjJcZJf8qT8xMZHc3Fzefvttjh075vQ40sFMv5g681X/rFmz6Nu3LxUVFZSUlPD6669flMf4ZMNWPlq3mfypd5N27dUX5WtK25k99QPk5OSwYcMG5syZwyOPPEJtbS1paWlMnDjxonz9U6dOU7JpOwP7pylSh5k+ona0yq+PsvS9IibeeZtCdZhCbUUoFMLrNftU/pKhUMUEHSrEBIUqJihUMUGhigkKVUxQqGKCQhUTFKqYoFDFBIUqJihUMUGhigkKVUxQqGKCQhUTFKqYoFDFBIUqJihUMUGhigkKVUxQqGKCQhUTFKqYoFDFBIUqJihUMeF/Dh+KBpEX/dYAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 206.997x264.88 with 1 Axes>"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "qc_dag = QuantumCircuit(num_qubits+1, 1)\n",
    "mycirc_dag = CA_dag(qc_dag,obs_list,4)\n",
    "mycirc_dag = transpile(mycirc_dag, backend)\n",
    "mycirc_dag.draw('mpl')"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the expectation value, we have the $V^{\\dagger}$ and $U^{\\dagger}$, therefore we can have a look at how to produce this circuits. So complex conjugate of the unitary operators are inverse of a quantum circuit for unitary operator. Below, we see an example of an ansatz and also it's inverse."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAKoAAAB7CAYAAADkFBsIAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAIsElEQVR4nO3dfUxV9x3H8TcPRkBAQeJDWVelXB+gQNDJtD5cbI116qJz6nxc1EQMIzq1zibFh21AtU6Tus7wxxoHnWOmolFTdR2rwnTOODedZW6Fqp3zIXZQLaCYwYX94bjxCtMLvXDOj/t5JfcPz304X+Mnv3PP8eR+Apqbm5sRsblAqwcQ8YaCKkZQUMUICqoYQUEVIyioYgQFVYygoIoRFFQxgoIqRlBQxQgKqhhBQRUjKKhiBAVVjKCgihEUVDGCgipGUFDFCAqqGEFBFSMoqGIEBVWMoKCKERRUMYKCKkYItnoAuztwDm7csWbfsVEw62sde+/Hx6H2M9/O442IfjD0Jd9/roL6FDfuwGUL/sG/rNrP4O51q6fwHR36xQgKqhhBQRUjKKhiBJ1M+bFX89P5+z//SFBQDwIDgxgQNZgFL2fjTJlj9WitKKh+buGkjSyctAGXq5FDp3/GlqIFxMemEhsTb/VoHnToFwCCgoL5xteX42pq5PLNC1aP04qCKgA0NP6H90/nA/CVmCEWT9Oagurnij7MY+bGPkx/PZRffLCBtXPeIe6ZZADyfjWfP5QfdL92c8FMzn9y3JI5bR3UpqYmtm/fjsPhICQkhJSUFMrKyhg6dCgZGRlWj9em4tx0zh7M9Xq71Ra8nM3BnLsU/7CKtGFTufBIEL83Yyfv/nYz9x/UcvKjA0SERZMa3wn/P+oFWwd12bJl5OTksGLFCo4dO8bcuXOZP38+V65cYeTIkVaP161EhEWxds47nP3HUU6XHwIgKrwfs8avYdehVRR9mEvG9O2WzWfboBYVFVFYWMjhw4dZt24dEydOJDs7mzFjxtDY2OgO6u3bt5k8eTJhYWGkpKRw/vx5iyc3V2RYNN8ev5bdv3mdpqYmAF4ZtYQbVZXMHLuKyLBoy2azbVC3bNnClClTcDqdHtvj4+Pp0aMHSUlJAGRmZjJs2DCqq6vJyspi9uzZuFwuK0buFr41/vt8XnOLkj+/6972TN94yy9X2TKo169fp7y8nDlzWl94vnbtGomJifTs2ZPa2lqOHDnCpk2bCA0NJSMjA5fLxZkzZ566j4CAAK8eZWWl7Z7/7KE88jP6eDxuVpxq9+eUlZV6PWdH5t6RWcrCSRs8tvUKieTAjz/nlVFL2j1vR2b2li0v+F+//vD+tAEDBnhsr6+vp6ysjKlTpwJQWVlJ3759iYmJcb8mKSmJS5cuMXbs2K4b+DFpM7JJm+kZgOLcdGuG6SZsuaK2BK+iosJj+7Zt27h16xYjRowA4N69e0RGRnq8JjIykrq6uqfuo7m52auH05num79UBzid6V7P2Zlzr59XwAuDx3XKzN6y5YoaFxdHcnIyb7zxBtHR0cTGxlJcXMzRo0cB3CdSvXr1ora21uO9NTU1hIeHd/nM0rlsuaIGBgayb98+EhMTyczMZOnSpcTExJCVlUVwcDDJyQ8vSDscDqqqqqiurna/t7y8nISEBKtGl05iyxUVYMiQIZw4ccJj2+LFixk+fDihoaEAREREMG3aNHJycti6dSt79uwhICCA0aNHWzEyALM3lLZru3jHtkFty7lz51qFMD8/n0WLFhEVFYXD4WD//v0EBQVZNKF0Flse+ttSV1dHRUWF+0SqRf/+/SkpKaG+vp6LFy+Smppq0YRm++BPBSzdNpS/Xi4D4L3Sn7B61zi2FC2k0dUAwIbd01m9y7uTKl8zJqjh4eG4XC5Wrlxp9Sjd1hznD0h53sndun9z4fIJ3so6xeCBye4bU3KXvW/ZbEYd+sV3qr64ydZfP1wt4wYmM/TZNPdzH//rLClx6QCMcEzi+Pkiy+/6N2ZFFd/q3SuGN5eX8FbWKe4/qMHV1OB+rq7+LmEhD69P9wrpTV29Rb/A8QitqH6q5n41Pz2QSV39XW7f+ZTYGAd9e8cCEB7ah6ovbgBw/0EN4SF9LJz0Ia2ofur4+SJeTJzJjsxSEgeNpV/Uc+7nhjw7iotXHp5U/aXydwx/zrrLfS20ovqp1PiXeHPvdzn9t4OtnosK70dS3ARW7xpHvz5fZdb41V0+3+MUVD8VH5vKz1/9yP3n318sZu+JrcTGOEh53sm8ia8xb+JrHu/ZsHs60ZEDu3pUQEGV/5mQPJsJybOf+BpdnrKx2Cgz9x3Rz3dz2GG/Ac3tuddKxCI66xcjKKhiBAVVjKCgihEUVDGCgipGUFDFCAqqGEFBFSMoqGIEBVWMoKCKERRUMYJu83sKtUu3j9qlLaJ2aXvQoV+MoKCKERRUMYK+o/oxlfaKMVTaK0ZRaa8YQaW9Ymsq7fUBE0t7TaPSXh8wsbTXtHbpFirt7SBvS3s3b95MQkICgYGBFBcXWzy12VTa2wHelvY6HA527txJWlpaWx8j7WTX0l5bXkdtKe1ds2ZNq+ceLe0FWLRoEQB5eXldOmN3sCOztNW2ltJeu7Hlivq00l5ffD9Vu3Tn8Kt26UdLe1uapKF1aa9ddbd26fXzCqwewZ5B9ba098vw9tc23y6x7n5UpzOd4tyO/Sroub3W3I/qdKbTnO/7XzK15aHf29Je8R+2XFHBu9JegIaGBlwuF01NTTQ0NPDgwQN69uzZru8/Yn+2DWpb2irtXb58OYWFhQCcPHkSgKtXrzJo0KCuHg9Qu3RnseWhvy3/r7S3oKCA5uZmj4dVITXZo6W9ja4GVr09hm9mh3Oj6hP3a1Ta6wWV9na+ltLeoMBgfrTkIOOTPFtS1IoiXe5Jpb0BAQFERfS3cLrWjFlRxbeeVNprR1pR/dSTSnvtSCuqn3pSaa8daUX1U08q7QXI+eVcyj89xY2qSr6Tvp4XX5jRtQM+RkH1U08r7d24+L1W71Fpr1hOpb2GU2mvPfar0l4xgs76xQgKqhhBQRUjKKhiBAVVjKCgihEUVDGCgipGUFDFCAqqGEFBFSMoqGIEBVWMoKCKERRUMYKCKkZQUMUICqoY4b9iaKWhMifL8gAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 206.852x144.48 with 1 Axes>"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def custom_ansatz(num_qubits):\n",
    "    \n",
    "    circ = QuantumCircuit(num_qubits) # building the quantum circuit\n",
    "    theta = ParameterVector('a', length = num_qubits) # initializing the parametere for the ansatz\n",
    "    \n",
    "    for idx in range(num_qubits):\n",
    "        circ.h(idx)\n",
    "\n",
    "    #A very minimal variational circuit.\n",
    "    for idx in range(num_qubits):\n",
    "        #print(idx, element)\n",
    "        circ.ry(theta[idx], idx)\n",
    "    \n",
    "    return circ, len(theta)\n",
    "\n",
    "ansatz_circuit, num_params = custom_ansatz(2)\n",
    "ansatz_circuit.draw('mpl')"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "After transpiling tha ansatz on the backend, we obtain the following circuit."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAyYAAACoCAYAAAAVZfoHAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAmRklEQVR4nO3de3wU9b3/8dfmAgmESyBKMMglJCAEEkkQDAIBBFRKj9afoKhY4RRooByh0mN/Pznto8cD3tBSL1A5VqVVahVrpQKtoBBQBAURCFITCCQEAkIuhJALye7+/hgJLLltkt2dSfb9fDz2oZn57mc/2fmwmc/OzHdsTqfTiYiIiIiIiIkCzE5AREREREREjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJhOjYmIiIiIiJguyOwEpHGqLsLpb6G0AAKC4Jq+0DHS7KzMY3dAei7k5EOADfp2g36Rxv+LAJwohAPH4WIVRHSAxN4QEmx2VuIPLhTAdxlQVQEhHaDbAGgTanZW5jlXBl8dg/Nl0K4tDOkFXcPMzkqsosoO+48bn9kBNoiNhNhuYNPfc79iczqdTrOTqIvD4eD555/nlVde4fjx4/Tv358XXniB2bNnk5KSwqpVq8xO0adO7IeMLWCvdF3eKQrifwht/ewD/ts8eHMHnC93XR4RBjNGQ1S4OXmJNZSUw+pPIfO06/LgQJh8I6TcYEpa4geqKuDgRjhz2HW5LQD63Ax9kv1rZ8vugHV7Yfu34Lhqj2NYNEwZZvy7FP918ASs+RwuVLguv7YjzBwNkZ3MyUt8z9Kncs2cOZMnnniCOXPmsHHjRqZOncq0adPIysoiKSnJ7PR86uRBOPRRzaYE4NxJ2PMX42iKvzh6BlZtMXY+r5Z/AV7aBGfP+z4vsYaLVbDyEzh8uua6Sju8vwc+zfB9XtL6OR3w9fs1m5JL67J2wNHPfZ+Xmd7fA2n/qtmUAHyRBW9+Btb9ilS8LfM0/CENSitqrjtzHl7cBIUXfJ+XmMOyjcmaNWtYvXo169atY9GiRYwdO5bHH3+c5ORkqqqqqhuT06dPM3HiRNq1a0dCQgJ79+41OXPPc9jhcFo9A5xQWggn032Wkuk27DP+yNX2t8zphPJK+Pgbn6clFrE32zgdoL59nfVfGw2MiCedPQpFufWPOboLKst8k4/Zzp5v+EuAfceN03HFP3241/i7Xdff8wsVsOWQz9MSk1i2MXnyySe5/fbbSUlJcVkeExNDcHAwgwcPBiA1NZUbbriB/Px85s2bxz333IPdbjcjZa8pOAYXSxsed2K/11OxhIIS4xuW+nY6ncCXR41vx8X/fH4YGjpTpqzSuD5JxJNOHqDB4nPa4dS/fJKO6b7IaniMDdh1xOupiAWdPgfZ+fX/PQejPhwOn6QkJrNkY5Kbm0t6ejpTpkypsS4nJ4e4uDjatm3L+fPnWb9+Pb/61a8IDQ1l9uzZ2O12du7c2eBr2Gy2FvOY8eBP3XrfTh8vMj1XXzwGJd7i1vtRZYfwa6JMz1cP3z/SM3Ib/EMHMOtnvzA9Vz1a1+PztL0N7mU5nA7+e/HTpufqi8fyFW/gsNd/aNLusLNm7QbTc9XD949hIye49fe8ogradQg3PV89mv5wl2UbE4DISNfppsrKykhLS6s+jSszM5OuXbsSERFRPWbw4MF8803rOofnQvm5Bsc4nU63xrUGF8vc+z2dTgeV5brQxB9VlBbhdOPrtYul/vFvRnynpKwQh6P+Q7UBtgC/+byuKDvX8JX+TqcxTvxOhZufwQ5HFZUVutDEH1hyuuBLjUZGRgaTJk2qXv7MM8+Ql5dHYmIiABcuXKBjx44uz+3YsSMlJSUNvoaFJyOroaoCtq0wrjWpi81mY9htvXAuaTm/V1M5nbD078a5y3X9tjYbDLgugIrSYp/mJtbw0QHY0MCpjQE22LFhFR1C/Gt2P/Gu3K/hX5sbHvfq2qWsCV/q9XzMlvUdvLCp/jEBgUE8/dg0Nr40zTdJiWU4HPDEB1BYz+nqNmBIryDs/jTDjx+z5BGT6Oho4uPjWbp0KX/84x/5+OOPSU1N5bXXXgOoPmLSvn17zp93/Ua8uLiYsLDWNW9uUFu4fkg9A2wQGAw9EnyWkqlsNpgwqIGzJZwwbqCvMhKrSY6B0OD6T/W/OQY6hPgsJfETkQOhTXvqLb5rY6Gdn0xn3uca41HX22HDmBI2LsqXWYlVBATA+LgGBtlgrP6e+w1LNiYBAQG8++67xMXFkZqayowZM4iIiGDevHkEBQURHx8PQGxsLGfPniU///J0Hunp6Qwc2PoquO9o6HbpvgtXfcIHBsONd0OoH83zPSwaJhllUOMPns0G05KNGzOJf+oQCj8dByFtXJdfOqNkUA+4279mHBcfCWoDiVOgTburVnxfe+HXw8A7fJ6WaWw2+Pcr7it19ed11w7w07EQaMm9EfGFEbEw4fvm5Or6CLDB9BHQO6LG06SVsvQNFq82ffp09u3bx/79l8/RuPvuu+nZsydPPfUUb775JkuXLiUzM5PAwNZ3tyanE4pOwIl9cOr7qfNiU6D7IP+9m3BeEXyWeXk6ylsHGh9yupuwAJRehC+zjPsogHGn6eQY3U1YvK/qovE5/a/vT2O6JhaiBkHXPsaNFv2N3QEHco1/jwdPGMseSIYbe+nmimI4UQifZcCO7+8BNHGQ8Xkd3t7cvMS3WtTH4+7du2vcWHHlypUcPHiQ8PBwXnjhBd57771W2ZSAsSMV3gMG/eDysl43+W9TAtC9M9xz0+WffzhETYlc1q6N6x3efzwS+kWqKRHvC2rjenptwp0Q0dc/mxIwjojc2BNmjbm87KZoNSVyWVQ4TB1++edJCWpK/JElL36vTUlJCRkZGcydO9dlebdu3di0qYEr60RERERExNJaTGMSFhbW6m6cKCIiIiIiBj89qCwiIiIiIlaixkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREynxkREREREREwXZHYCrd23n8D577wXf/fbno/Z4VroP87zcf2Nt7e91amOmqel14+2v4iINJYaEy87/x0U5XovvjdjS/N4e9tL66b6ERERf6NTuURERERExHRqTERERERExHRqTKRFczjg1LnLP585Dw6nefmI9Zwru/z/xwug0m5eLiIiIlI3XWMiLY7dAQdy4fNMOHoGLl6xo7lkHYQEQ2w3uKUf9IuEAJt5uYo5ThTCZxmQfgKKr2hMntto1ENUOAyLhpuijXoRERER86kxkRYl6zv4807jyEhdyiuNxuVALvTsCvcnQ2Qn3+Uo5ikph/d2w97susc4nMaRk+MFsH4f/CjJaFJsamBFRERMpVO5WpAHl/Zm85433V7emjid8FE6vLip/qbkajn58OwG+CLLe7lZhT/XB0D2WXhqff1NydXKK41G941P/eMUL3+vERERsTY1JtIi/PMAbNgHTbl8xO6ANZ/DriMeT0ss4ng+rPjYOGLSFPty4PVtRq2IiIiIOdSYiOX9Kw/+caD+McsfMB71eecLyCvyWFpiERWV8Pp2qKiqe4w79fHNSdh00LO5iYiIiPvUmIilVVTC2zs9E+vSkRPN2tW6fPg1FFzwTKyPDsDJQs/EEhERkcaxdGPicDhYtmwZsbGxhISEkJCQQFpaGv3792f27Nlmpyc+sPsoFJV6Lt7xAvg2z3PxxFzny2HHYc/FczhhyyHPxRMRERH3WboxmTlzJk888QRz5sxh48aNTJ06lWnTppGVlUVSUpLZ6flcYGAwVY7KGsur7JUEBbbOOU8/y2wZMa3AH+tj1xHPXxeyNxsuVHg2plX4Y42IiEjLYdnGZM2aNaxevZp169axaNEixo4dy+OPP05ycjJVVVXVjcmvf/1rBg4cSEBAAGvXrjU5a++KDO/NybOuXw+XVZRQWHKa7l2jTcrKe0rK4WSR5+Nmnmqdp3P5W30AZJzyfMwqh3F/nNbIH2tERERaDss2Jk8++SS33347KSkpLstjYmIIDg5m8ODBAMTGxvK73/2OYcOGmZGmT00c+jAbdq3iQNZ27A4750sLWfHBI/TuFkfMdUPMTs/jjhd4J25FVeOmHG4p/K0+nE7I9VKN5OR7J67Z/K1GRESkZbHkDRZzc3NJT09n4cKFNdbl5OQQFxdH27ZtAXjwwQcBWLJkSaNew+aju6kt++kWEvqO8UisWxMfoKKylBffn8fpomxC24QRH53CEzM/JDDQc5syLW0rN00b67F4TRWXMpPxs/7gsqyhmZXqWr/gLdefhyaPI/ebLc3IrmGe3Pbu8FV9uMvbdRQYHMLPXi9zWeap+nj+pVf5wauzmpFd83mjfnxZI1b5HDHbpmeNw7O++pvTEjzypt4TqZvqo3VyOt07VcWyjQlAZGSky/KysjLS0tKYNGmSGWlZwqThs5g03NwdJp/x4odSa/3A86f68OY2tNksezC52fypRkREpGWxZGMSEREBQEZGhksT8swzz5CXl0diYmKzX8Pdzq25dr8NRbk+eSmPSUkZg3Ol+RdhHDgOf9jmuuzqb7YvufRNeF3rr7Zz+8dcF9703NzREre9J3m7jpxO+OU7rvcv8VR9zE+dyaZVM5uXYDO19PqxyueI2TYvM/7rq785LcGlf4d6T6Q2qg//ZsnGJDo6mvj4eJYuXUqXLl2Iiopi7dq1bNiwAcAvZ+TyRz26eCduUCB06+Sd2OI7NhtEhUOWFy5Uv95LtSciIiJ1s+T5CgEBAbz77rvExcWRmprKjBkziIiIYN68eQQFBREfH292iuIDndtBl/aej9snAgItWfnSWH2v9XxMG9DnGs/HFRERkfpZ8ogJQL9+/diyxfXi5OnTpzNgwABCQ0Orl1VWVmK323E4HFRWVlJeXk7btm1b7TUE/sRmg+QYWL/Ps3GTYzwbT8wzvC9sPgiePOA/MMpoikVERMS3WtT3xrt3765xGtesWbMIDQ1l+/bt3H///YSGhpKdnW1ShuJpyTEQ4sH7vnUNg/jrPRdPzBXRAeJ7ejbm2AGejecN//zyDWY80599R9IAeGfrsyx4eSRPrnmAKrtxA8XFr01mwcsjq5+zfO0cThdm84tXbuXnK0az+LXJlJbXPW/2x1+9xSMvjWDxa5O5UF7c7HEiIiINaTGNSUlJCRkZGTUufH/jjTdwOp0uj969e5uTZCOdPXeCFR8s4FDOLh55aQQLV4xi5bqaUyRfaeW6hSxcMYqXP3jEIzl4Op6nhYXAjzx4SdF9NxvXmFhFfe//qYJjTPlNNx5dOYbHVk0EoKSsiH1HtrLvyFZKyopqjXn1TuvVr3Hi7GHmPH8jr/9jcZPyaso4b7p7KIS28Uys5BiI6eaZWN42JeUXJPRNoajkDF8f2cLyeZ/Sp3s8n6X/DYD/mfmhy/jisgLCw7rxy2lv8vzcbYyIu5OPdr9Ra+wqeyUf7vw9z6duY3zidNbvfKVZ40RERNzRYhqTsLAw7HY78+fPNzsVj9mTsYmkfhPo1rkXz875hN/O3U5RyXcczTtQ6/jM3K8ov3iB387dTlXVRb49/mWzXt/T8bxlWDQM7V3/mAVvNTzj0sRBEGuhnU533v+k2Ak8l7qVp2d/BED60U/ZvOdPbNrzR/Znbasx/pJLO621vUZURAxz71zerLwaM87bOoXCg8kQUM/Zm+7UR1Q43Nn8Cf887uy5kyz6/VgWvDySF/46t8b6b49/QUL0GAASY8dzKGdnjTG5ZzKIioilTXAIXTt2ByAgIIiAgNq79NwzGfSJHExgYJARM7tmzMaMExERcUeLaUxaun1HtvKjX4Xz6MoxPLCkF79+4y72HdlKfHQKXTpG0iY4BIDAenYWvsn+nMTY8QD17gTU9lrNiWc2mw2mJcPQPk2PMT4O7rDYnAnuvP9fH9nCwhWjeG/bb6uXOa+4oqKhndambGN3n2Ol+onrAT8e2fRJDa7vAqnjPHvaoKd0ah/B07M2sXzep5SWF5N7JtNlfUlZEe1COgLQPqQTJWWFNWLsydjE0H4Tq38uqyhhw85VjBtyf62v6U7MxowTERFxhxoTHxncZzT9rx/Gc6lGMzL/Ry9TUVlKaNuw6jFZJ/dz7sJZenUbWGuMkvIi2rW9vBNwvo6dgNpeqznxrCAwAB5Ihmk3N27nsVMozBoDk2/06v0am6Sh979Lx+68/lgGy+Zs4avMzWSd3M+gPiOZkPQQE5IeIj56dMM7rU3Yxu4+x2r1k9ATFt3RuKl+bcCtA+E/JhqnDVpRcWk+//2ne3h05RjSj31KfvFJl/VhoZ0p/f7ajtLyYsJCOteI8U32Dgb2GgEY9wZY9s5MZtyxhLDQmmOvjnmhopj2zRwnIiLiDsvOytXa5BVk0b1LNABnzh2nqOQM0dclVK8vLi3gpb/9jMUPvlNnjLCQzpRWXN4JqGun4urXiugU1ax4VmGzGbMwDYyCnYdhRyYUltY+tltHuKWfcRqY2d+CFxSfYslb97ks69IhkkF9RtX7/rcJagu0BeDmAZM5djqd6OviSeg7pnpMfnEeL/w1lZKyIk4XHqu509qEbezuc6xYP907w4LbjJtzbs+AI9/VPi4kGG6KhltiIdLi97T5ZO8aRsTdxW03PcyTax7g6jnI+l1/E+t2rODesf/JV5mbGdDrZpf1dnsVDqeD4CDjQpzV//wVcb1vYUjMuOoxZ8+dcPmc6HFNP46dSsfusLM3czMDet7cqHEiIiJNocbER7JPHaRXZBx2hx2bLYA9GR+RFDsBMHYcnvrzg8ya/CxdOkZWLysuzSe8w+WLIgb2Smb9zldISZjK3szNTBz6cK3jrn6txsRrCTqEwIRBxulZBRcgtwBKyo3GpVMo9Ohq/NcqunSM5LnUrTWWZ+Z+Ve/7X1p+nnYhHQA4eOwz7hxZ8/qqhnZaG9rGzakzq9ZPYADc2Mt4lF006uNsCdgdRkMSFW40rgEt5HjxkJhxPP32Q+w4+Lda14eHXcvg6NEseHkk13buyd2jFrisP5SziwE9hwPGqX9/2fo0A3uN4LP09xmTcC8/HJHK028/xNOzNhHw/ZsSFBjMHcNn8fMVowgLDef/3b8GwO1xIiIiTaHGxEeOnT7IwF7JVFZVUFTyHXsPf8yUlEUApO1/l4zjX/Lq+scA+Pc7nqRT2DX8ZcvT/HzK/1bHiO2RSHBwCAtXjCK6ewI39BzGibOHa4y7+rXyi/Mov3jBrXgtic1mTP/bNazhsVZU2/tfUHyKjV/+gQdufZwDR7ez+p//RXBQWwb1Hlm9c3mlhnZaG9rGpwqPNbnOWkL9hLaB2EiINTuRZoiJGsL/Puo6Ica5C2t5e8tTREXEktA3hfvGPsZ9Yx9zGbP4tcl06didwMAgRg66G4CITtex8amLLuMcDge9ug2sbjYumZA0nQlJ0xs9TkREpKlsTqfTk/cmk6vsfhuKcmsuT9v3DikJU+t83vb97xHWLtzldAtfjAPo3AOG3tfgMGlAXdve27btN3Za50x+joS+KTXWnzh7mKf+/CCj46cQGd7b4/VzieqoecyqH0/R9jdsXmb8d/wic/Owkksz5C1/wNw8xJpUH/5NR0xMUl9TAjAq/v+4FcfT46TlGx1/D6Pj76lzfVREDC/Ob9wMWqofERER8TY1Jl7W4VqzM2i8lpizFfn7++jvv39ztfT3r6XnLyIivqfGxMv6u3/mi7Qy2vbSHKofERHxNy1kXhoREREREWnN1JiIiIiIiIjp1JiIiIiIiIjp1JiIiIiIiIjp1JiIiIiIiIjp1JiIiIiIiIjp1JiIiIiIiIjp1JiIiIiIiIjp1JiIiIiIiIjp1JiIiIiIiIjp1JiIiIiIiIjp1JiIiIiIiIjp1JiIiIiIiIjpgsxOQERExGzffgLnv/Ne/N1vez5mh2uh/zjPx/U33t72Vqc6ar6WXENW2/5qTERExO+d/w6Kcr0X35uxpXm8ve2l9VMNeY5O5RIREREREdOpMREREREREdPpVC4RERHxmgsVsP845OTDqXOXl/95J1zfBeKvh46h5uUn5nI44F95kHUGcgsuL/9DmlEf/bpDr65gs5mXo/iOGhMRERHxuHOlsH4ffJUNVfaa63cdMR5/3Q0JPeEHCRDRwfd5ijnsDtieAWmHoLC05voDucZjw37oEQ4TBxtNrLRuakxERETEo/YchbVfQlllw2MdTtibDem5cGci3BKrb8dbu++K4c0dxlE0d+QWwmvb4MaeMHUYtGvr3fzEPLrGRERExMMeXNqbzXvedHt5a/LJN/CnHe41JVeqtBvNzN+/BqfTK6lZhj/Xx8lC+N1H7jclV/o6B17aDCXlns/LSvy5PtSYiIiIiEfsOQbr9jYvxiffwLZvPZKOWMz5Mlj5iXHdUVOdLIJX04xTwaT1UWMiIiIizXauFNZ+0fC45Q8Yj/r8/Ws4XeyRtMQinE5490s438DRDnfq49hZ2HLIc7mJdVi6MXE4HCxbtozY2FhCQkJISEggLS2N/v37M3v2bLPTExERke9t3N/407fqUmWHdV95JpZYw+HTxuxsnvKP/a3/lC5/ZOnGZObMmTzxxBPMmTOHjRs3MnXqVKZNm0ZWVhZJSUlmpyciIiJAaQXsPubZmN+cgPwSz8YU83ya4dl4VQ7YecSzMcV8lm1M1qxZw+rVq1m3bh2LFi1i7NixPP744yQnJ1NVVUVSUhIVFRU8/PDDREVF0blzZ8aNG8ehQzq2JyIi5goMDKbKUfPwQZW9kqDAYBMy8q4DubVPCdwcTozZulojf6uPiiqjRjztq2Oej2kF/lYfV7JsY/Lkk09y++23k5KS4rI8JiaG4OBgBg8eTFVVFTExMXzxxRfk5+czfvx47r33XpMyFhERMUSG9+bk2cMuy8oqSigsOU33rtEmZeU9x5sww5JbcQsaHtMS+Vt9nCgwpoX2tFPn4GKV5+Oazd/q40qWbExyc3NJT09nypQpNdbl5OQQFxdH27Ztad++PYsXLyYqKorAwEDmz5/PgQMHKC/XSYciImKeiUMfZsOuVRzI2o7dYed8aSErPniE3t3iiLluiNnpeVzeuYbHNClukXfims3f6uOUl+rD4TTuidLa+Ft9XMmSN1jMzTWO90VGRrosLysrIy0tjUmTJtX6vB07dtC7d29CQkIafA2b7t4kIiLfW/bTLST0HeOxeLcmPkBFZSkvvj+P00XZhLYJIz46hSdmfkhgoGf+9KalbeWmaWM9Equ57v3NLiL7DnNZ1tDMSrWtX/CW68+Hs7Kx2Xo3L7kGeHrbu8MX9eEuX9TRjbc9Qsr05S7LmlIfULNGht18C3mZO5qenAe0xM+PS3z1OeJ08+ZElmxMIiIiAMjIyHBpQp555hny8vJITEys8ZzCwkLmzZvHkiVLfJaniIhIXSYNn8Wk4bPMTsMn7JXeOVOh6mKZV+JagerDM6oqW2eN+FN9XMnmdLeF8SGHw8GQIUPIy8tj2bJlREVFsXbtWjZs2EBOTg47d+5k+PDh1ePLysqYOHEio0ePVmMiIiKNtvttKPLCxbne1LkHDL3P7CwM7+yCHYcbHgeXvwm/+pvv2sRfDzNHNz0vd7TEbe9JvqijI6fhxc3ujW1Mfdhs8NRUaGvy1+wtuYas9DkCFr3GJCAggHfffZe4uDhSU1OZMWMGERERzJs3j6CgIOLj46vHVlVVMXXqVGJjY9WUiIiImOD6rl6K28U7ccW3oroYTYSndetoflMinmXZzdmvXz+2bNnismz69OkMGDCA0NDQ6mU/+clPcDgcrFq1ytcpioiICDC4B6z9EuwOz8a9sadn44k5QoJh4HVw8IRn4w7p5dl4Yj5LHjGpy+7du11urJidnc3q1av55JNP6Ny5M2FhYYSFhZGTk2NiliIiIv4lLASGeLiJ6N8druno2ZhinltiPRsvwAY3x3g2ppivxTQmJSUlZGRkuFz43qtXL5xOJ2VlZZSUlFQ/evbUVywiIuIZZ8+dYMUHCziUs4tHXhrBwhWjWLluYb3PWbluIQtXjOLlDx7xSA6ejucNkxI8d1pNoA3utNisqPVtg1MFx5jym248unIMj62aCEBJWRH7jmxl35GtlJQV1Rrzn1++wYxn+rPvSBpnz50kdXkik/5vCHa7cXOOE2cPM+f5G3n9H4ublFdTxnnLgOuMh6eMj4NOoQ2PM9uV27jKXsl/vJjMDx8P48QV9ylZ/NpkFrw8svrn5WvncLowm1+8cis/XzGaxa9NprT8fJ2v8fFXb/HISyNY/NpkLpTXPX+yu+PM1GIak7CwMOx2O/Pnzzc7FRER8SN7MjaR1G8C3Tr34tk5n/DbudspKvmOo3kHah2fmfsV5Rcv8Nu526mqusi3x79s1ut7Op63dAmDu5IaHrfgrYYvbL49Hq4L90xenuDONkiKncBzqVt5evZHAKQf/ZTNe/7Epj1/ZH/WtjpjT0n5BQl9U+jYrgvPzP6YAT1vrl4XFRHD3DuXNyuvxozzJpsN7h0OoW3qH+dOfUSFw8RBnsvN2y5t48CAIH7z8N8YNfgel/X/M/NDl5+LywoID+vGL6e9yfNztzEi7k4+2v1GrbGr7JV8uPP3PJ+6jfGJ01m/85VmjTNbi2lMREREvG3fka386FfhPLpyDA8s6cWv37iLfUe2Eh+dQpeOkbQJNu6TFRgQREBAYK0xvsn+nMTY8QAkxo7nUPZOt1+rOfGsIDmm+TuMI2KNb8OtxJ1t8PWRLSxcMYr3tv22epmTyxOfnj13kkW/H8uCl0fywl/n1nh+m+AQOrRrXDfmbm1YpYY6t4M5Y41rTprqmg4weywE1f7Pz1QNbWObzUZ4h271xsg9k0FURCxtgkPo2rE7AAH1fN7knsmgT+RgAgOD6t227o4zmxoTERGR7w3uM5r+1w/juVSjGZn/o5epqCwltG1Y9Zisk/s5d+EsvboNrDVGSXkR7doaF0e0D+nE+bJCt1+rOfGsYlICTLkJghu54xgYAD/4/rlWuwdyQ9ugS8fuvP5YBsvmbOGrzM1kndzPoD4jmZD0EBOSHiI+ejSd2kfw9KxNLJ/3KaXlxeSeyfR6Xo0d5wu9I2D+BIjs1Pjn3tAd/mOCdU/h8sQ23pOxiaH9Jlb/XFZRwoadqxg35P5ax5eUFdEu5PK2LamrBtwcZzbLzsolIiLia3kFWXTvEg3AmXPHKSo5Q/R1CdXri0sLeOlvP2Pxg+/UGSMspDOlFcb52xcqigkL7ezWa0V0impWPCu5pZ9x8frf98KBXHDUc8c0G8bYOxOhe2dfZVi7guJTLHnL9aYOXTpEMqjPqHq3QZugtkBbAG4eMJljp9OJvi7e5W7g+cV5vPDXVErKijhdeIz84pPNztfd2rBaDUWFw6I7YPNB2J4BFyrqHx/RASbEwbBo6zWtVyouzW/2Nv4me0f1jRWdTifL3pnJjDuW1L1tQztTWn5527Zv5jizqTERERH5Xvapg/SKjMPusGOzBbAn4yOSYicAYLdX8dSfH2TW5Gfp0jGyellxab7L6RkDeyWzfucrpCRMZW/mZiYOfbjWcVe/VmPitQQRHWDGaCgqhb3ZkJMPp4rgot04mtKtk3Gfkht7GmOtoEvHSJ5L3VpjeWbuV/Vug9Ly87QLMX6Jg8c+486RNa+H/WTvGkbE3cVtNz3Mk2seABp3f+vm1JoVaygo0LiW6NY4OHAcss7A8QK4UG40H53bGfXRrzv0izRm4bI6T2xjh9NBcJBxIc7qf/6KuN63MCRmXPWYs+dOuHyJ0eOafhw7lY7dYWdv5ubqa5TcHWc1OpVLRETke8dOH6R3tzgqqyooKvmOvYc/pl+PoQCk7X+XjONf8ur6x3h05Ri+OfY5pwqP1ZgxKbZHIsHBISxcMQqbLYAbeg6rddzVr5VfnOd2vJakczsYOwB+PBIemwz/dSf8cjLMGGVcS2KVpqQ+tW2DguJTvPWxcWPnA0e3M3d5Eo+8NIKuHa9jQM/hNWIMiRnH2m3P8es37qL84oUa66vslfznK+PJytvHL1+9jUM5u1zWN6fWrFxDwYGQ2BvuuQkW3gaL74TH/w3mjYd/SzRO32oJTQk0vI0BnvjTVPZkfsQzb/+YHekfuKw7lLOrunbOnjvJX7Y+zWfp7/PoyjH8fcdKAJ5++yEcjss3DAoKDOaO4bP4+YpRfLR7NZNvntOocVZjczqdjWvnREREWpndb0NRbs3lafveISVhap3P277/PcLahbt8o+mLcQCde8DQ+xocJg2oa9v7wrb9a3l7y1PMmfwcCX1Taqw/cfYwT/35QUbHTyEyvLfHawhUR55QXw01tI3BmC64TXAoU1IW0bVDd64Nr/22Fw6HgxXrHuFnd71Ybz7ujgPrbX81JiIi4vfM3DltKqvtULRULXHbe5LqqPlacg1ZbfvrGhMREfF7Ha41O4PGa4k5W5G/v4/+/vt7Qkt+D62Wu46YiIiIiIiI6XTxu4iIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImE6NiYiIiIiImO7/A37B00u9rCN9AAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 1049.8x204.68 with 1 Axes>"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Ansatz_gate = ansatz_circuit.to_gate(label=\"ansatz\")\n",
    "CAnsatz_gate = Ansatz_gate.control(1) \n",
    "\n",
    "new_qc = QuantumCircuit(3)\n",
    "new_qc.append(CAnsatz_gate, [0, 1, 2])\n",
    "\n",
    "new_qc = transpile(new_qc, backend)\n",
    "new_qc.draw('mpl')"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "we can easily get the inverse of the circuit by invoking the inverse function in the qiskit."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAswAAACoCAYAAAD0K7h/AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAmRUlEQVR4nO3deXhU9b3H8Xc2sgIBIltYQwJCIJEAIgoEFCnl2itawSLgRVqgAamoVHuvttrHAoqAFAUu1NsKV9EqWqWCCghEELECsgSUhJ2EsCSQhKxklvvHXJZhkpNtMkvm83qeeWDO+Z3ffDP55ZzPnDmLn9VqtSIiIiIiIhXyd3cBIiIiIiKeTIFZRERERMSAArOIiIiIiAEFZhERERERAwrMIiIiIiIGFJhFRERERAwoMIuIiIiIGFBgFhERERExoMAsIiIiImJAgVlERERExIACs4iIiIiIAQVmEREREREDCswiIiIiIgYUmEVEREREDCgwi4iIiIgYUGAWERERETGgwCwiIiIiYkCBWURERETEgAKziIiIiIgBBWYREREREQMKzCIiIiIiBhSYRUREREQMKDCLiIiIiBhQYBYRERERMaDALCIiIiJiQIFZRERERMSAArOIiIiIiAEFZhERERERAwrMIiIiIiIGFJhFRERERAwoMIuIiIiIGFBgFhERERExoMAsIiIiImJAgVlERERExIACs4iIiIiIAQVmEREREREDge4uQKQuzBZIy4RTueDvB11aQdfWtv+LAGRdggOn4YoJohpDUicICXJ3Ve5Tkg/nfoTyUmgUDq1vheAId1flPkVlsOcEXCqyjYvEDtCqqburEk9hscChM3Aix/a8cxR0bwv+Prq70WqF/DOQexwsZoiIgpZdIcAH1ql+VqvV6u4iKmOxWFi4cCHLly/n9OnTdOvWjcWLFzNlyhSSk5NZsWKFu0sUNzqcDW/vgMul9tOjIuCxwRDdzD11iWcoLIWV2yHjnP30oAC47zZIvtUtZbmNxQQ/bITsgzfN8IP2vSFuiG+FAKsVNqTZHmaL/bz4aJhwl29/sBI4fsG2Dskrtp8eGQYTB0GnKPfU5S4lebBvLRSet58eGAy33mv78N2QefTqcdKkSbz00ktMnTqVzz77jDFjxjB27FiOHTtGnz593F2euNHxC7Biiy0U3Sy3CN7YCDmXXV+XeIYrJli2GY6cc5xXboZ/7Ibt6a6vy53SPqsgLANY4fQeOPyly0tyq40H4bP9jmEZ4GCWbf1S0TzxDWcuwdIvIb/EcV5+MSzZBNl5Li/Lba4Uw66/Q+EFx3mmMkj7FC4ccX1druSxgXn16tWsXLmStWvXMmvWLIYOHcpzzz3HgAEDMJlM1wLzuXPnGD58OGFhYSQmJvL999+7uXJxhfX7wGKFir4esVqhtBy+POTyssRDfH/SdiiG0ddn6/bagrUvyM+G84eN22Ttg+JLrqnH3YrL4IsDxm2OXYBDWa6pRzzP5wfAZLZtT25mxTavqjHUkGTuhbLLGK5UM1Irfr8aCo8NzHPnzmXEiBEkJyfbTY+NjSUoKIhevXoBkJKSwq233kpubi7Tp0/noYcewmw2u6NkcZGLhbav2Y3+Lq3Ad8dtexPF93xzBKo6jL2k3Hb8uy84k0bVbwhwpqI90A3QnpNV7z32A7496pJyxMMUlcGBzKq3MftP2z58+YKs/VW3Kb5k+3DeUHlkYM7MzCQtLY3Ro0c7zDt16hTx8fEEBwdz+fJl1q1bxx/+8AdCQ0OZMmUKZrOZnTt3Vvkafn5+enjpo2fSXdUaRyYzNLsl2u316uH6R1p6puHG7qrJj//W7bW64vHhu59WuSPBbDbxlyWr3F6rKx7P/n4OVqtxYrYCW77Z5/Za9XD9o2Nsz2rtKbVYoX1MD7fX64pHWWE1VqjAqJEPu73WmjxqwmMDM0Dr1q3tppeUlJCamnrtcIyMjAxatGhBVNT1I+979erFoUP6Lr4hu1KSX612VquF8lIdyOyLyorzsFqqPgD1SnH1xpK3KyrNx3h/mW0nQlE1/7a83ZWSAvz8jDd/FouZsqI81xQkHqWsBuuFMh/5mykuLahWO9u6pmHyyMvKXQ3A6enpjBw58tr0efPmkZ2dTVJSEgBFRUU0adLEbtkmTZpQWFj1RyEPvjiIVMFqhTn/tJ3UV9lv0c8Purf1p6y4en/k0rBsOADrq/gK0d8PdqxfQeOQhn+1nfNHYP/Hxm38/QP40xszeOPjGS6pyZ0uXIbZa43b+PsHMPPRZNb8SdsKX7Twczida7yN6dgCCi/6xoHuP2yo+rCMwBD4Nu1z/D0yWdadR+5hjomJISEhgTlz5rBq1Sq+/PJLUlJS+Otf/wpwbQ9zeHg4ly/b70EsKCggIsKHLyrqA/z84N6eVewvs8LdPVxVkXiaAbEQGmR82O4dsdA4xGUluVVUDIS3wPANadoGmka7rCS3uqUx3Nah8vl+frax0bez62oSzzIsvopjmK22Nr6ifRL4BRi36diPBhuWwUMDs7+/Px988AHx8fGkpKTw2GOPERUVxfTp0wkMDCQhIQGAuLg4cnJyyM3NvbZsWloaPXooKTV0t8fASNswcMgAfn4wdgDEtXJ5WeIhGofCr++GkEb2068estazHTzoQ1em9PeH3j+HsMibZvz/+xFxCySOuv7++IKxd1xfR9z8c0cEw7R7dB1mX5bQHh4wWEf8vK9tPeIrIqIg8f4KAvH//+20S4ROt7u8LJfy6BuX3GzChAns27eP/fuvfy/w4IMP0qFDB15++WXefvtt5syZQ0ZGBgEBVXwUkgYhOw++zrh+Td17esCdcdBCXzIIUHwFvjtmu+4yQO+Otr3Pca18KxxeZTbB+XQ4uN72PCoG2sTDLbHg74OrTIsVfjwDO4/arngAMLqfbc9ysMKyAOcLbNuY1B9tz5NvhYFxcEsT4+UaqrIiOHMAjm63PW/T0xaWm7Zxb12u4JF7mCuza9cuhxuWLFu2jIMHD9KsWTMWL17Mhx9+qLDsQ9pEwkP9rj//WW+FZbkurJH9Hf3+Y6Dt1um+GJYBAgKhzQ1fwN32ILTq5pthGWzHsfeIhkmDr0+7q6vCslzXson9nuYH+vhuWAYIDofOd1x/Hj/CN8IyeOhJfxUpLCwkPT2dadOm2U1v1aoVGzdudFNVIiIiItLQeU1gjoiI0A1JRERERMTlvOqQDBERERERV1NgFhERERExoMAsIiIiImJAgVlERERExIACs4iIiIiIAQVmEREREREDCswiIiIiIgYUmEVEREREDCgwi4iIiIgYUGAWERERETGgwCwiIiIiYkCBWURERETEgAKziIiIiIgBBWYREREREQMKzCIiIiIiBhSYRUREREQMKDCLiIiIiBhQYBYRERERMaDALCIiIiJiQIFZRERERMSAArOIiIiIiIFAdxfgrQ5vhsvn3V2F+zRuCd3udncVIiK+wRe2Odqu1I47x8au91z7eu4cIwrMtXT5PORlursKERHxBdrmSGXcOTZ8aUzqkAwREREREQMKzCIiIiIiBhSYRUTEZ10xwanc688LS91Xi3geixXOFVx/fq7ANk18j45hFhERn1JUBv86Bt8dg+x8sN4QgJ7/ECLDIKE93BUHrZq6r05xD4sFDp2BHRlw9DyUma7Pm/tPCA6ELi1t46N7W/DXrkefoMAsIiI+wWKFr9Phn3tte5Yrk1cMXx22Pe7oAvcnQWgjl5UpbnQqF979xvZBqjJlJlugPnQG2kTCIwOgfXOXlShuos9FXmT8nE5s2v12taeLiIhNWTms2AIf7jIOyzfbeRTmrYPsvHorzWP4+jYm9UdY9IVxWL5Zdh689jlsO1xvZXkUXx4jCswiItKgXTHB8i3wY3btlr9UDG9sgrM1CFLiXbb+CP/YXbvjky1W2wex1B+dX5d4DgVmERFp0D7dC8cuGLdZNM72qExRGazcBiazU0sTD3D8Anyy27hNVeMD4OPdcCLHeXWJZ1FgFhGRBuv4Bed9XZ6dDxsPOqcv8QwmM7y7E5xx4QsrsPobfahqqDw6MFssFubPn09cXBwhISEkJiaSmppKt27dmDJlirvLExERD7fpoHPC0FWpP9pfNUG8277TcL6g6nbVdb4ADvjQ3e98iUcH5kmTJvHSSy8xdepUPvvsM8aMGcPYsWM5duwYffr0cXd5LhcQEITJUu4w3WQuJzAgyA0ViYh4rouFcCjLuX2WlsOeE87t01P44jbm63Tn97m9Hvr0FL44Rq7y2MC8evVqVq5cydq1a5k1axZDhw7lueeeY8CAAZhMpmuB+YUXXqBHjx74+/uzZs0aN1ddv1o368SZnCN200rKCrlUeI42LWLcVJWIiGfKOOfcvctXpZ+th049gK9tY66Y4Hg9HHN8/AKUN9DDMnxtjNzIYwPz3LlzGTFiBMnJyXbTY2NjCQoKolevXgDExcXx5z//mdtvv90dZbrU8L4TWf/tCg4c24bZYuZy8SWWfvIEnVrFE9u2t7vLExHxKJkX66ff07lVt/FGvraNybpkf9MaZ7FY4cwl5/frCXxtjNzII29ckpmZSVpaGk8++aTDvFOnThEfH09wcDAA48ePB2D27Nk1eg0/P7861Tj/11tI7DKkTn3U1D1J4ygrL+b1f0znXN5JQhtFkBCTzEuTPiUgwLW/ytTUrfQbO9Slr2nkibdta726/l6lYdL4sLfxVd94P/5t5kfE9n3AblpVVzqobP7Md67//+zFEvz8wupYXc24Ypvj7m2Mq7crsf1+zr89Yf/NtDPGB8A9P32Ao7s+rn1xNeDKPNLQxoi1Bp+YPDYwA7Ru3dpueklJCampqYwcOdIdZXmEkf0nM7L/ZHeXISLi8fyorw8EDfeDhk9tY+r1A6PGSEPjkYE5KioKgPT0dLtwPG/ePLKzs0lKSqrza9TkU0VFdr0HeT58Jmxy8hCsy+rj6MDaufrpvq6/V2mYND7sbZpv+7ehvx9//xa+sT/c0mFP4FVX9xxWNv9GUZEhLn/vfGGb4+rtSsY5WLLJfpozxgfAF+s+okvL2tdWE74wNq5yZ/bwyMAcExNDQkICc+bMoXnz5kRHR7NmzRrWr18P4JNXyBARkZpp16ye+m1eP/2Ka9XX+PADouupb3Efjzzpz9/fnw8++ID4+HhSUlJ47LHHiIqKYvr06QQGBpKQkODuEkVExMN1aVVP/bpoz6HUr9BG0LYegm10Mwhp2FdY80keuYcZoGvXrmzZssVu2oQJE+jevTuhoaHXppWXl2M2m7FYLJSXl1NaWkpwcHCDP5lFRESMtW4KnW+xXebLWQL9oV9n5/Un7nVnLKz5zrl9Dohzbn/iGTxyD3Nldu3a5XA4xuTJkwkNDWXbtm088sgjhIaGcvLkSTdVKCIinmRod+f2d3sMRIQ4t8+aWrb2SZ5cOoglnzzhMG/VhheZujCRp5cNYU3qQgAKS/JYteHFa/+vyLz3JjJjcX+KSvIxm028/O4EZi4ZyHubXwYgK+cIUxfext8+f77Sut7f+iozlwxk7upxmMyON7eoaTtX6NvZub/PxiHQt5Pz+hPP4TWBubCwkPT0dIcT/t566y2sVqvdo1OnTu4psgJGK7azF08w+o+teHrZEJ5dMRywrcz2Hd3KvqNbK12xffHdWzw2rxv7jqZW+BrVWbEZ1VWbdiIiniihPSR2cE5fTUPhPjdfajYjcw+lV4p4bdo2TKYrHD7tuHt06n0LWJCylYeSnwJg0563OXB8G6s2vMiF/MrPDvvdI+8QHtqUHYfW0qFldxZN307aie1cLDhLdFQs0+5fVOmyeYUX2Ht0C4umb6dzmwS+Tvu4Tu1cJSQIxjjxNg4P94dgLzkcIyc/i/9686c88cadPLl0EMvWOl7K19m8OXt4TWCOiIjAbDYzY8YMd5dSbdVZsfWJu5cFKVt5ZcoGANKOb2fT7v9l4+5V7D/2VaV9j07+LYldkit8japWbNWpqybtREQ82eh+cEtj4zYz3zG+AkKgP0y4C8IaObe2mjp08huS4oYBkBQ3jB9O7nRo8+b6Z3lm+TCOZO21m371Mnvf/fg5Ty1LZtqf+7Jx1yqH5X+44TUSuwzlcGbV6/7Dp/9FYsyQ63WdcqyrJu1cKaE9DOxq3Kaq8QEwuBv0bOe8uurb7vSN3JM0nlenbua1advIKzzP8ewD9fZ63p49vCYwe6PqrNj2Ht3Ck0sH8eFXr12bZr3hZq45+WeY9d9DmblkIIs/mlar16jtMrXpW0TE00SEwPRh0KpJ7ZZvFAi/GgKx9XQSYU0UluYRFmz7QcJDmnK5xP6WcqMG/oalM3fzmweXseQT2w6mYUnj6dV5EBOGv8AtTdvRK2YwC1NSef3xnaz7doXja5TkERZy/TUKi6u+bZ3DMiUVL1Pddq72YF8YWIdjjwd1hVEeegGvfUe38sAfmvH0siGMm92RF94adW36nfH30yjIdkxKgH8g/v4BTuv/Zt6ePRSY61FVK7bmTdrwt2fTmT91C3syNnHszH56dh7IvX0e5d4+j5IQM5im4VG8Mnkji6Zvp7i0gMwLGTV6jdrUVZe+RUQ8UWQYPPVTW7CpiS4t4ZmRcGub+qmrpiJCIikuKwCgqKyAiNBIu/lNwmzXvGt3y/X0FxEayaPDX7z2/4zM3TyzfBi/XXEPJ88ddHyN0EiKS22vUVxaQPhNr1FhXTctExFS8TLVbedq/n7w837w6F0QXoNvESKCYeIg27L+HnqtgV6dB9Ot/e0sSNlKQkwyMx5YgtVqpay8mNDgCACOndlPflEOHVv1cEr/FfH27OGxV8nwJhcLzjL7nV/YTWveuDU9Ow8yXLE1CgwGbLf4vqP7fZw4l0ZM2wS7W1zmFmSz+KMUCkvyOHfpBLkFZ+z6qGrlWZHqLlObvkVEPFVwoC3Y3BEL29Nh9wm4YnJs5wd0awN3xUF8O88KQj06DmDdzuUkJ47h+4xNDO870W5+UWkB4SFNyC/KwWyu4IcD3t86j6dGv0lU02gmvuK4W7V7xwF8f+RLbu1wO/uObmFo77EObXLys4hqGn3tedf2/Vi7YykPD32GPRmb6N7xjhq18wR+fpDUyfa7//Yo7MiAnMKK20ZFwF1dbSeBhge7tMxKVZZFJo74E22axwBwIf80UU2jOZK1l5i2iQAUFF/kjY8f5/nx7zut/4p4e/ZQYHaC5k1asyBlq8P0jMw9hiu24tLLhIXYDqw7eOJr7h/oeHz25u9Xc2f8KH7SbyJzV48D7O9wU9XK02w2UVCcS7PGrQyXqW47ERFvF93MdnLWQ/3gbD6cuQRlJgjwh5ZNPPs6unHtkggKCuHJpYOIaZPIrR1sZ6y98fEMHh/1On/59LccP5uG1WrhlyNfrrCPu3o+wAtv3U+XtrcREep4IeIBPX7Gq3//kJlLBnL7rSNp0cRx9/rc1ePstnvNIlrSK2YwM5cMpGVkBx4cNLNG7TxJeDDc3cN2hZW8Yjh9ES6XAH7QJATat7CdAOppV6+tLIvsSPuEjq3jMVvM+PnZDizYnb6BPnH3/v8VUcYz+b5Xad6k9bVlKsoENem/IWYPBeZ6VNGK7WLBWT777n8Yd89zHDi+jZVf/J6gwGB6dhpI9w79HfroHXs3r7z3KDsOflzt17jR2Usn+PuWV3hq9F8Ml8nKOVKtdiIiDUWAvy0ce9td2abf/2eHaY+Peh2AmQ8tr3L5EbdPYsTtk+ymNQ2P4tW/T2T2pHWEhzblPx+xP8MtK+cIb67/HYMTRlNQlEtCTLJDv78Y+iy/GPrstefVbeep/PygWbjt4c1OnDtIj44DKDeVkVd4ntyCbDKydjM6eRZb9/2d9NPf8eY62+/jlz+dS49OAyrMDjXpv/RKUYPLHn5Wq9U9N+X2cu68d/tX+9fw3paXmXrfAhK7OK6MsnKO8PK74xmcMJrWzToREdaM3rF3G/a5bf+H1Wp3VWQ76PuLqtu5ytWzlxeNc28d4pk0Puxtmm/7d9gs99Yh1efObY6reNp2xVvUZmyk7nuf5MQxlc6vaSao7fLelD20h9kLDU54iMEJD1U6Pzoqltdn1Oys0kEJP69rWSIiIuIFjMIy1D0TVHd5b8oeCsy11LiluytwL1//+UVEXMkX1rm+8DPWB19639z5syow11K32n1LISIiUmPa5khlNDZcQ9dhFhERERExoMAsIiIiImJAgVlERERExIACs4iIiIiIAQVmEREREREDCswiIiIiIgYUmEVEREREDCgwi4iIiIgYUGAWERERETGgwCwiIiIiYkCBWURERETEgAKziIiIiIgBBWYREREREQOB7i7AWx3eDJfPu7sK92ncErrd7e4qvIcvjBeNidpx59jY9Z7rX1PjRES8kQJzLV0+D3mZ7q5CvIXGi1TGnWNDY1JEpHp0SIaIiIiIiAEFZhERERERAzokQ7ySxQLp5+DoOTh98fr0N1OhXXPo2go63wJ+fu6rUdzrfAEczILTudenvbEJ2kZCpyjo2Q4aaQ0oIiLVoM2FeBWLBb7OgK0/Qm6h4/y0TNvjc6BNJNwbD707Kjj7khM58Nk+OHzWcd6Rc7bHV4chtBEMiIXhPSEkyPV1ioiI91BgFq+RWwhv74DjF6rXPjsPVn0N35+Eh/tDREi9liduZrbAun2w5RBYq9G+5ApsPmQbH+MGQGyrei9RRES8lI5h9iLj53Ri0+63qz29ITmXD4u+qH5YvtGBTFi8EQpKnF+XJ/Hl8WG2wN+22QJwdcLyjS4VwbIv4cDpeinNo/jyGBERqQsFZvF4RWWwbDNcLq19H+cLYMVWMJmdVpZ4kDXf2Q7FqS2zFd7aDqdyq24rIiK+R4FZPN4/dkFesXGbReNsDyOZF2FDmvPqEs9wKAu+OWLcpjrjw2yB1d9AuT5UiYjITTw6MFssFubPn09cXBwhISEkJiaSmppKt27dmDJlirvLExc4kQO7Tjivv00Hqw7f4j0sVvhol/P6O5sPX6c7rz8REWkYPDowT5o0iZdeeompU6fy2WefMWbMGMaOHcuxY8fo06ePu8sTF9ju5PBisVa9N1K8x+FsyKngail18XWGbZyIiIhc5bGBefXq1axcuZK1a9cya9Yshg4dynPPPceAAQMwmUz06dOHsrIyJk6cSHR0NJGRkdx999388MMP7i693gQEBGGylDtMN5nLCQxoeNfFMplh7ynn97v7hPP79AS+Nj6gfn6XFy7bX7u5IfHFMSIi4gweG5jnzp3LiBEjSE5OtpseGxtLUFAQvXr1wmQyERsby7/+9S9yc3MZNmwYDz/8sJsqrn+tm3XiTI797tGSskIuFZ6jTYsYN1VVf87m189JejmXofiK8/t1N18bH1B/J+ndeDOchsQXx4iIiDN4ZGDOzMwkLS2N0aNHO8w7deoU8fHxBAcHEx4ezvPPP090dDQBAQHMmDGDAwcOUFpah8speLDhfSey/tsVHDi2DbPFzOXiSyz95Ak6tYontm1vd5fndGfz67HvvPrr2118bXyYLXChoH76zs6rn37dzdfGiIiIs3jkjUsyM23Xh2rdurXd9JKSElJTUxk5cmSFy+3YsYNOnToRElL1HSr86njrt/m/3kJilyF16qOm7kkaR1l5Ma//Yzrn8k4S2iiChJhkXpr0KQEBrv1VpqZupd/YofX6Gj2HTuaeX66wm1bVlQ4qmz/zHfvnQ++5l1Npm+pQXc24Yry4e3y4YkzcKDA4jOn/U2Q3zVnj4y//8zfG9J9Uh+qqz5XrEnePEXD9OBERqYzVWv0TVjwyMEdFRQGQnp5uF47nzZtHdnY2SUlJDstcunSJ6dOnM3v2bJfV6Q4j+09mZP/J7i7DJUzl9fdNgelKw7yLiS+ND4vJdlyN1Wqt8wfgm5muNMxvqcC3xoiIiLP4WWsSr13EYrHQu3dvsrOzmT9/PtHR0axZs4b169dz6tQpdu7cSf/+/a+1LykpYfjw4QwePNhlgXnXe5BXhxsleLvIdtD3F/X7GqdyYeHn1Wt7dc/hzXsKKzP7IQgPrl1dteEL48UVY+Jms9faTtKrSk3Hx4N9YXC32tdVE74wNm7kjnEiIlJXHnkMs7+/Px988AHx8fGkpKTw2GOPERUVxfTp0wkMDCQhIeFaW5PJxJgxY4iLi2vwe5d9TdtICKiHEdoiwrVhWepP++be1a+IiHgnjzwkA6Br165s2bLFbtqECRPo3r07oaGh16b96le/wmKxsGLFipu7EC8XGAC92jn/0nK3dXBuf+I+vTvCnpPO7bNZOHRs4dw+RUTEu3nkHubK7Nq1y+6GJSdPnmTlypVs3ryZyMhIIiIiiIiI4NSperh4r7jFwK7O7c8PuDPOuX2K+/SIhmZhzu3zrjjw96o1o4iI1Dev2SwUFhaSnp5ud8Jfx44dsVqtlJSUUFhYeO3RoYPn7EJctvZJnlw6iCWfPOEw7+zFE4z+YyueXjaEZ1cMB6CwJI99R7ey7+hWCkvyKuzzi+/e4rF53dh3NJWc/DOkLEpi5H+GYDabAMjKOcLUhbfxt8+fr1VdtWlXX7q0hIT2zutvSHfbIRnuZPSertrwIlMXJvL0siGsSV0I2MbEqg0vXvt/Rea9N5EZi/tTVJLPkazvmbygF+PndLo2vzpj4v2trzJzyUDmrh6Hyex4c4uatnOFAH8Y5cSbfkZFwCAXHbvsDDn5WfzXmz/liTfu5Mmlg1i29sl6f01vWXeIiDiT1wTmiIgIzGYzM2bMcHcp1ZaRuYfSK0W8Nm0bJtMVDp/+zqFNn7h7WZCylVembAAg7fh2Nu3+XzbuXsX+Y19V2vfo5N+S2CWZJmHNmTflS7p3uOPavOioWKbdv6hOddWkXX3y84PR/SCiimOOZ75T9QldrZrATxOM29S36rynU+9bwIKUrTyU/BQAm/a8zYHj21i14UUu5Fd+dtjvHnmH8NCmtG0Ry+IZO4lq2u7avKrGRF7hBfYe3cKi6dvp3CaBr9M+rlM7V0rsAH06Gbepzvjw94OxAyDYYw9Uc7Q7fSP3JI3n1ambeW3aNvIKz3M8+0C9vZ43rTtERJzJawKzNzp08huS4oYBkBQ3jB9O7nRos/foFp5cOogPv3rt2jQr1y9ckpN/hln/PZSZSway+KNpDss3CgqhcVgzp9dVk3b1rXEo/PpuCGtU+z6ah9v6aOTmMFSd9/TN9c/yzPJhHMnaazfdD9ul07778XOeWpbMtD/3ZeOuVQ7Lh4U0JrRReI3qOnz6XyTGDLle16mKf9fVbedqv7gDurauul1l/IBHBti+0fBE+45u5YE/NOPpZUMYN7sjL7w16tr0O+Pvp1GQ7drzAf6B+PsHOK3/m3nbukNExFkUmOtRYWkeYcFNAAgPacrlkkt285s3acPfnk1n/tQt7MnYxLEz++nZeSD39nmUe/s8SkLMYJqGR/HK5I0smr6d4tICMi9k1HtdNW3nCu2aw2+GQ3TNPhsAENcKnhhuO5nL3ap6T0cN/A1LZ+7mNw8uY8kntm9ThiWNp1fnQUwY/gK3NG1Hr5jBLExJ5fXHd7LuW+ec7FpYkkdYyPW6CisbE9Vs52pBATB5CNwZW/NlG4fYlu3b2dlVOU+vzoPp1v52FqRsJSEmmRkPLMFqtVJWXkxosO0Yo2Nn9pNflEPHVj2c0n9FvHHdISLiDF705aPnulhwltnv2F9YtHnj1vTsPIjiMtu9e4vKCogIjbRr0ygwGLAda3BH9/s4cS6NmLYJdnf9yi3IZvFHKRSW5HHu0glyC87Uud6IkEjDumrazlVaN4WnRsDmQ5B6GAqruLdEiwi4pwcMiLUd2uEJqnpPm4TZrmfW7pbrZyZGhEby6PAXr/3/wLFt/O/GP2KylHPy3EHn1BUaSU5+FgDFpQVEhETWqZ07BAXAmP62QzTW74OTuVW379cZ/u02z7nMYGXrkokj/kSb5jEAXMg/TVTTaI5k7SWmbSIABcUXeePjx3l+/PtO678i3rruEBGpKwVmJ2jepDULUrY6TM/I3MO6nctJThzD9xmbGN53ot384tLLhIU0BuDgia+5f6Dj8dmbv1/NnfGj+Em/icxdPQ6o2X1mzGYTBcW5NGvc6tq0Hh0HONRV3XbuFuAP9/aEod0hLROOXoDMXLhcZvtavWmY7Rq6ca3g1jaed7WDqt7TotICwkOakF+Uc+0kzpu9v3UeT41+k6im0Ux8pXaX/MjJz7ILRV3b92PtjqU8PPQZ9mRsonvHO2rUzpN0a2N7nM6FtCzbvzmFYLZAaBC0bQadomyXFwzzkKB8VWXrkh1pn9CxdTxmixk/P9ug3p2+gT5x92I2m3j53fFMvu9Vmje5flxKRX/TNem/oa07RETqwsPiRMMS1y6JoKAQnlw6CD8/f27tcDsXC87yzpe2G6wcOL6NaYv68MQbd9KiSVu6d+jv0Efv2LtZ89UCXnhrFKVXihzmm8zlPLN8GMey9/G7N3/CD6e+tZt/9tIJhysjVFRXddt5isAAuK0j/LwvPPETeP7f4bl/h8eHwf1JtsuNeVpYhsrf0zc+tn1Y+sunv+WJN+7i93/9Gb8c+XKFfdzV8wFeeOt+Fn7wKyJCHY9ROZ93mmeWD+PE2TSeWT6MsxdPOLSxffi6rllES3rFDGbmkoEcPbOXO+NH1aidJ2rfwnaS55Sh8F8/g9/fD7NG2o5VvjPO88KykRPnDtKpVTzlpjLyCs+TW5BNRtZuurbrS+r+D0g//R1vrnuWp5cN4dCJb4CK//Zr0n9DXXeIiNSGR94a2xu483a2X+1fw3tbXmbqfQtI7JLsMD8r5wgvvzuewQmjad2sExFhzegde7dhn9v2f1itdlfp9rY1U9/jZfk/Z/HDqZ3MnrSO8NCmDvNvHBM/6TuRf2xfzH/85I+GfRYU5Var3VUaE7VT27GRuu99khPHVDq/pn/TtV1e6w4R8QUKzLXkzsDsCbTRqxlfGC8aE7XjC2PjRhonIuKNdAxzLTX20MtPuYqv//w15Qvvly/8jPXB1943X/t5RaRh0B5mEREREREDHnhalIiIiIiI51BgFhERERExoMAsIiIiImJAgVlERERExIACs4iIiIiIAQVmEREREREDCswiIiIiIgYUmEVEREREDCgwi4iIiIgYUGAWERERETGgwCwiIiIiYkCBWURERETEgAKziIiIiIgBBWYREREREQMKzCIiIiIiBhSYRUREREQMKDCLiIiIiBhQYBYRERERMfB/N4gV0XTEto4AAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 929.397x204.68 with 1 Axes>"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "new_qc_inverse = transpile(new_qc, backend)\n",
    "new_qc_inverse = new_qc_inverse.inverse()\n",
    "new_qc_inverse.draw('mpl')"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Similar procedure goes for $U^{\\dagger}$. Building $U$ circuit and inversing the circuit, we obtain $U^{\\dagger}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "def U_b(vector: np.array, to_gate:True):\n",
    "    \n",
    "    n = int(np.log2(len(vector)))\n",
    "    qc = QuantumCircuit(n)\n",
    "    qc.isometry(vector / np.linalg.norm(vector), list(range(n)), None)\n",
    "    qc = transpile(qc, backend)\n",
    "    \n",
    "    return qc.to_gate(label=\"$U_b$\") if to_gate else qc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "vector = np.ones(4)\n",
    "my_circ1 = U_b(vector,to_gate=False)\n",
    "my_circ1 = transpile(my_circ1, backend)\n",
    "my_circ1.draw('mpl')"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$U^{\\dagger}$ circuit is:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mycirc1_inv = my_circ1.inverse()\n",
    "mycirc1_inv.draw('mpl')"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Having dealt with the operators used to calculate the expectation valie, we combine them and write the main program. The explanation of step by by step code is given as follows:\n",
    "\n",
    "- First, we determine the number of qubits based on the vector.\n",
    "- The matrix is decomposed into a string of Pauli operatos along with the coefficients.\n",
    "- To begin, we select the ansatz input, which can be either a user-created (Hardware-efficient) version or one of the predefined versions found in the qiskit library.\n",
    "\n",
    "- The expectation values are then computed.\n",
    "- We compute the cost function by comparing the optimization routines requested by the user, such as SPSA versus SciPy, and dispatching to the appropriate optimizer. The chosen optimizer starts at x0 and calls 'cost_loc' every time the cost function needs to be evaluated. The outcome is a SciPy OptimizerResult object.\n",
    "- Finally, we have the result, the cost function, and the quantum solution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [],
   "source": [
    "# The entrypoint for our Runtime Program\n",
    "def main(backend, user_messenger,\n",
    "         matrix,\n",
    "         vector,\n",
    "         ansatz,\n",
    "         ansatz_config,\n",
    "         ansatz_circuit,\n",
    "         num_params,\n",
    "         x0,\n",
    "         optimizer,\n",
    "         optimizer_config,\n",
    "         shots,\n",
    "        ):\n",
    "    \n",
    "    \n",
    "    #determine number of qubits\n",
    "    num_qubits = int(np.log2(len(vector))) # number of system qubits\n",
    "    tot_qubits = num_qubits + 1  # Addition of an ancillary qubit\n",
    "    ancilla_idx = num_qubits  # Index of the ancillary qubit for the Hadamard test(last position)\n",
    "    \n",
    "    #Convert the matrix to Hamiltonian\n",
    "    pauli_op = SparsePauliOp.from_operator(matrix) \n",
    "    coeffs = pauli_op.coeffs\n",
    "    obs_list = pauli_op.paulis\n",
    "    \n",
    "    if ansatz=='custom':\n",
    "        \n",
    "        ansatz_circuit = ansatz_circuit\n",
    "        num_params = num_params\n",
    "        \n",
    "    else: \n",
    "            \n",
    "        # We grab the requested ansatz circuit class from the Qiskit circuit library\n",
    "        # n_local module and configure it using the number of qubits and options\n",
    "        # passed in the ansatz_config.\n",
    "        ansatz_instance = getattr(lib_local, ansatz)\n",
    "        ansatz_circuit = ansatz_instance(num_qubits, **ansatz_config)\n",
    "\n",
    "        # Get the number of parameters in the ansatz circuit.\n",
    "        num_params = ansatz_circuit.num_parameters\n",
    "    \n",
    "    # Use a given initial state, if any, or do random initial state.\n",
    "    if x0:\n",
    "        x0 = np.asarray(x0, dtype=float)\n",
    "        if x0.shape[0] != num_params:\n",
    "            raise ValueError('Number of params in x0 ({}) does not match number \\\n",
    "                              of ansatz parameters ({})'. format(x0.shape[0],\n",
    "                                                                 num_params))\n",
    "    else:\n",
    "        x0 = 2*np.pi*np.random.rand(num_params)  \n",
    "        #x0 = [float(np.random.randint(0,10))/1000 for i in range(0, num_params)]\n",
    "           \n",
    "    \n",
    "    # Here we define a callback function that will stream the optimizer parameter vector\n",
    "    # back to the user after each iteration.  This uses the `user_messenger` object.\n",
    "    # Here we convert to a list so that the return is user readable locally, but\n",
    "    # this is not required.\n",
    "        \n",
    "    def callback(xk):\n",
    "        user_messenger.publish(list(xk))\n",
    "\n",
    "    \n",
    "    def CA(qc, obs_list, idx):\n",
    "\n",
    "        operators = obs_list[idx]   \n",
    "\n",
    "        for j in reversed(range(len(operators))):\n",
    "            if operators[j] == Pauli('I'):\n",
    "                None\n",
    "            elif operators[j] == Pauli('X'):\n",
    "                qc.cx(ancilla_idx, j)\n",
    "            elif operators[j] == Pauli('Y'):\n",
    "                qc.cy(ancilla_idx, j)\n",
    "            elif operators[j] == Pauli('Z'):\n",
    "                qc.cz(ancilla_idx, j)\n",
    "            \n",
    "    def CA_dag(qc, obs_list, idx):\n",
    "\n",
    "        operators = obs_list[idx]   \n",
    "\n",
    "        for j in range(len(operators)):\n",
    "            if operators[j] == Pauli('I'):\n",
    "                None\n",
    "            elif operators[j] == Pauli('X'):\n",
    "                qc.cx(ancilla_idx, j)\n",
    "            elif operators[j] == Pauli('Y'):\n",
    "                qc.cy(ancilla_idx, j)\n",
    "            elif operators[j] == Pauli('Z'):\n",
    "                qc.cz(ancilla_idx, j)\n",
    "                \n",
    "    def local_hadamard_test(weights, l=None, lp=None, j=None, part=None):\n",
    "    \n",
    "        qc = QuantumCircuit(num_qubits+1, 1)\n",
    "\n",
    "        # First Hadamard gate applied to the ancillary qubit.\n",
    "        qc.h(ancilla_idx)\n",
    "\n",
    "        # For estimating the imaginary part of the coefficient \"mu\", we must add a \"-i\"\n",
    "        # phase gate.\n",
    "        if part == \"Im\" or part == \"im\":\n",
    "            qc.p(-np.pi / 2, ancilla_idx)\n",
    "\n",
    "        qubit_list = []\n",
    "\n",
    "        for s in range(num_qubits):\n",
    "            qubit_list.append(s)\n",
    "            \n",
    "        # Variational circuit generating a guess for the solution vector |x>\n",
    "        qc.append(ansatz_circuit.bind_parameters(weights), qubit_list)\n",
    "\n",
    "        qc.barrier()\n",
    "        # Controlled application of the unitary component A_l of the problem matrix A.\n",
    "        CA(qc, obs_list, l)\n",
    "\n",
    "        nb = int(np.log2(len(vector)))\n",
    "        vector_circuit = QuantumCircuit(nb)\n",
    "        vector_circuit.isometry(vector / np.linalg.norm(vector), list(range(nb)), None)\n",
    "        vector_circuit = vector_circuit.inverse()\n",
    "\n",
    "        qc.append(vector_circuit, qubit_list)\n",
    "\n",
    "        qc.barrier()\n",
    "\n",
    "        # Controlled Z operator at position j. If j = -1, apply the identity.\n",
    "        if j != -1:\n",
    "            qc.cz(ancilla_idx, j)\n",
    "\n",
    "        qc.barrier()\n",
    "\n",
    "        # Unitary U_b associated to the problem vector |b>.\n",
    "        qc.isometry(vector / np.linalg.norm(vector), list(range(nb)), None)\n",
    "\n",
    "        # Controlled application of Adjoint(A_lp).\n",
    "        CA_dag(qc, obs_list, lp)\n",
    "        #CA(qc, obs_list, lp)\n",
    "\n",
    "        qc.barrier()\n",
    "        # # Second Hadamard gate applied to the ancillary qubit.\n",
    "        qc.h(ancilla_idx)\n",
    "\n",
    "        qc = transpile(qc, basis_gates=['id', 'rz', 'sx', 'x', 'cx'])\n",
    "        qc.barrier()\n",
    "        qc.measure(ancilla_idx, 0)\n",
    "\n",
    "        job = backend.run(qc, shots=shots)\n",
    "        result = job.result().get_counts()\n",
    "\n",
    "\n",
    "        if len(result)==1:\n",
    "            for s in result.keys():\n",
    "                if s=='0':\n",
    "                    expectation = result['0']/shots\n",
    "                if s=='1':\n",
    "                    expectation = -result['1']/shots\n",
    "        else:        \n",
    "            expectation = (result['0']-result['1'])/shots\n",
    "\n",
    "        return expectation\n",
    "    \n",
    "    def mu(weights, l=None, lp=None, j=None):\n",
    "        \"\"\"Generates the coefficients to compute the \"local\" cost function C_L.\"\"\"\n",
    "\n",
    "        mu_real = local_hadamard_test(weights, l=l, lp=lp, j=j, part=\"Re\")\n",
    "        mu_imag = local_hadamard_test(weights, l=l, lp=lp, j=j, part=\"Im\")\n",
    "        \n",
    "        return mu_real + 1.0j * mu_imag\n",
    "    \n",
    "    def psi_norm(weights):\n",
    "        \"\"\"Returns the normalization constant <psi|psi>, where |psi> = A |x>.\"\"\"\n",
    "        norm = 0.0\n",
    "\n",
    "        for l in range(0, len(coeffs)):\n",
    "            for lp in range(0, len(coeffs)):\n",
    "                norm = norm + coeffs[l] * np.conj(coeffs[lp]) * mu(weights, l, lp, -1)\n",
    "\n",
    "        return abs(norm)\n",
    "\n",
    "    \n",
    "    def cost_loc(weights):\n",
    "        \"\"\"Local version of the cost function. Tends to zero when A|x> is proportional to |b>.\"\"\"\n",
    "        mu_sum = 0.0\n",
    "\n",
    "        for l in range(0, len(coeffs)):\n",
    "            for lp in range(0, len(coeffs)):\n",
    "                for j in range(0, num_qubits):\n",
    "                    mu_sum = mu_sum + coeffs[l] * np.conj(coeffs[lp]) * mu(weights, l, lp, j)\n",
    "\n",
    "        mu_sum = abs(mu_sum)\n",
    "        pn = psi_norm(weights)\n",
    "\n",
    "        # Cost function C_L\n",
    "        cf = 0.5 - 0.5 * mu_sum / (num_qubits * pn)\n",
    "\n",
    "        cost_function.append(cf)\n",
    "\n",
    "        return cf\n",
    "    \n",
    "    # Here is where we actually perform the computation.  We begin by seeing what\n",
    "    # optimization routine the user has requested, eg. SPSA verses SciPy ones,\n",
    "    # and dispatch to the correct optimizer.  The selected optimizer starts at\n",
    "    # x0 and calls 'cost_loc' everytime the optimizer needs to evaluate the cost\n",
    "    # function.  The result is returned as a SciPy OptimizerResult object.\n",
    "    # Additionally, after every iteration, we use the 'callback' function to\n",
    "    # publish the interm results back to the user. This is important to do\n",
    "    # so that if the Program terminates unexpectedly, the user can start where they\n",
    "    # left off.\n",
    "    \n",
    "    cost_function = []\n",
    "    \n",
    "    # Since SPSA is not in SciPy need if statement\n",
    "    if optimizer == 'SPSA':\n",
    "        opt = SPSA(maxiter=500)\n",
    "        res = opt.minimize(cost_loc,x0)\n",
    "    # All other SciPy optimizers here\n",
    "    else:\n",
    "        res = minimize(cost_loc, x0, method=optimizer,options=optimizer_config, callback=callback)\n",
    "    # Return result. res is a subclass of dict.\n",
    "    \n",
    "    circ = ansatz_circuit.bind_parameters(res.x)\n",
    "    circ.measure_all()\n",
    "    \n",
    "    t_circ = transpile(circ, backend)\n",
    "    job = backend.run(t_circ, shots=shots)\n",
    "    result = job.result()\n",
    "    \n",
    "    sorted_keys = sorted(result.get_counts())\n",
    "    solution_probabilities = np.zeros(2**num_qubits)\n",
    "    i=0\n",
    "    for key in sorted_keys:\n",
    "        solution_probabilities[i] = result.get_counts()[key]\n",
    "        i=i+1\n",
    "\n",
    "    o = solution_probabilities/shots\n",
    "    \n",
    "    # quantum probabilities to actual classical solution\n",
    "    quantum_solution = np.sqrt(np.real(o))*np.linalg.norm(vector)\n",
    "    \n",
    "    \n",
    "    return res, cost_function, quantum_solution,"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Giving the inputs for backend, matrx and vector."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [],
   "source": [
    "backend = Aer.get_backend('qasm_simulator')\n",
    "msg = UserMessenger()\n",
    "\n",
    "matrix = np.eye(2)\n",
    "vector = np.ones(2)\n",
    "\n",
    "A_inv = np.linalg.inv(matrix)\n",
    "classical_solution = np.dot(A_inv, vector)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Executing the main program to obtain the result, cost function and quantum solution"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "c:\\Users\\mummanen\\Anaconda3\\lib\\site-packages\\scipy\\optimize\\_minimize.py:548: RuntimeWarning: Method COBYLA does not support callback.\n",
      "  warn('Method %s does not support callback.' % method, RuntimeWarning)\n"
     ]
    }
   ],
   "source": [
    "# Execute the main routine for our sparseMatrix and vector\n",
    "# Try another qiskit ansatz: EfficientSU2, TwoLocal, NLocal, PauliTwoDesign, UCCSD, PUCCD, SUCCD\n",
    "# Or another scipy optimizer: SPSA, Nelder-Mead, Powell, CG, BFGS, Newton-CG, L-BFGS-B, TNC, COBYLA, SLSQP etc.\n",
    "\n",
    "\n",
    "result, cost_fn, quantum_solution = main(backend,msg,\n",
    "                                         matrix,\n",
    "                                         vector,\n",
    "                                         ansatz='TwoLocal',\n",
    "                                         ansatz_config={'reps': 3, 'entanglement': 'linear', 'rotation_blocks': ['ry', 'rz'], 'entanglement_blocks': 'cx'},\n",
    "                                         ansatz_circuit = None,\n",
    "                                         num_params = None,\n",
    "                                         x0=None,\n",
    "                                         optimizer='COBYLA',\n",
    "                                         optimizer_config={'maxiter': 200},\n",
    "                                         shots = 1024,\n",
    "                                         )\n"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Visualization of the results (both cost function and quantum solution)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAfUAAAFXCAYAAAC7nNf0AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAABUsklEQVR4nO3deXxU1fk/8M/MnSXLZIOEQFgCBCJYVDZRxCAuaBXXAgIqaMWFblq11qUuiDRgka8bamltraJVFBCLCvJDUJRNQSJE9n0JhJCFzL7e3x+Te2cmmT2ZDLn5vF8vXiVz5945RyjPPOc85xyVKIoiiIiIqM1TJ7sBRERE1DIY1ImIiBSCQZ2IiEghGNSJiIgUgkGdiIhIIRjUiYiIFEKT7AY0V1WVsUWfl5OThtpaS4s+82zEfioL+6ks7KeytHQ/8/IyQl5jpt6IRiMkuwmtgv1UFvZTWdhPZWnNfjKoExERKQSDOhERkUIwqBMRESkEgzoREZFCMKgTEREpBIM6ERGRQjCoExERKQSDOhERkUIkLKh7PB4888wzmDBhAiZPnozDhw8HXP/yyy8xduxYjBs3Dh9//HFU9xAREVFoCdsmdtWqVXA4HFi4cCHKysowe/ZsvPnmmwAAt9uNuXPnYvHixUhLS8N1112HK6+8Eps3bw55DxEREYWXsEx9y5YtKCkpAQAMHDgQ5eXl8jVBEPDFF18gIyMDdXV1AID09PSw91DsjhqP4P0d70IUxWQ3hYiIWkHCMnWTyQSDwSD/LAgCXC4XNBrvR2o0GqxcuRIzZszAZZddBo1GE/GeYHJy0lp8X91wm+W3JTN++Dvm/TAP1w0YjeKOxU2uK6WfkbCfysJ+Kgv72bISFtQNBgPMZrP8s8fjaRKcr776alx11VV4/PHHsXTp0qjuaaylT/jJy8to8ZPfkuVE3SkAwJGTJ5Hj6RJwTUn9DIf9VBb2U1nYz/ifF0rCht8HDx6MtWvXAgDKyspQXOzLFE0mE+644w44HA6o1WqkpqZCrVaHvYdiZ3Z5vyDZ3Y4kt4SIiFpDwjL10aNHY926dZg4cSJEUURpaSmWLVsGi8WCCRMm4IYbbsDtt98OjUaDc845BzfeeCNUKlWTeyh+ZqcU1G1JbgkREbWGhAV1tVqNGTNmBLxWVFQk/37ChAmYMGFCk/sa39Oavjq8EqcOHMek3r9OWhtaksVpAgA43PYkt4SIiFoDN5/x83rZq3joy4cUUy3uy9Q5/E5E1B4wqPtRqwR4RA9cHleym9IipKDOTJ2IqH1gUPejVXtnIxweZWS2Joe32tLOoE5E1C4wqPvRCjoAgFMhw9Vypq6QLylERBQeg7ofnbohqCtg+N3hdsjB3O5i9TsRUXvAoO5Hq9YCAJwKyGwtTt8mPiyUIyJqHxjU/WgFb1B3KCAImv2COgvliIjaBwZ1P1p5+N2Z5JYAHtEDt8cd9/3+Qd3uYVAnImoPGNT96ARp+D35Qf3yhZfg91/dH/f95oaNZwBljDwQEVFkCdtRri2SM/UkB0GP6MHOmh04XH8ITrdTnhaIhbTvO8BCOSKi9oKZuh+pUC7ZS8CkdeUWlwXlp7fF9QwzC+WIiNodBnU/Ukac7B3lbC6r/PuNJzbE9Qxp4xkAcHBOnYioXWBQ9yOtU0/2HLTNb7h844n1cT2D1e9ERO0Pg7ofzVmyTt3q9mXq35/YENcBMxx+JyJqfxjU/ejkderJrX73z9SrbdXYV7c35mcEVr8zUyciag8Y1P1I1e+uJC9pk6rVM3VZAOIbgvfP1G0M6kRE7QKDuh9dw4Euya5+t7m9Qf3SriMBAJviKJZjpk5E1P4wqPuR935P8vC7taH6/YK8gcjWZ8dVAc9COSKi9odB3c/Zsk5dmlNP1aZiWOeLcaT+EE6YKmJ6BgvliIjaHwZ1P9Lwe7Ln1G0N1e8pQiouKrgEQOxD8NLwe5Y+m5k6EVE7waDuRyuvU092oZw3CKdoUnBxl+EAYi+WMzlN0Kg1yNBmJH3dPRERtQ4GdT9atXcr/LNlnXqKkIIL8gYhRUjBphMbY3qGxWmGQWuATtDJhXdERKRsDOp+tMLZtaNciiYVOkGHwflDsaO6HGfsdVE/w+w0I11rgF7Qc/idiKidYFD3oztL1qlLe7+naFIAABd1uRgiRGw++X3UzzA7TUjXpkMn6JP+JYWIiFoHg7of6UAXR7I3n3H7MnUAuKiLt1gulqVt3kw9HTpBJ5/6RkREysbz1P341qkneU5dWtImeDP1CzsPg1qljrpYzuVxwea2IV1rgAoquEU33B43BLWQsDYTEVHyMVP3I1e/J32dunf4Xd8w/J6hy8QvOp6HrZVbAvaFD8XSsEZdytQBMFsnImoHGNT96M6W89Sl4feGTB3wzqs7PA6UVW2NeL85IKjrAfiG9ImISLkY1P1oz5rz1L2ZemrDnDoAXNwwr/59FPPqvqBugP4sqegnIqLEY1D3oz1LzlOXTlWTqt8B4CJpE5qKyPPqJqcRAJAWkKlz+J2ISOkY1P341qmfLUvafJl6fnpn9Mzshe9PboLb4w57v5SpGxrWqQPM1ImI2gMGdT86tTSnnuygboMKKnndvOTigktQ7ziDXTU7w94v7fue7hfUmakTESkfg7of7dlynrrLihRNClQqVcDrF3WObh94FsoREbVPDOp+zpbz1G1uW0Dlu2Rgp8EAEEWm7gvqLJQjImo/GNT9CCoBKqjgPAuG3/3n0yW5qbkAgDP22rD3+w+/s1COiKj9YFD3o1KpoBN0Z0H1u02eC/eXpc8GANTaIgX1psPvPNSFiEj5GNQb0Qm6s6L6PTVIpp6iSUGqJhV1EU5r8w/qKXKmzuF3IiKlS9je7x6PB9OnT8fu3buh0+kwc+ZMFBYWytc/++wzvPPOOxAEAcXFxZg+fTrUajVuvvlmZGRkAAC6deuGWbNmJaqJQWkFbbMydbvbjg93vY8J59wWsM48Ft7h9+D3ZumzURdh+F1ap+4//O7wMFMnIlK6hAX1VatWweFwYOHChSgrK8Ps2bPx5ptvAgBsNhtefvllLFu2DKmpqXj44YexZs0aXHrppQCABQsWJKpZEXmH3+PP1D/b/yke/eaP0At6TOx3e8z3i6LoLZQLkqkDQI4+BxXmirDPCLZO3e5iUCciUrqEDb9v2bIFJSUlAICBAweivLxcvqbT6fDhhx8iNdUbuFwuF/R6PXbt2gWr1Yq7774bU6ZMQVlZWaKaF5JO0DWr+r3aehoAcMpyKq77pYK2YNXvgDdTr7efCbsBjf82sTzQhYio/UhYpm4ymWAwGOSfBUGAy+WCRqOBWq1Gbq63knvBggWwWCwYMWIE9uzZg6lTp2L8+PE4dOgQ7r33XqxYsQIaTehm5uSkQaNpuSNFtWotrB4r8vIy4rpf1Hm/EDjU5rieUWv1HiaTmWYIen9+Zh7EEyJ0mR50SM0O+gynyrsmvWdBZ+wwe9+jS1U1eV68fWxr2E9lYT+Vhf1sWQkL6gaDAWazWf7Z4/EEBGePx4M5c+bg4MGDeO2116BSqdCrVy8UFhbKv8/OzkZVVRW6dOkS8nNqay0t2m6doEOdtQ5VVca47j9Z683Uj9eejOsZJ81VAAC1WxP0/jSV9y/GvuNH0StLG/QZteYzEFQC6mscsJm8GX11fX3A8/LyMuLuY1vCfioL+6ks7Gf8zwslYcPvgwcPxtq1awEAZWVlKC4uDrj+zDPPwG6344033pCH4RctWoTZs2cDACorK2EymZCXl5eoJgblnVOP/+hVU8Ma8UjLzkKxBtn33Z+0rK0uzPPNTjPStQaoVCr5THYuaSMiUr6EZeqjR4/GunXrMHHiRIiiiNLSUixbtgwWiwUDBgzAokWLMHToUNx5550AgClTpmDcuHF44oknMGnSJKhUKpSWloYdek+E5q5TNzrqAQC19pq47rcHOaHNX44+BwDCLmszO01I16YDANepExG1IwmLmGq1GjNmzAh4raioSP79rl27gt43d+7cRDUpKlpB26wtVaXd3MJl0uFIJ7TpQxXKpWR7nx9mWZvZaUZ2Q0avbzgUxsbqdyIixePmM43oBB3cohse0RPX/UaHd96kxhZfpm5zeYvcUiNk6rVhgrqlYfgdANepExG1IwzqjUhLwOJdqy7NqdfZayGKYsz3y3PqQvg59TO2uqDX3R43LC6LPPzO89SJiNoPBvVG5KAeZxA0NWTqdrcdFlfslfm2hiNSw20+A4TO1C0u38YzALhOnYioHWFQb0Q6fjXeM9Wl4Xcgvnl1u0sK6uHn1M+EKJTz3/cdAKvfiYjaEQb1RnyZenzD71KhHADUxFEBL2fqIQrlImXq/seuAn6FcgzqRESKx6DeSHPm1F0eV8CQezyZum+devCgnqnLAhB9ps4lbURE7QeDeiNSUI9n+N0/SweA2jgq4G2u8HPqglrwntQW4gtDk+F3FsoREbUbDOqNSHPq8Qy/S/PpGrV3+X+4ZWeh2OTq99DHtmbrs0NuPiN9sUhrGH4X1AIElcBCOSKidoBBvZHmZOrScrYCQzcA8WXqdre0Tj14pg4A2fqckJvPNM7UAW+2zkydiEj5GNQbkYK6K65M3btFbI+MHgDi2//d2jD8Lg2bB5Otz4bVZZWH6v2FDurM1ImIlI5BvRFfph57UDc5vJl694agHm4r11Bs7vAHugDeTB0IXiwnrZM3aH2n+OgEvVxVT0REysWg3ohWaJhTj2v43RtQu8uZenMK5cLMqaeEXtbG4XciovaLQb0ROVOPIwhKmXqBoSvUKnVchXLRzalnAwh+UpsvqBvk13SCjoVyRETtAIN6I/KcehzD79KceoYuEzn6nLgy9ejm1BuOXw0yZ292SZvP+DJ1HTN1IqJ2gUG9EV+mHsecekP1e4YuA9kpOXEVyslL2qLK1KMbfk9hoRwRUbvAoN6IvE69GUvaDFoDsvU5qLXXxHxSW6RtYgHfnHr4oO4//K6H3W2P69Q4IiJqOxjUG2nONrHS5jMGXQY6pHSAy+NqsstcJDaXFSlCClQqVcj3hJ9TDz78LkKM+zhZIiJqGxjUG2nOgS7ScrIMbYY8710T47y6zWUPW/kORJhTd5qhgiqg0E4vTylwCJ6ISMkY1Btp3o5yUqZuQIeUDgBiX6tuc1uhDzP0DgA58vB7XZNrZqcZ6VpDQKYvHepiZ7EcEZGiMag3Iq9Tb8aSNoM2Q573jj1Tt0XM1LPCFMqZHEYYdIaA15ipExG1Dwzqjfjm1F0x32t0GpGmSYOgFpAjZeoxVsDbXNawa9QBIE2TBp1aF3L43X8+HYCc+XOtOhGRsjGoN+IL6vFk6kYYdN7tWXOkOXV7jJm62xa28h0AVCqV9/jVMMPv/nQ8fpWIqF1gUG9EWtIWTwA0OowwNATUeDJ1URRhc9mgjzD87n1+05PaPKIHFlewTN37RcXO/d+JiBSNQb2R5ixpMztNyNBlAvBl6rFsFevwOCBCjJipA5AzdY/okV+zuCwAgHRNYFD3Fcpx+J2ISMkY1BuJN6i7PC5YXJYmmXosW8VKu8lFmlMHvF8aPKJHXkYHBN94BmjefvZERNR2MKg34lunHlsANPttEQv4LTuLYfg9mhPaJFlBNqAJtvEMAOjVzNSJiNoDBvVGpCVtsZ6nLu0mJ2XJBm0GNGpNTEva5C1io8nUg2wVG2zfdwDyHH08a++JiKjtYFBvJN5M3dQoU1epVMjWNy1mC8cmn9AWb6YefPhdLpRzsVCOiEjJGNQbiXdO3eS377ukQ0qHuObUoxl+zwmyVazZ4dvRzh8L5YiI2gcG9UbiXacuH+bilyV7T2qrDahQD8faMPyeKkQefs8OslVsqOF3nZqFckRE7QGDeiO+deoxZuqNht8B77y3R/TA6KiP6hmxZOrBzlQPPfzekKl7mKkTESkZg3oj8Wbq8vC71j+oS8vaoptXl4bHoymUk05q8392qOp3345yDOpERErGoN5I3HPqzqZz6r7AG928upSpS5l1ONKzz0Qx/J6i4TaxRETtAYN6I75T2uJb0uY/py4dvxrtrnLWGDafkebUa4MMv6eF2Pvdxup3IiJFY1BvRK1SQ1AJMa/pDjanHnumHsPmM7osAI0z9RDD72pm6kRE7QGDehA6QRfzOnVjiCVtQPBzz4OxuRsK5aKoftcKWhi0GY3m1ENsPiOtU2ehHBGRojGoB6FV62I+T10qlMvwK5STh8ijLZRzSYVykefUAW8FvH+mLo0W+BfrASyUIyJqLzSJerDH48H06dOxe/du6HQ6zJw5E4WFhfL1zz77DO+88w4EQUBxcTGmT58OAGHvaS06QRt79buz6cYvsR7qYo0hUwe8XxoOnjkg/xxy73cWyhERtQsJy9RXrVoFh8OBhQsX4pFHHsHs2bPlazabDS+//DLeffddfPjhhzCZTFizZk3Ye1qTRq2NOQCaHE2z5FiPX41lTh3wZupmp0ku6jM7zVBB1aTQjge6EBG1DwkL6lu2bEFJSQkAYODAgSgvL5ev6XQ6fPjhh0hN9QYfl8sFvV4f9p7WpFPr4Ipx+N3oNCJNkwZBLcivxZqp+zafiTJT1wfuKmd2mpGmTYdaFfjHym1iiYjah4QNv5tMJhgMvqFoQRDgcrmg0WigVquRm5sLAFiwYAEsFgtGjBiB5cuXh7wnlJycNGg0Qsjr8UjR6VFvr0deXkbkNzewus3I0GcE3JMrGqATdDC5o3yWxg0A6NopF3kdIr+/S3YnAIA63YG83AzYPBZk6A1NPktna9imVnAHXIulf20Z+6ks7KeysJ8tK2FB3WAwwGw2yz97PJ6A4OzxeDBnzhwcPHgQr732GlQqVcR7gqmttbRou/PyMiCIGtiddlRVGaO+74ytHgatock9OfoOOGWqiupZdWbve8xnXKhyR36/XvTOnR84cQwdxa6otxmRpk1r8lk2l3cqwWg1y9fy8jJi6l9bxX4qC/upLOxn/M8LJWHD74MHD8batWsBAGVlZSguLg64/swzz8But+ONN96Qh+Ej3dNaNGptzOepmxzGgOVskpyUnICT1MKJfU498KQ2s9PcZN93wLdLHgvliIiULWGZ+ujRo7Fu3TpMnDgRoiiitLQUy5Ytg8ViwYABA7Bo0SIMHToUd955JwBgypQpQe9JBp2ghSuGoO72uGFxWQKWs0my9TnYXbMLbo87YL49mFjWqXufnQ3AO6cuiiLMTlOTynfAu6GOTq3jnDoRkcIlLKir1WrMmDEj4LWioiL597t27Qp6X+N7kkGr1sHhdkAURahUqojvD7acTZKT0gEiRNQ7zsiFc6HEmqnnyMev1sLqskKEGDSoA95iOWbqRETKxs1ngtCqtRAhwi26o3p/sOVskpwYtoq1u23QqXVNqtdDyfLL1ENtPCPRCzrY3dz7nYhIyRjUg5AOdYk2szUGOaFNkhPDoS5Wly3q5WyA7wtDna025MYzEp2g5/A7EZHCMagHoVN7C8uinVeXt4gNUSgHRJep21zWqIfegcBMPdS+7xKdoOPwOxGRwjGoB6GVqsWjDOrBjl2V+DagiZyp22LN1P3m1H1BvWkbAO8Z7dz7nYhI2RjUg9CqvfWD0Z7UFuzYVUksx6/a3FakCtFn6gZtBgSV0JCphx9+1wspsDNTJyJSNAb1ILRqKVOPMqjLmXrToN4hhjl1m8sOfQzD7yqVCtn67IY59cjD7yyUIyJSNgb1IKTNWqKeUw9TKOfbnz18UBdFETa3FSkxZOqA96S2Wrt/oVzo4XenxwmP6Inp+URE1HYwqAehVUvV782fU+8Q5aEuUsCNZU4d8J2pHrn6nbvKEREpHYN6EFJQj/ZMdXmNeJDNZ7Ll6vfwmbp0QltqDMPvgHckwOlxospyytuGMJk6ABbLEREpGIN6ENoYs1rfkrbMJtdSNalI1aRGzNStDfPd0W4RK5GWtR0zHQMQbvjd+2WBxXJERMrFoB6EriFTj/ZM9XDD74A3m45UKGdv2CJWr9FH20wAvmVtx41SUI80/M5MnYhIqSLu/X78+HG89957OHPmDERRlF+fNWtWQhuWTL516tFlteYwS9oA71r1Y8ajYZ8h7/seZ6Z+3BQ+qEvD76yAJyJSrohB/Y9//COGDh2KoUOHRnW4iRLIc+rRbhPbkKmHGvrO0edgR3U5XB4XNOrg/8mlE9pinVOXtoqtMB0P2wYpU+fwOxGRckUM6i6XC4899lhrtOWs4VunHu2SNhPSNGkhA3a2vPNbHXJTc4O+xyqf0BZfpi6NKoTb+x3g8DsRkZJFnFMfMmQIVq9eDYej/WR4OkGaU492SVt9yAwZ8C1rqwtTAS8Ni0vD5NFqfJxrqiYt6Pv0amn4vf38ORIRtTcRM/UVK1bgvffeC3hNpVJh586dCWtUssmZegzbxIaaTwd8G9DUhKmAl5a0xb5OPUf+fZomDYJaCPo+qQDP4WGmTkSkVBGD+nfffdca7Tir+NapR39KW+f0LiGvS9l0nT1cUPdm6rGvU8+Wf58WYugd4Dp1IqL2IGJQt1qtmDdvHjZs2AC3242LL74YDz74INLSgg/zKoF0nno0Qd3tccPisiAjyL7vkpwoMnWrlKnHWP0uzdcDoZfUAb5COZuLQZ2ISKkizqnPmDEDVqsVpaWleOGFF+B0OvHss8+2RtuSRhfD8Ltv3/fQAdWXqYeeU7dJm880I1MPN6/PQjkiIuWLmKn//PPP+N///if//Mwzz+C6665LaKOSTVqnHs02sSZH+INUAL9zz8MVykmbz8R4oIte0CNNkwaLyxKy8h3wL5RjUCciUqqImbooiqivr5d/rq+vhyAEL8ZSCuk89WgOdDE6Q28RK5Ey9bCFcu745tQB37K2cEFdF+OGOkRE1PZEzNTvuusujBs3DldccQVEUcSaNWtw3333tUbbkkaqfo9mSZspwhaxgG9OPdzwuzXO6nfAWwF/wlwRdrRAGtbn8DsRkXJFDOpjx47Feeedhx9++AEejwevvfYazjnnnNZoW9LEktWaImwRC/iK2WrCDL/L28TGkalnp2QDiJSp6xs+h0GdiEipQg6/r1mzBgCwdOlS7NixA+np6cjIyMDOnTuxdOnS1mpfUvi2iY1i+D2KTN07750eoVDO2vDeOIJ6w0hA2KCu5oEuRERKFzJT3759Oy6//HJs2rQp6PWbb745UW1KOimoR5Op+w5zCT2nDnh3lQt3/Kq9IYOOZ049W55TD//FAgDs3HyGiEixQgb1Bx54AABw/fXXY8SIEQHXVq5cmdhWJZlU/R7NnLrR4S0iDLekDfAOwR86czDk9XjXqQNRZuoxnhFPRERtT8ig/sUXX8DhcODVV1+VAzzgPeBl/vz5uPrqq1ulgcngG6qOplDOm6kbwmw+A3iL5cqd2+BwO+QA68+3Tj32oC4tmYs0BQBw+J2ISMlCBnWz2Ywff/wRZrM5YAheEAQ89NBDrdK4ZJFOW4tmnboxis1nAP8NaOrQKa1Tk+u+vd+bM6ceLqh7n8sDXYiIlCtkUB8/fjzGjx+PDRs2oLi4GB07doTVasWpU6dQWFjYmm1sdbEMVUtL2iLNqUuBt9ZWEzSoS5vCpMRRKHd90U34qWorru09JuR7fOep22J+PhERtQ0RN5/Zu3cv7rnnHgBATU0Npk2bhoULFya8YcnkW6fuivjeaKrfASAvLQ8AcMpSGfS61WWFVq0NecpaOLmpuXjp8nnokNIx5Hs4/E5EpHwRg/pHH32E999/HwDQtWtXLFmypMlRrEojnaceW/V7+Dn1/LTOAIBKy8mg120uW1zz6dGS1qlz+J2ISLkiBnWn0wmdzlfYpdVqE9qgs4FGXqcexZx6Q6Yebj4bgHw0a6U5eKZuc1vjGnqPll7gOnUiIqWLuKPcVVddhTvvvBPXXnstVCoVvvzyS1xxxRWt0bakkavfo9km1mlCqiZVLq4LJT8tH0D4TD01gZk6C+WIiJQvYlB/9NFHsWLFCvzwww/QaDSYMmUKrrrqqtZoW9JIATradeqRlrMBQH66d/j9VIigbnfbkKPvEEMrY6NRa6CCipk6EZGCRQzqAFBUVITc3FyIoggA+OGHH3DhhRcmtGHJpFKpoFProjxP3RRxORsA5KV2ggoqnDQHD+pWlw1d0hOXqatUKugFPavfiYgULGJQf+6557BmzRp0795dfk2lUuHdd99NaMOSTaPWwhnlKW3SfHk4WkGLjqkdwwy/W+Naox4LnaDn8DsRkYJFDOrr1q3DihUrkJKS2IBzttEJ2oiFcm6PGxaXJeJyNkmntM44Un+4yetOtxNu0Z3Q6nfAu1adw+9ERMoVsfq9e/fu8rB7LDweD5555hlMmDABkydPxuHDTYOZ1WrFxIkTsX//fvm1m2++GZMnT8bkyZPxxBNPxPy5LUWr1kXM1KNdzibJT8uHyWmE2WkOeF06oS2lYdlZougFPfd+JyJSsIiZelZWFsaMGYNBgwYFLG2bNWtW2PtWrVoFh8OBhQsXoqysDLNnz8abb74pX9++fTueffZZVFb6lnjZ7d4scsGCBTF3pKXpBF3E6vdol7NJ5GVtlpPonVUkvy6dcZ7oTF0v6Jt8oSAiIuWIGNRLSkpQUlIS84O3bNki3zdw4ECUl5cHXHc4HHj99dfx5z//WX5t165dsFqtuPvuu+FyufDwww9j4MCBMX92S9CoNRGzWlOUx65KpA1oTpkrA4O6nKkndopDL+jDHv9KRERtW8SgftFFF8X1YJPJBIPBl8EKggCXywWNxvuRQ4YMaXJPSkoKpk6divHjx+PQoUO49957sWLFCvmeYHJy0qDRxL61ajh5eRlI1aXAYjEjLy/00Pp+uxsA0CmrQ9j3SYryvXvmWzVnAt5frfK2P8eQGdVz4pWmT4XdaJc/I5GfdTZhP5WF/VQW9rNlRQzqd9xxB1QqFURRhMvlwunTp9G/f38sXrw47H0GgwFms2+o1+PxhA3OANCrVy8UFhZCpVKhV69eyM7ORlVVFbp0CV1dXltridSFmOTlZaCqygi1qIHD7URVlTHke4+e8k4dCC592PdJ0jzZAIA9Jw6gqpPv/RVVp72/cQlRPSdealEDu8uOqiqj3E+lYz+Vhf1UFvYz/ueFEjGor169OuDnbdu2yXvBhzN48GCsWbMG1113HcrKylBcXBzxnkWLFmHPnj2YPn06KisrYTKZkJeXF/G+RNCpI1e/R3uYi8S3/3vgVrFWl3ftuL4Vht/dohtujzuhn0NERMkRsfq9sfPPPx8///xzxPeNHj0aOp0OEydOxKxZs/DEE09g2bJlYU94GzduHIxGIyZNmoSHHnoIpaWlEbP7RNGotREPdJGq3w3RVr+nN2wV22gDGmlDmMSvU5eOX+WyNiIiJYoYMefNmxfw8969e9GxY+gjPiVqtRozZswIeK2oqKjJ+/wr3XU6HebOnRvx2a1BJ+jg8rggiiJUKlXQ9xgd9QBiWdIWPFO3uRoK5RJe/e790sC16kREyhQyqFssFqSlpTV5fdiwYRgzZkxCG3U20EontXmccobbmMnhzdSjXdKWoklBtj67yf7vtoZMPTXhw+/M1ImIlCxkUL/99tvxySef4PTp05g+fXorNunsIAVyh8cRMqgbnd459WgzdcCbrTcefrc2ZOr6VtgmFmBQJyJSqpBB3Wq14k9/+hO+/fZbeVMYf5E2n2nrAs5UD3GEvEkulIs+qHdK74zdtbtgd9uhbxRkW2OdOgDuKkdEpFAhg/rbb7+NTZs2YcuWLRg2bFhrtumsIJ2p7vS4Qr5Hqn6PLVP3Fcv1yPSuW2+tOXUWyhERKVvIoN6lSxfcfPPN6NevH/r169eabToraAW/TD0Eufo9yjl1wL9YzhfUpSVtqYkefldLmTqDOhGREkVc0tYeAzrgy9TDLWszxbikDfBf1uargG+tTD1F0zDcH2GpHhERtU0xr1NvL3xz6qEPdTE6jEjVpEKjjn4tvX+mLpGq3/UJPqVNLpRrGBkgIiJliRjU161b1+S1lStXJqQxZxOd4FvSForJaYypSA7wO9TFL6hLQTbxc+ocficiUrKQKeYXX3wBh8OBV199FQ888ID8utPpxD/+8Q9cffXVrdLAZNHKhXKhh6qNDiMMuujn0wGgc3pDpu4//C6vU0/w5jNqqVCOw+9EREoUMqibzWb8+OOPMJvN2LRpk/y6IAh46KGHWqVxySRtPuMIM/xucpjkM9Kj1akhqJ+0nJBfs8pz6q2zTt3hYaZORKREIYP6+PHjMX78eGzYsAHDhw+XX298pKpSydXvITJ1t8cNi8scU+U74K2UT9caGhXKtc7wu/SlgevUiYiUKeKcutVqxZw5c2A2m3HttdfiyiuvxJIlS1qjbUnlW6cePFOPZzmbJD8tP7BQTsrUE14ox3XqRERKFjGov/7667jhhhvwxRdf4Pzzz8fq1avx3nvvtUbbkkrbEABDrVOXj12NYTmbJD+9M6qtp+Fq2NhG3lEu0YVyala/ExEpWVRL2vr164evv/4aV1xxBdLT0+F0hp5nVgptwzI1R4hMXV6jHmP1O+DN1EWIqLKcAuCdU9eoNTEtjYuHfKAL16kTESlSxKCem5uL559/Htu3b0dJSQlmz56NgoKC1mhbUkXO1GM7dtVf47XqNrcNKQmufAe4pI2ISOkiBvW5c+fivPPOw3vvvYe0tDR07979rDnzPJEizan7dpOLY069oWJeOlfd5rLKu70lkp5BnYhI0SKO96anp8NsNuPFF1+Ey+XCRRddFPScdaXxP089GPkwlziH3wHgpNm7rM3utrdKpu47FY7D70REShQxqP/tb3/D4cOHMXbsWIiiiCVLluDo0aN46qmnWqN9SSMtaQu1/Mscx77vknx5Axrv8LvVZUW2PjuOVsbGd546C+WIiJQoYlBft24dli5dCrXaO1I/atQo3HDDDQlvWLJF2lHOFMexqxLfnLo0/G5DSnrrZepcp05EpEwR59TdbjdcLlfAz4IgJLRRZwNpTXeo89Sl4ff0ONepA779321uK1KExO4mB3CdOhGR0kXM1G+44QZMmTIFY8aMAQB8/vnnuP766xPesGST59RDZLVSoVw8mXqWPht6QY9K80m4PC64PK6EbxELsFCOiEjpIgb1adOm4dxzz8WGDRsgiiKmTZuGUaNGtULTkkve+z3E8Lu0pC2edeoqlQr5aZ1RaamUD3Np3Uydw+9EREoUNqifOXMGbrcbI0eOxMiRI7Fp0yb07du3tdqWVNI6dVeEJW3xZOqAt1jux8rNsDgtABK/mxwA6AVp73dm6kREShRyTn3Hjh0YM2YMysvL5dfWr1+Pm266Cbt27WqVxiWTLsIpbVKhXDx7vwPeYjm36EaF6RiAxJ/QBvgvaWNQJyJSopBB/YUXXsDcuXMxcuRI+bWHHnoIpaWlmD17dqs0LpnkHeVCVb83Y0kbAOSne4vlDtcfAtA6w++CWoCgEhjUiYgUKmRQr6+vx0UXXdTk9ZKSEtTW1ia0UWcD33nqoQ90SdWkxr1fu7SsTQ7qrZCpA95snUvaiIiUKWRQd7lc8Hg8TV73eDzt4kAXafjdFWJJm8lpjKtITtI0qCd+Th3wFstxTp2ISJlCBvULL7wQ8+bNa/L6G2+8gQEDBiS0UWcDafg9dPW7Ma593yW+4ffDAFpn+B3wFstx+J2ISJlCjh0//PDDuO+++7B06VL069cPer0eO3bsQIcOHfDmm2+2ZhuTQj7QJdQ6dYcJnRsOZolHJzlTPwig9TJ1Dr8TESlXyKBuMBjw/vvvY+PGjdi5cyfUajVuv/12DB06tDXblzQaeZ1606kGt8cNi8scd+U7APkLwTHjUQBAipD4U9oA7/B7nb2uVT6LiIhaV9gqL5VKheHDh2P48OGt1Z6zhk6Q5tSbBnVzM9eoA0CHlA7QqDXynH3rzakzUyciUqqIe7+3V9KBLsECYHP2fZeoVWp0Ss2Xf2696ncWyhERKRWDegjhzlP37SaX2azPkIrlACC1FTN1u9sOURRb5fOIiKj1MKiHIKgFqFXqEJm6tO97/Jk64FvWBrRu9bsIMeRSPSIiarsY1MPQqXVB59R9u8k1L6h38gvq+lYcfge4VSwRkRIxqIehFXRBq99Njobh92ZsPgMEDr+nCK03/A4ANpetVT6PiIhaD4N6GFq1Jug6dZOz4TCXZlS/AwhY557a2pm6i5k6EZHSMKiHoVXrgu4oJ53Q1pwlbQCQn+Zf/d66mTqH34mIlCdhQd3j8eCZZ57BhAkTMHnyZBw+fLjJe6xWKyZOnIj9+/dHfU9r0gm6oAVlLbGkDQgslNO31uYz6oagzkydiEhxEhbUV61aBYfDgYULF+KRRx5pclzr9u3bcfvtt+Po0aNR39PatGpt0Op3UwtsPgMA+el+1e+tlKmnaJipExEpVcKC+pYtW1BSUgIAGDhwIMrLywOuOxwOvP766+jdu3fU97Q2rVob9Dx135K25gX13NQ8qKAC0Hpz6iyUIyJSrvgOA4+CyWSCweAbnhYEAS6XCxqN9yOHDBkS8z3B5OSkQaMRWrDlQF6eN1in6lPgNDvlnyUuwZvl9uzSBXlZzQvsndI7odJcie6dO0HbsDVtIuVkeNtrd9mR1615bW8rGv/5KRX7qSzsp7K0Vj8TFtQNBgPMZrP8s8fjCRuc472nttbSvIY2kpeXgaoq75y52iPA6XbKP0tO19cAAOz1QJXD2OQZschP64Jaay3qamwAEp89Swm63W1v0i8l8v/zVDL2U1nYT2Vp6X6G+4KQsOH3wYMHY+3atQCAsrIyFBcXJ+SeRPKuUw89p97cJW0A8Ozw5/F/l7/W7OdES65+Z6EcEZHiJCxTHz16NNatW4eJEydCFEWUlpZi2bJlsFgsmDBhQtT3JJNWrYVH9MDtcUNQ+4b4jQ4jUjWp0Kib/5+vpNtlzX5GLLijHBGRciUsqKvVasyYMSPgtaKioibvW7BgQdh7kkkrn6nuQKraV51uchqbXSSXLNJ2tMzUiYiUh5vPhKFryGob7/9udBibve97sugajpRl9TsRkfIwqIfhO1M9MKibHKZmH7uaLHruKEdEpFgM6mFoG+bM/dequz1uWFzmZh+7miwslCMiUi4G9TC0DcPvTr/hd3ML7SaXLCyUIyJSLgb1MKT5Z/+T2lpq3/dkYaEcEZFyMaiHIe3w5n+mum/f97Y5py4f6MJMnYhIcRjUw5CWtPln6vJZ6m01UxdY/U5EpFQM6mFI1e/+c+rGFjpLPVlYKEdEpFwM6mHogg2/Oxq2iG3jmTqH34mIlIdBPQxtkEI5afi9zc6pc506EZFiMaiHIc+p+61TNzUMv7fVHeX0AqvfiYiUikE9DN86dZf8mjSn3mb3fufwOxGRYjGoh6GTDnQJGH6Xjl1tm5m6NPzO6nciIuVhUA/Dl6n7bz5TD6Atz6k3ZOocficiUhwG9TC04TL1Nlr9rlapoVVrOfxORKRADOphSEHd5Tenbmrj69QB7xA8M3UiIuVhUA9DGqp2eJpm6m1173cASBH0zNSJiBSIQT0MTZBtYo0OI1I1qdA0HMvaFukEPQvliIgUiEE9DLn6PeBAF2ObXc4m0Qk6Dr8TESkQg3oYUvW7yx2493tbXc4m0XP4nYhIkRjUw5DOUw+YU3eY2uxyNgkL5YiIlIlBPQzfnLo3U3d73LC4zG12OZukY0pHGB1GnLHXJbspRETUghjUw/Cd0ubN1M0Nle9teTkbAAzqNAQAUHZqa5JbQkRELYlBPQzplDZXQ6GcEpazAcCgfG9Q33pqS5JbQkRELYlBPQx5nXrD8LvR0baPXZUMbsjUf2RQJyJSFAb1MKS16NLe79JZ6m19Tj0/vTO6ZnTFj5WbIYpisptDREQthEE9DKn63elpnKm37Tl1ABjWdRhOWSpxwlyR7KYQEVELYVAPQz6lrWH43eRo24e5+BvWdRgA4MdKDsETESkFg3oYvh3lAoff2/qcOuAL6mWnfkxyS4iIqKUwqIehEQL3fpdOaGvrO8oBwJAuQ6CCihXwREQKwqAeRqg59ba+9zsAZKVkoU92X2w99SM8oifZzSEiohbAoB6GdJ66s9E6dSVk6oB3vbrJacS+2r3JbgoREbUABvUwVCoVtGotHA3D70ZHPQBlzKkDvp3lfjy1OcktISKilsCgHoFWrW2aqSug+h3wbULDYjkiImVgUI9AK+h8QV1B69QB4NzcAdCqtSyWIyJSCAb1CLRqra/6XSF7v0v0gh4Dcs9D+entPF+diEgBGNQj0Kl18jp1o8OINE2avH2sEgzqNAROjxM/n96e7KYQEVEzJSw6eTweTJ8+Hbt374ZOp8PMmTNRWFgoX1+9ejVef/11aDQajB07FrfeeisA4Oabb0ZGhnd4u1u3bpg1a1aimhgVjaD17SjnNComS5d4i+X+ia2ntmBw/tBkN4eIiJohYUF91apVcDgcWLhwIcrKyjB79my8+eabAACn04lZs2Zh0aJFSE1NxaRJk3D55ZcjM9NbVb5gwYJENStmOrUWVqcFgDdTV8p8ukQK5D9WbsHU85LcGCIiapaEDb9v2bIFJSUlAICBAweivLxcvrZ//3706NEDWVlZ0Ol0GDJkCDZv3oxdu3bBarXi7rvvxpQpU1BWVpao5kVNq9b5TmlzmGBQWFAvyu6DDF0mK+CJiBQgYZm6yWSCweAbqhYEAS6XCxqNBiaTSR5iB4D09HSYTCakpKRg6tSpGD9+PA4dOoR7770XK1asgEYTupk5OWnQaIQWbXtenq9tafoUOI1OdOiYBovLjA7p2QHX2zKpHxd2HYrVB1dDl+FBVkpWklvV8pTy5xUJ+6ks7KeytFY/ExbUDQYDzGaz/LPH45GDc+NrZrMZGRkZ6NWrFwoLC6FSqdCrVy9kZ2ejqqoKXbp0Cfk5tbWWFm13Xl4GqqqMvhc8ajjcDhyqOAEA0CM18Hob5d/PAdkDsRqrsWrnWozsNiq5DWthTf48FYr9VBb2U1laup/hviAkbPh98ODBWLt2LQCgrKwMxcXF8rWioiIcPnwYdXV1cDgc2Lx5MwYNGoRFixZh9uzZAIDKykqYTCbk5eUlqolR0TWsU1fSvu+NDcr3bkKzlcewEhG1aQnL1EePHo1169Zh4sSJEEURpaWlWLZsGSwWCyZMmIDHH38cU6dOhSiKGDt2LPLz8zFu3Dg88cQTmDRpElQqFUpLS8MOvbcGbcOhLrX2WgBQ3Jw6AAzqNBgA8CM3oSEiatMSFjHVajVmzJgR8FpRUZH8+yuuuAJXXHFFwHWdToe5c+cmqklx0TasSa+11QBQzm5y/rqkFyA/rTN3liMiauO4+UwEWqEhU28I6krZ992fSqXCoPwhOGk+gROmimQ3h4iI4sSgHoF0pnqNgoM64DvcZSuXthERtVkM6hFIZ6rX2KoBKOfY1cYGyUGdQ/BERG0Vg3oEWsEb1OXhd50yM/WBnQYBYLEcEVFbxqAegbbJ8LvyCuUAIEufjZ6ZvfDz6W0QRTHZzSEiojgwqEegayeZOgAMyD0fNbYaVJiOJ7spREQUBwb1CHyZurLn1AFgQK73RJfyah7DSkTUFjGoR+ArlFN29TvgF9RPb0tyS4iIKB4M6hHIhXINO8opcfMZyXm5FwAAyk8zUyciaosY1COQ1qmfsdcBANIVnKl3Tu+CjikdmakTEbVRDOoRSDvKAUCaJg0adXL3ok8klUqFX+Sej8P1h1BvP5O0dlRbq3HlRyX4/MCypLWBiKgtYlCPQOsXxJWcpUukefWfq8uT1oaVh5Zj++mf8N6O/yStDUREbRGDegRS9Tug7Pl0yXm55wOIXCzncDuwdO9iuDyuFm/Dd8e9R/ZuOrExIc8nIlIqBvUIdH7D70o8drWxAXJQD18s93b5P3Hf//s1/r39Hy36+aIoYt3xbwEAJqeR8/tERDFgUI9AWtIGABkK3U3OX1F2H6QIKRGD+ldH/h8A4KM9H7bo5x88sx8V5uPyzn3rK9a16POJiJSMQT0C/6Cu5N3kJBq1Bv07notdNTvgcDuCvsfqsmJjxXoAwLaqMuyq2dlin/9dQ5Z+z3n3AwDWN/xMRESRMahHEDD83g4ydQAYkHsBnB4n9tTuDnp9Y8V62Nw29M4qAgB8vLvlsvV1DfPpt54zCT0ze2HjiQ1we9wt9nwiIiVjUI/Av1CuPcypA5F3lvv66GoAwPRL/ooMXSYW7/kIHtHT7M8VRRHfHf8W+WmdUZTdB5cUXIp6xxnsSGIlPhFRW8KgHoF0oAvQPqrfAb9lbSHm1b8++hVShBRc1v1y3Fh0MyrMx+XitubYU7sbVdZTGNG1BCqVCsMLRgAA1ld81+xnExG1BwzqEWj859TbwTp1AOjf8RdQQRW0WO6k+QR21uzA8IIRSNWkYnzxRADAxy1QMCctZSvpehkA4JKulwIA1jGoExFFhUE9Al07W6cOeL+89M4uQnn19iZnq0tD76O6XwkAuLjgEnTP6IFl+z+FxWlp1udKQX1E1xIAQPeMHuiRUYhNFetbZHifiEjpGNQj0Ar+mXr7COqAdxOaM/Y6HDUeCXj966NfAQAu7+EN6mqVGuOKb4XZacLyg5/F/Xke0YP1x79FN0N3FGb2lF8fXjACtfZa7KzeEfeziYjaCwb1CHTtsFAOCL4JjUf04Juja9AlvQDn5PSTXx9fPAlA84bgd1T/jFp7LS7tNhIqlUp+/ZIC7xD8Bg7BExFFxKAeQXucUweCV8Bvr/oJ1bZqjOp+RUDg7ZPTF4M7DcHXR1ej0nwy6PM+2v0B/l3+z5Cf993xbwAAIwpKAl6X5tW5CQ0RUWQM6hG0x+p3APhFkD3g1zQMvY/qfkWT948/ZyI8ogdL9i4KeN0jevD0uifw+6/ux+NrH8HqI6uCfp5UPX9p15EBr/fIKERXQzdsqPiuyfw+EREFYlCPoD2uUweA/LR8dErLDxh+//roaqigwshulzd5/819xkGj1gQMwdtcNty/8m7M/+l19MrqDbVKjcfXPgKryxpwr8vjwvqKdeiZ2QtdM7oFXJOWtlXbqrG7dlcL95KISFkY1CNob3u/+xuQex6OmY6i1lYDk8OIH05uwgV5A9ExtWOT93ZM7Ygre4xG+elt2FH9M+pstZjw2S34dP8SDC8YgS/HrsG95/8Gh+oP4tUf/y/g3u1VP8HoqEdJt8uCtkOaV+d6dSKi8BjUI9AGnNLWfubUAWBAR+8Q/M/V5VhX8R2cHqdc9R7Mred4C+Ze+/ElXP/J1dhQsQ43Ff0KC6//BNkpOXjswifRJb0Ar/34EvbX7ZXv+67CO/QuLWVrTJpX33Cc8+pEROEwqEegU7fPJW1AYLHc1/J8euigPrrwl8jUZWHx3o+wp3Y3pl3we8y/+t9I0aQA8E5fzLx0NhweB/689hF5jlza771xkZykV2ZvdE7vgnUV33JenYgoDAb1CKRMPU2TBkEtJLk1reu8PG+mvr1qG9Yc+QrpWgOG5g8L+f4UTQom9b8DKqjw/IhZmDGiFGpV4F+x63vfhCt7jMa3x77GJ/sWwel2YmPFBhTnnIP89M5Bn6tSqXBJwQictlZhn1+GT0REgRjUI5Dm1NPb0XI2Sc/M3kjTpGPN0VU4cGY/SrqODNiMJ5jpw2di2117cP8Fvwt6XaVSYVbJi0gRUvD0d09g7bE1sLjMIYfeJcPbwLy6xWnBnctvw9j/3djs3fWIiOLBoB6BoBKggqpdLWeTCGoB53b8BU5bTwMARoWZT/e/Jz8tP+x7emb1wkNDHkWV9RR+/5X33PTGS9kak4bmz9ZNaKwuK6Ysn4TlBz/Dt8e+xmNrH+ZUARG1Ogb1CFQqFVI1qcjUZSa7KUkhzasDwdenx+u3gx5An+y+qLZVAwAuCTGfLinK7oO81E5Yd/zsW69ud9vx6xW3Y+2xNbim57UY1GkwFu7+Lxbs+E+ym0ZE7QyDehRmlbyIxy96OtnNSAppu9jCzJ7onVXUYs/VC3r87bKXAAC/6Hhe0GVy/rzz6pei0nIS++v2tVg7msvhduCeL6dg9ZFVuKrH1Xjrmnfx1jXvIkefgye/fRRbK7cku4lE1I4wqEdhUv87cEWPq5LdjKQYkn8hAGB04TUt/uxLu47EW1e/g5cvnxfV+6XldJM+H3tWBEun24n7Vv4aXx5ajsu6XY5///I96AU9umf0wN9H/xtOjxNTv5yCamt1sz9LFEV8vPtD/Pmbh1Brq2mB1hOREjGoU1i/yB2A/93yJZ686JmEPP/GPrfggk6DonrvhHNuw4ODH8GR+sMY88lovFH2WtKOZHV5XPjdV/fii4PLcGnXkXjn2g/kpXuA9wvIn4c9iWOmo/jNqqlwe9xxf5bJYcRvV92L3311H/7z879wzaLLsatmZ0t0g1pQna0W/9o+Hx/t/gA/ndoKs9Oc8M+sNJ/Eqz/+Hw6eOZDwzwpHFEXsqtmJ5Qc/b7JjZEs4Un8YL/4wG18c+Aw2l63Fn68kKvFsm6CMUVWVsUWfl5eX0eLPPBu15X5+c3QNfrvqXlRZT+HKHqPx2pXzkZuaG/S9/v10eVzYdGID1hz5Cv07notf9R0fcDBNMHW2WnyybzGOGo+gwnQcJ8wVOGGqwAlzBexuOy7ucgk+uH4x0rXpTe71iB7c8fmtWHVkJR4Z+hgeG/YX+ZrT7cQx01GcNJ9Avw79kZPSIejnb6sqw70r78LBMwcwJH8oLux8Mf7+0zwYtBl4c/RbuKbntU36KdlTsxs6QYfCzJ4R+5ksZ+x1sLlscHgccLodcHiccHqcUEON3NRcdEjpGLDiojl/b0VRxPqK7/BTVRluKLoJ3TN6tEgfRFHEp/uW4C/fPYYq66mAa90zeqBvTjGKc/qhf4dz0a9DfxR36Bf074u/SP10eVz41/b5eOH7UpicRmToMvHSqNdwY59bWqRPAHDMeBTfn9wIvZCCXlm90TOzF9K0afJ1t8eNH05uwvKDn2P5wc9wqP4gACA/rTMeGPwQ7jj3LqRqUsN+RqR+Ot1OzN/2Bub8UCp/WcjQZeK6Xtfjlr5jUdJ1VMQVOZEcOLMfyw98jtzUXIzsNgpdDAXNel4wLf3vbV5e6MLthAV1j8eD6dOnY/fu3dDpdJg5cyYKCwvl66tXr8brr78OjUaDsWPH4tZbb414TzAM6vFp6/08ZTmFP3x1P9Yc/Qr5aZ3xt8tewrkdf4FsfTYydJny+vj0bAEf//gplh/8DCsPLUetvVZ+xujCazB31KvonN6lyfNFUcSy/UvxxLePBvxDrYIKndLyUWAowIDc8/HcJX8NeyZAra0Goz++DEeMh/GrvuNRZTmFw/WHcMx0VB5lEFQCLim4FNf1vh7X9roeBYauEEUR/9z2Jp7b8DScHif+MOghPD7sKWgFLZbuXYwH1/wWNpcNT170DB4Y/DA6dcpEVZURJocRi/d+jAU7/oNtVWUAgLzUTriw80UY1uViXNh5GM7PGwi9oG+JP4a4mJ1mfLpvCRbs+A+2VP4Q8f05+hzkpuYhNy0PF3Q5D6O7jsHwghHQqDVRfV6lpRILd72P93e+K2e0GrUGE8+5HQ8OeQSFmT3j7stR4xE89s3DWHVkJVKEFDww+GF0SO2IvbW7sadmN/bU7kalJfDkQhVUKMzsiX4dz8U5Of3QM7MXemb1Qs/MXuhiKIBapQ77/8+NJzbgsW8exs6an5Gtz8bEfnfg3Z/fhsVlxl2/mIoZI2YFjBpJRFHET1VbsatmJ3JSOqCD/KsjsvTZOGGqwLqKb7Hu+LdYV/EdjtQfavKMzuld0CurNzqm5GJDxXdyoWu61oAre4xG5/TOeG/Hu7C4zFEF93D93FL5Ax75+kHsqC5Hbmou/nThEzhafwRL9y3GcdMxAECHlA4Y0/sm3FB0E0YUlEQd4OvtZ/Dp/k+wcNd/8f3JjQHX+mYXo6TbZSjpNgrDOl8Mo+MMDtUfwuH6QzhSfxiH6w/hjL0ORdl9cG7HAejf8Rc4t+O5yAhTTK2IoL5y5UqsXr0as2fPRllZGebPn48333wTAOB0OnHddddh0aJFSE1NxaRJk/D3v/8dW7duDXlPKAzq8VFCPz2iB2+UvYbSTc/B5XHJr6tVamTqMpGlz8YpS6X8Db9zehdc22sMLut2Bf5d/k+sPbYG2fpslJbMwdi+t8rZbIXpOB5b+zC+PLQcKUIKHhzyCEq6jkKBoQD5aZ1jzgy2V/2EMUtGw+a2ye0ozOyJwsyeyEvthA0V3+HHU74agUGdBsOgzcC3x79Bbmoe5l05v0lNx7aqMkz5YhIqzMdxS5+xeKTkIfx94z+xZO8iWFxmCCoBowuvgV5IwfcnN+KEuUK+Vy/o0b/DuTgvbyDOyz0f5+Wdj3M7DoiYVTXm9rhx0nwCR41HcMR4GEeNR3DceAzpOgP6ZPdFn+y+KMrug/y0zlCpVNh+ehsW/Pw2Fu/9GEZHPVTwFj92SusEraCDTq2DRq2BTtDBLbpRbT2N09bTOG2twmlrFaqt1RDh/ecqNzUP1/e+ETf1+RUu7nKJvDGU2+NGta0ap61VOHjmAD7a/QFWHloOt+hGqiYVNxTdjEGdBuNf2/+BfXV7IagE3HrOJDw45JGAQlCXx4UqyymcslTCLbqRnZKDHH0OMnVZENQCXB4X3tr+d8zeNBMWlwUl3UZhzmUvBS0mPWOvw66aXdhdsxO7anZgZ/UO7Kz5GTVBaiN0ah16ZBaib24f5OsL0D2zEIUZheie0QNZ+iy8tOVFLNz9XwDA7f2n4KmLn0PH1I7YW7sH93x5J3bW/IwBuefjrav/g97ZfQB4v9Qs2r0QC3e/H3LaRq1SB0xlZemzMbxgBIZ3GQGP6MHBMwdwsP4ADp05gGPGoxAholNaPq7peR2u6zUGI7qOlL9IVFur8WbZa3hr+3w5uE8+9y70yuqNLoYCFKQXoHN6AdK0aUH/Haq3n0Hpphl4u/wtiBBxR/878fTw5+TRLI/oweaTP2DpvkX43/6lOGWpBOD98vfLXmNwfe8bMbL75fIXV7fHjVOWShwzHcWR+sNYeWg5lh/8HDa3DSqocGm3yzCu762os9fh22NfY33FOlhcsU+d9MgoRJ+cvihI74ouhgJ0SS9Al/Qu6GLoihHFQ1Fb3XLTEkkJ6rNmzcL555+PMWPGAABKSkrw7bfePb537dqFOXPm4F//+hcAoLS0FIMGDUJZWVnIe0JhUI+PkvpZdupHLN7zEWrttThjr0OdvU7+345pHTC6+7W4ttcYXNBpkJzBi6KI//z8Lzy3/mlYXGZc1+sGvHDZ/+Gz/Z/irxufg8lpxIiCEswd9Yr8j2NzHDUegdVpRffMHkGDZ4XpOJYf/BxfHPwM649/C7foxshul+P1q/4Rct3/Kcsp/HrF7fjh5Cb5te4ZPXBH/zsxqf8d8giEKIo4bjqGH05uavj1PXZW/wyHxyHfJ6gE9MgshFqlhtPtDBgOd3tcAFRQq9RQqbz/q4YKRqcx4MtUKAZtBjqmdsThhsyvS3oBbus/Gbf1nxzTELjL48IuSxne2fI+Pj/wqbx/Ql5qJ3RM7dgk8EvOy70Ad5x7J37Vdxyy9NkAvP/Qf7p/Cf5v89+wp3Y3BJWAi7oMR72jHifNJ1BtPd3kOYA3y87SZ0Gj1uC09TQ6pHTAc5eU4tZzJsU0xSGKIqqsVdhbuxuH6w/h0JmDOFR/QP69/4hSY+flXoAXRs7F0M6BuztaXVY89d1jWLDjP0jXGvDAoIewufJ7rD6yCm7RDZ1ah1/2GoOR3UbB6DCixlaNGls1qm3VqLFWo0NKB1zS9VKMKCjBuR0HhNxB0+6245SlEl0N3ZrsGOmvcXBvLFufjQx9BhwuJ9yiG26PC27RA5vLCofHgb7ZxXhx1CsYXjAi5Ge4PW5sPLEenx34FJ8fWIaT5hMAvEP0/Tr0R6X5JCrMx5v8PS3K7oMJ59yGccUT0C2je8A1h9uBH09twdqja/BT1VZ0TM1FYWZP9MgoRGFmLxRm9USGNgN7a3djZ80O/Fxd7v2iVv1zk+kXya2/uBXzLnsrZD9ilZSg/pe//AVXX301LrvMe/LWqFGjsGrVKmg0GmzevBnvvfceXn75ZQDAK6+8goKCApSVlYW8JxSXyw2Npn1t30ot50DtAfz6019j7eG10Kg1cHlcyE7JxoujX8Tdg+5Oylx0jbUG+2v2Y0jBkLD/aAKA3WXH46seR4WpAncPvBtX9b4qqu2MHW4HdlbtxNaTW/HjiR+x9eRW7KneA0ElQCtooRN08i9BJUCECFEU4RE9EOH9X4POgF7ZvdAzu6f8q0dWD9TZ6rCneg92n96N3dXeX0fOHMFlhZfhviH34Zd9fhn10HkoLo8L3xz6Bh/v+BhLdy2Fw+1Ap/ROyDfko1N6J3RK64TOhs4YUzwGg7sMDvkct8eNxTsX4/m1z6P8VDnStenoktEFXQxd5P/VqDWotdaixlaDGqv3V52tDqN7j8bfRv8NuWnB6zma44ztDA6fOYxDdYfkX8fqj+HynpfjviH3hf0z/mD7B7jvs/tgcpgAABcWXIi7Bt6FiQMmokNq8NqNRKq2VGP90fU4bjyO4/XHcdx4HMfqj+G48TisTisEtQCNWgNB5f1fraDFTefchEcveRR6TfTTRB7Rg03HNmHRjkVYvHMxjpw5goKMAvTI6hHw68KCCzGs67CE/H/baDf6+tfQ1+P1x3Ft32txffH1Lf55wSQ0U7/gggtw3XXXAQBGjhyJtWu9B3fs2rULc+fOxT//+U8A3kx98ODB2Lp1a8h7QmGmHh/208cjevDWtr9j9vd/xRU9rsJfS/4WcVe8sw3/PJtHFEVYXJaIBWytpbn9PHBmP9YcWYURXUeiX4f+LdiylpXIP0+36G72l8eW0ppz6glb0jZ48GA5IJeVlaG4uFi+VlRUhMOHD6Ourg4OhwObN2/GoEGDwt5DlChqlRr3XfBb7LvnKN665p02F9Cp+VQq1VkT0FtC76wiTD3v/rM6oCeSSqU6awJ6a0tYr0ePHo1169Zh4sSJEEURpaWlWLZsGSwWCyZMmIDHH38cU6dOhSiKGDt2LPLz84PeQ9RaIg11ExGd7bhOvREOYyoL+6ks7KeysJ/xPy8UpiZEREQKwaBORESkEAzqRERECsGgTkREpBAM6kRERArBoE5ERKQQDOpEREQKwaBORESkEAzqRERECsGgTkREpBBtfptYIiIi8mKmTkREpBAM6kRERArBoE5ERKQQDOpEREQKwaBORESkEAzqRERECqFJdgPOFh6PB9OnT8fu3buh0+kwc+ZMFBYWJrtZLeqnn37Ciy++iAULFuDw4cN4/PHHoVKp0LdvXzz77LNQq9v2dzyn04knn3wSx48fh8PhwG9+8xv06dNHcf10u9146qmncPDgQQiCgFmzZkEURcX1U1JdXY1f/epX+Pe//w2NRqPIft58883IyMgAAHTr1g3Tpk1TZD/nz5+P1atXw+l0YtKkSRg2bJji+rlkyRJ88sknAAC73Y6dO3fiv//9L0pLS1unnyKJoiiKX375pfjYY4+JoiiKW7duFadNm5bkFrWsf/zjH+L1118vjh8/XhRFUbz//vvFjRs3iqIoik8//bS4cuXKZDavRSxatEicOXOmKIqiWFNTI1522WWK7Of/+3//T3z88cdFURTFjRs3itOmTVNkP0VRFB0Oh/jb3/5WvPrqq8V9+/Ypsp82m0286aabAl5TYj83btwo3n///aLb7RZNJpP46quvKrKf/qZPny5++OGHrdrPtv2VqAVt2bIFJSUlAICBAweivLw8yS1qWT169MBrr70m//zzzz9j2LBhAICRI0di/fr1yWpai/nlL3+JBx98UP5ZEARF9vOqq67C888/DwCoqKhAbm6uIvsJAC+88AImTpyITp06AVDm39tdu3bBarXi7rvvxpQpU1BWVqbIfn733XcoLi7G7373O0ybNg2jRo1SZD8l27dvx759+zBhwoRW7SeDegOTyQSDwSD/LAgCXC5XElvUsq655hpoNL7ZFlEUoVKpAADp6ekwGo3JalqLSU9Ph8FggMlkwgMPPIA//vGPiuwnAGg0Gjz22GN4/vnncc011yiyn0uWLEGHDh3kL9uAMv/epqSkYOrUqfjXv/6F5557Dn/6058U2c/a2lqUl5fjlVdeUXQ/JfPnz8fvfvc7AK3795ZBvYHBYIDZbJZ/9ng8AUFQafznc8xmMzIzM5PYmpZz4sQJTJkyBTfddBNuuOEGxfYT8GaxX375JZ5++mnY7Xb5daX0c/HixVi/fj0mT56MnTt34rHHHkNNTY18XSn97NWrF2688UaoVCr06tUL2dnZqK6ulq8rpZ/Z2dm49NJLodPp0Lt3b+j1+oDgppR+AkB9fT0OHDiAiy++GEDr/nvLoN5g8ODBWLt2LQCgrKwMxcXFSW5RYp177rnYtGkTAGDt2rUYOnRoklvUfKdPn8bdd9+NRx99FOPGjQOgzH4uXboU8+fPBwCkpqZCpVJhwIABiuvn+++/j/feew8LFixA//798cILL2DkyJGK6+eiRYswe/ZsAEBlZSVMJhNGjBihuH4OGTIE3377LURRRGVlJaxWK4YPH664fgLADz/8gEsuuUT+uTX/HeKBLg2k6vc9e/ZAFEWUlpaiqKgo2c1qUceOHcPDDz+Mjz76CAcPHsTTTz8Np9OJ3r17Y+bMmRAEIdlNbJaZM2di+fLl6N27t/zaX/7yF8ycOVNR/bRYLHjiiSdw+vRpuFwu3HvvvSgqKlLcn6e/yZMnY/r06VCr1Yrrp8PhwBNPPIGKigqoVCr86U9/Qk5OjuL6CQB/+9vfsGnTJoiiiIceegjdunVTZD/feustaDQa3HXXXQDQqv/eMqgTEREpBIffiYiIFIJBnYiISCEY1ImIiBSCQZ2IiEghGNSJiIgUgkGd6CxhNpvx3HPPYfTo0bjxxhtx2223YcOGDRHvW7NmDd5++20AwAcffIAPPvgg6s/8y1/+gu3bt8fc1uZ8ZiyOHj2KJ598MiHPJlIi5W6ZRtSGiKKIadOmoX///vj888+h0+mwY8cO3HfffZg7dy4uuuiikPf6n1MwadKkmD73r3/9a1ztbc5nxqKiogJHjx5N2POJlIZBnegs8P3336OiogLvvvuuvEf0ueeei9/85jd44403cNFFF2Hy5Mno168fNm/eDLvdjieffBKdO3fGhx9+CAAoKChARUUFAOAPf/gDRowYgSuvvBLbtm1Dbm4uxo4diwULFuDkyZOYPXs2hg0bhsmTJ+P3v/89du/ejcWLFwMAbDYbjh49im+++QYHDx7ESy+9BJvNhvr6ejzxxBPo2bNnyM9cs2YNXn75ZXg8HnTv3h0zZsxAbm4urrjiCtx444347rvvYLVa8cILL2DAgAEB/w3efvttfPLJJ1Cr1Tj//PMxY8YMzJw5E8eOHcNzzz2HZ599Fv/4xz+wfPlyuN1uXHrppXj00Udx/Phx/OY3v0Hv3r2xb98+FBQUYM6cOUhPT8eTTz6JvXv3AgBuu+023HrrrYn/wyRKIg6/E50Ftm/fjgEDBsgBXXLhhRcGDI+bTCZ88sknmDt3Lh5//HH06NEDEydOxMSJEzF27NiAe0+fPo2RI0di6dKlsNvtWLVqFf773//iD3/4A955552A906ZMgWffvopli5dir59++Lhhx9GXl4e3nvvPcycOROffPIJZs6ciVdeeQV9+vQJ+pnV1dV45pln8Prrr2PZsmUYPHgwZsyYIV/Pzs7GokWLMHHiRHmbW4nb7cb8+fOxePFiLFmyBE6nE5WVlXjqqacwYMAAPPvss1i7di3Ky8uxaNEiLF26FJWVlfjf//4HANizZw9uu+02fP755ygqKsK8efOwdetWnDlzRt5Wd/Pmzc37QyJqAxjUic4CKpUKbre7yetOpzMg0EuZZv/+/ZGXl4fdu3eHfe7IkSMBAF27dpUPlygoKEB9fX3Q97/yyivQarW45557AABz5szB3r178frrr+Ptt98OOPSosW3btuH8889Ht27dAAATJkzAxo0b5evSaWt9+/ZFXV1dwL2CIGDQoEEYN24c5s2bh1//+tfIz88PeM+GDRuwbds2/OpXv8Itt9yC8vJy7Nu3DwDQs2dPeYri5ptvxsaNG9G3b18cPHgQU6dOxYoVK/DnP/857H8rIiVgUCc6C1xwwQUoLy+H0+kMeL2srCxgmNp/v+hoThLU6XRB7w1mxYoVWLNmDUpLS+XXbrvtNmzbtg0DBgzAtGnTwt7v8XgCfhZFMeD4Yr1eDwBNRiMkb7zxBqZPnw5RFHHPPffg+++/D7judrtx55134tNPP8Wnn36Kjz/+WG5T42OFBUFATk4OPv/8c9xxxx04ePAgbrnllpBfZoiUgkGd6CwwdOhQ9OnTB6WlpXJgLy8vx5tvvonf/va38vu++OILAN7h+vr6ehQXF0MQhIDgGY+dO3fihRdewLx585CamgoAqKurw6FDh/Dggw9i5MiR+Oqrr+TRhGCfecEFF+Cnn37CsWPHAAALFy4MW+Dnr6amBtdddx2Ki4vx4IMPYsSIEdi9e3fA51x88cX49NNPYTab4XK58Lvf/Q5ffvklAO+BGTt37gTgPbJVau+jjz6KUaNG4amnnkJaWhpOnDjRrP9ORGc7FsoRnSXmzZuHl156Cddffz0EQUBWVhbmzJkTEBiPHj2KW265BQDw0ksvQRAEXHjhhXjssceQm5sb92fPmTMHLpcLDz74oBy4n376aYwbNw5jxoyBRqPBxRdfDJvNBovFEvQzc3NzMWPGDPz+97+H0+lEQUFB1NX1HTp0wIQJEzBu3DikpqaiV69eGDt2LOx2O4xGIx599FHMmTMHu3btwq233gq3242SkhLccsstOH78OLKysvDqq6/iyJEjOOecczBz5kxotVqsXLkSY8aMgV6vx4033ohzzjkn7v9GRG0BT2kjaiOkSvVos9/24tixY5gyZQpWr16d7KYQJR2H34mIiBSCmToREZFCMFMnIiJSCAZ1IiIihWBQJyIiUggGdSIiIoVgUCciIlIIBnUiIiKF+P//a/1ndWpDvQAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 576x396 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.style.use(\"seaborn\")\n",
    "plt.plot(cost_fn, \"g\")\n",
    "plt.ylabel(\"Cost function\")\n",
    "plt.xlabel(\"Optimization steps\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAbAAAAETCAYAAABXxlNqAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAo70lEQVR4nO3df1xU9Z4/8BfMAIKjKC2yd69JNxQrzYt4TUlJ+ZWpoAmrAyhgW1hbahayUioOaEri3Vvir7TVTFshla6iloaiGOpVVEz0+mO1UFgV/IEyQ/ycz/cPvs4yF4chmhE++no+Hj0enXNmznnPBz/zOudzzpxjI4QQICIikoxtWxdARETUGgwwIiKSEgOMiIikxAAjIiIpMcCIiEhKDDAiIpISA8yM+vp6rFu3DqGhoRg7dixGjRqF1NRU1NTUAAASEhLwX//1Xxbd5t69e7FgwYJWvbe4uBj9+/e3SB1paWlITk42+7o5c+agsLAQADB79mwcOnTIItsnImoOA8wMjUaDkydPYv369di2bRu2bNmCn376CbNnz7baNgMCAjBnzhyrrd/SDh06hPs/J/zoo4/w4osvtnFF1F7V1dVh9erVCAkJQUhICEaPHo2kpCTcuXPHatvcvHkzvvrqK6ut/9eKiorCd9991+xrKioqEB0dbZgeO3Ys7t27Z+3SpKNs6wLas+LiYmRlZeGHH36ASqUCADg5OSEpKQknTpxo8votW7YgIyMDtbW1uHv3LmJjYxEZGYmysjLMmjXL0EmHDRuGGTNmmJyfmZmJ3bt347PPPkNZWRnmzZuHy5cvw9bWFuHh4YiOjkZBQYHhSLCsrAwvvvgiFi5caPKz1NXVYf78+Thx4gTs7OzQvXt3LFq0CB07dkR2djaWLVsGvV6Pjh074oMPPkC/fv2M3u/v749PP/0Uzz//vNF0dnY2SktLMXPmTCxevBhLlizBxIkT8corr5hcb1paGkpKSlBWVoaSkhK4ubkhNTUV3bp1s8jfjdqv+Ph46PV6bNy4Ec7OzqitrcUXX3yB8PBwbN261dDPLOn48ePo1auXxddrTXfv3sXp06cN09u2bWvDatovBlgzzpw5g549ezbpVK6urhgxYoTRPJ1Oh82bN2P16tXo2rUrCgoK8NprryEyMhJff/01unfvjrVr16KyshKzZ89GRUWFyfmNJSUl4amnnsKKFStQUVGBiIgIDBs2DF9++SWmT5+OQYMGQafTISAgAIWFhejSpcsDP0tBQQGOHj2KXbt2wcbGBqmpqTh//jycnZ0xb948pKen48knn8Thw4fx9ttvm91DvO+9995DVlYWlixZYgg3ALh06VKz683Pz8df//pXqFQqvPXWW0hPT8f06dNbtE2S048//ohjx47h+++/h6OjIwDAzs4OsbGxOHHiBDZt2oTY2FiTO0vPP/88Vq1ahb1796Kqqgq//PILZs2ahaCgIJM7RadOncK+ffuQl5eHDh064Pbt27hz5w4SExMBNAyT35+OiopCnz59UFBQgNu3b2PChAm4efMmjh49il9++QWffPIJevfubfSZTO2EAsDy5cuxc+dOKBQK/OEPf8DcuXPh6upqeG9xcTFCQkJw8uTJJtMffPABqqqqMHbsWGRmZuK5557D4cOH4eLiYnK9UVFR8PLywokTJ3Dt2jX4+Phg/vz5sLV9dAfaHt1PZgG2trbQ6/Utem3Hjh2xatUqHDhwAJ988glWrVqFyspKAICvry/27NmD2NhYZGRkIC4uDp06dTI5v7FDhw5BrVYDADp16oQdO3bA3d0dKSkpqKiowKpVq5CUlITq6mrD9h7E09MTCoUC48ePxyeffIIRI0bA29sbR44cweDBg/Hkk08CAHx8fODi4mI4p9Va5tb7wgsvGHYMnnvuOdy9e/c3bY/avxMnTqBv376G8GpsyJAhDxzVaKykpASHDh3Chg0bkJWVhffeew9Lly41LM/Pz8enn36K7777Do6OjkhPT0dQUBD8/f0xefJkTJw40WyNJSUlSE9PR2pqKlJTU/HCCy8gMzMTvr6+2LhxY5PX398J/eabb/DVV1+hqKgIFRUV2Lp1Kw4ePIgtW7YgKysLvXr1QkJCQgtaqcGiRYvQoUMHbNu2DQqFwjDf3HqvXLmCDRs2YPv27cjNzcXRo0dbvE0ZMcCa0a9fP1y+fBlardZo/o0bNzBlyhRUVVUZ5l2/fh2vvvoqSkpKMGDAAMNe2P317N27F2q1GiUlJRg/fjwKCwtNzm9MqVTCxsbGMH316lVotVpMmjQJBw4cwNNPP4133nkH3bp1Q3O3tezcuTO2bduGWbNmQaFQYMaMGfjqq6+g1+uN1g8AQgjU1dU1WUfj9d+/iMUUc+vt0KGDYb6NjU2ztdPjwdy/gd///vdYvHix4Yg/PT0dOp3OsNwSO0VBQUEAYNjx8vX1BQD06NHjgesztROam5uL0NBQODk5AQCio6Nx5MgRs/3GHHPr9fPzg62tLVQqFdzd3R/5HUMGWDPc3NwQEhKCDz/80BBiWq0WGo0GXbp0MfoSLiwshIuLC95++20MHToUOTk5ABquYlyyZAlWrFiBwMBAzJ49Gz179sTFixdNzm/Mx8cHW7duBdBwYjcmJgY///wzTp8+jZkzZ+Lll1/G9evXceXKlWaPFnNycjB58mT0798f06ZNw6uvvorCwkL4+Pjghx9+wNWrVwEAhw8fxrVr1/DHP/7R6P2Nj57+9re/oayszLBMoVA0CbyWrpceH97e3jh9+jR++eUXAA07QfeH3o4cOWL0b+NBO0tnzpyBWq2GVqvFkCFD8MYbbxitvyU7Rf84v7a21mi5vb290bSdnV2zn8nUTug/7sDp9fomfcRcLQ9ibr2P244hz4GZMW/ePKxYsQLh4eFQKBSoqalBYGAgpk2bZvS6IUOGYMuWLXjllVdgY2ODF154AS4uLigqKkJMTAwSEhIQHBwMe3t79O7dG6NHj8bdu3cfOH/Hjh2G9SYmJkKj0SAkJARCCLz55pvo27cvpkyZgnHjxsHJyQlubm7w9vZGUVGRYc/xH7300kvIzc1FcHAwnJyc4OzsjPnz56N79+6YN28epk6divr6enTo0AGrVq1qMpQ5c+ZMaDQaZGRkoE+fPujTp49hWVBQEOLj46HRaAzzevbs2aL10uOjX79+GDRoEBISEpCcnIybN29i2rRp8PDwwPnz5w0XId3fWerXr5/RztKxY8fQt29fvPbaa6ivr0dSUhLq6+vNbrfxDlbXrl2Rm5sLIQR0Oh1ycnIwcODAVn+mJUuWQAiB+Ph4BAQE4Pz587h48SJ8fX2xdetWjB49Gk5OTtiwYQMGDhxoFJCdO3dGbW0t/ud//gc9e/bEzp07DcuUSiXq6+shhDAKrJas93Fiw8epENHDUl9fj7Vr12L79u0AgKqqKvTp0wc///wzYmJiMG7cOBw5cgQajQYODg7o06cPzp07h6SkJPzud7/D9OnTcefOHej1evj5+SE9PR0//PAD1q1bZ/LijN27dyMlJQXh4eGIjIxEfHw8Ll26BDc3N/Tu3RtCCMNFHPevoL19+zZ8fHxw/vx5AMDGjRtx9OhRo3NuQMNFHAkJCbhx44ZhJzQpKQlKpRJpaWnYs2cP9Ho93N3dodFo8M///M9G2/niiy+wfv16uLi44JVXXsGKFStw8uRJ1NfXIyYmBrdv38ZXX32FwYMH4/Dhw+jSpUuL1gugyfSjiAFGRG3u3r17KCws5G8I6VdhgBERkZR4EQcREUmJAUZERFJ66Fch1tXV484d0z+4bY+6dnWSqmbWa32uru3nakr2KetjvdbXmj710I/AlEqF+Re1M7LVzHofLzK2n2w1s972iUOIREQkJQYYERFJiQFGRERSYoAREZGUGGBERCQlBhgREUmpRQF26tQpREVFNZm/b98+hIWFQa1W4+uvv7Z4cURERKaY/SHzmjVrsH379iZPUa2trcWiRYuwZcsWODo6IiIiAn5+fkaPzCYiIrIWs0dgPXr0QFpaWpP5ly5dQo8ePeDs7Ax7e3sMGDAA+fn5VimSiIjoH5k9AhsxYgSKi4ubzNdqtUYPJ+zYsaPhqcXmmLtlyD88ib6daB+3DmrpswPka+P20b5Ay9uYqLFuKzq3dQntVunb96yy3lbfC1GlUkGn0xmmdTpdi5+2W1ZWYeYV7efLrL0x33YN4cU2br2WtjERta1WX4Xo4eGBoqIilJeXo6amBvn5+ejfv78layN6ZPHCKKLf7lcfgWVlZaGyshJqtRoJCQl4/fXXIYRAWFgY3NzcrFEj0SOFF0YRWUaLAqx79+6GvcGQkBDDfH9/f/j7+1unMqJH1P0Lo/7jP/7DaH7jC6MAGC6MGjlypNl1yjikKVvNstXbnlir7R7688CIHnfWuDCqJeft2pOWnadtP2Srt72x1nll3omDqJ34LRdGET2OGGBE7QQvjCL6dTiESNTGeGEUUevYCPHwf7Zpbjy0WzcOm5hSWmqZ34GxjU1raRu3J2b/3vyRrUkt+ZFti/oU29iklrbxr8UhRCIikhIDjIiIpMQAIyIiKTHAiIhISgwwIiKSEgOMiIikxAAjIiIpMcCIiEhKDDAiIpISA4yIiKTEACMiIikxwIiISEoMMCIikhIDjIiIpMQAIyIiKTHAiIhISgwwIiKSEgOMiIikxAAjIiIpMcCIiEhKDDAiIpISA4yIiKTEACMiIikxwIiISEoMMCIikhIDjIiIpMQAIyIiKTHAiIhISmYDTK/XIzExEWq1GlFRUSgqKjJavn37dowbNw5hYWH47//+b6sVSkRE1JjS3Auys7NRU1ODjIwMFBQUICUlBStXrjQsX7x4MXbs2AEnJyeMHj0ao0ePhrOzs1WLJiIiMhtgx48fh6+vLwDAy8sLhYWFRst79+6NiooKKJVKCCFgY2NjnUqJiIgaMRtgWq0WKpXKMK1QKFBXVwelsuGtvXr1QlhYGBwdHREUFITOnTub3aira6ffUPLjraVtxzZuPbYdkRzMBphKpYJOpzNM6/V6Q3idO3cO+/fvx969e+Hk5IT4+Hh8++23GDlyZLPrLCurMLNVfoGYYr7tGr6A2cat19I2bi29Xg+NRoPz58/D3t4eCxYsgLu7u2H59u3bsW7dOtja2iIsLAyRkZGt3hbRo8zsRRze3t7Izc0FABQUFMDT09OwrFOnTujQoQMcHBygUCjg4uKCe/fuWa9aokdA4/PKcXFxSElJMVq+ePFirFu3Dps2bcK6detw9+7dNqqUqH0zewQWFBSEvLw8hIeHQwiBhQsXIisrC5WVlVCr1VCr1YiMjISdnR169OiBcePGPYy6iaTF88pElmE2wGxtbZGcnGw0z8PDw/D/ERERiIiIsHxlRI8onlduX3he2fqs1XZmA4yILKttziuTKZY7r0ymWOu8Mu/EQfSQ8bwykWXwCIzoIeN5ZSLLYIARPWQ8r0xkGRxCJCIiKTHAiIhISgwwIiKSEgOMiIikxAAjIiIpMcCIiEhKDDAiIpISA4yIiKTEACMiIikxwIiISEoMMCIikhIDjIiIpMQAIyIiKTHAiIhISgwwIiKSEgOMiIikxAAjIiIpMcCIiEhKDDAiIpISA4yIiKTEACMiIikxwIiISEoMMCIikhIDjIiIpMQAIyIiKTHAiIhISgwwIiKSEgOMiIikxAAjIiIpMcCIiEhKSnMv0Ov10Gg0OH/+POzt7bFgwQK4u7sblv/4449ISUmBEAKurq5ITU2Fg4ODVYsmIiIyewSWnZ2NmpoaZGRkIC4uDikpKYZlQgjMnTsXixYtwqZNm+Dr64uSkhKrFkxERAS04Ajs+PHj8PX1BQB4eXmhsLDQsOynn35Cly5dsH79ely4cAHDhg3D008/bb1qiYiI/j+zAabVaqFSqQzTCoUCdXV1UCqVuHPnDk6ePIm5c+fC3d0db731Fvr27QsfH59m1+nq2um3V/6YamnbsY1bz9ptx2F5IsswG2AqlQo6nc4wrdfroVQ2vK1Lly5wd3dHz549AQC+vr4oLCw0G2BlZRVmtsovX1PMt13DFzDbuPVa2sat1XhYvqCgACkpKVi5ciWA/xuWX7p0Kdzd3bF582aUlJRwZIPoAcwGmLe3N3JycjBq1CgUFBTA09PTsOzJJ5+ETqdDUVER3N3dkZ+fj3/913+1asFEsrPGsDyPuFuPoxrWZ622MxtgQUFByMvLQ3h4OIQQWLhwIbKyslBZWQm1Wo2PPvoIcXFxEEKgf//+GD58uFUKJXpUWGNYviVHjfRglhvVIFOsNaphNsBsbW2RnJxsNM/Dw8Pw/z4+PtiyZcuv3jDR48oaw/JEjyP+kJnoIfP29kZubi4ANDssDwD5+fno1atXm9RJ1N6ZPQIjIsvisDyRZTDAiB4yDssTWQaHEImISEoMMCIikhIDjIiIpMQAIyIiKTHAiIhISgwwIiKSEgOMiIikxAAjIiIpMcCIiEhKDDAiIpISA4yIiKTEACMiIikxwIiISEoMMCIikhIDjIiIpMQAIyIiKTHAiIhISgwwIiKSEgOMiIikxAAjIiIpMcCIiEhKDDAiIpISA4yIiKTEACMiIikxwIiISEoMMCIikhIDjIiIpMQAIyIiKTHAiIhISgwwIiKSktkA0+v1SExMhFqtRlRUFIqKih74urlz52LJkiUWL5CIiOhBzAZYdnY2ampqkJGRgbi4OKSkpDR5TXp6Oi5cuGCVAomIiB7EbIAdP34cvr6+AAAvLy8UFhYaLT958iROnToFtVptnQqJHjEc1SCyDKW5F2i1WqhUKsO0QqFAXV0dlEolSktLsWzZMixbtgzffvttizfq6tqpddVSi9uObdx61m67xqMaBQUFSElJwcqVK41ec39UY+DAgVathUhmZgNMpVJBp9MZpvV6PZTKhrd99913uHPnDqZMmYKysjJUVVXh6aefRmhoaLPrLCurMLNVfvmaYr7tGr6A2cat19I2bq1fM6px+fLlVm+H6FFnNsC8vb2Rk5ODUaNGoaCgAJ6enoZl0dHRiI6OBgBkZmbi8uXLZsOL6HHHUY32haMa1mettjMbYEFBQcjLy0N4eDiEEFi4cCGysrJQWVnJ815ErdA2oxpkiuVGNcgUa41qmA0wW1tbJCcnG83z8PBo8joeeRG1DEc1iCzDbIARkWVxVIPIMhhgRA8ZRzWILIO3kiIiIikxwIiISEoMMCIikhIDjIiIpMQAIyIiKTHAiIhISgwwIiKSEgOMiIikxAAjIiIpMcCIiEhKDDAiIpISA4yIiKTEACMiIikxwIiISEoMMCIikhIDjIiIpMQAIyIiKTHAiIhISgwwIiKSEgOMiIikxAAjIiIpMcCIiEhKDDAiIpISA4yIiKTEACMiIikxwIiISEoMMCIikhIDjIiIpMQAIyIiKTHAiIhISgwwIiKSEgOMiIikpDT3Ar1eD41Gg/Pnz8Pe3h4LFiyAu7u7YfmOHTuwfv16KBQKeHp6QqPRwNaWuUhERNZlNmmys7NRU1ODjIwMxMXFISUlxbCsqqoKn3zyCb788kukp6dDq9UiJyfHqgUTyU6v1yMxMRFqtRpRUVEoKioyWr5jxw6MHz8e4eHhSExMhF6vb6NKido3s0dgx48fh6+vLwDAy8sLhYWFhmX29vZIT0+Ho6MjAKCurg4ODg5mN+rq2qm19T72Wtp2bOPWs3bbNd4pLCgoQEpKClauXAng/3YKs7Ky4OjoiPfffx85OTkICAiwak1EMjIbYFqtFiqVyjCtUChQV1cHpVIJW1tb/NM//RMAYMOGDaisrMSQIUPMbrSsrMLMK/jla4r5tmv4AmYbt15L27i1uFPYvnCn0Pqs1XZmA0ylUkGn0xmm9Xo9lEql0XRqaip++uknpKWlwcbGxiqFEj0q2mankEyx3E4hmWKtnUKz58C8vb2Rm5sLACgoKICnp6fR8sTERFRXV2PFihWGvUYiMq0lO4Uff/wx8vLyuFNI1AyzR2BBQUHIy8tDeHg4hBBYuHAhsrKyUFlZib59+2LLli3405/+hJiYGABAdHQ0goKCrF44kay8vb2Rk5ODUaNGmdwptLe3x4oVK3hFL1EzzAaYra0tkpOTjeZ5eHgY/v/cuXOWr4roEcadQiLLMBtgRGRZ3CkksgyOTxARkZQYYEREJCUGGBERSYkBRkREUmKAERGRlBhgREQkJQYYERFJiQFGRERSYoAREZGUGGBERCQlBhgREUmJAUZERFJigBERkZQYYEREJCUGGBERSYkBRkREUmKAERGRlBhgREQkJQYYERFJiQFGRERSYoAREZGUGGBERCQlBhgREUmJAUZERFJigBERkZQYYEREJCUGGBERSYkBRkREUmKAERGRlBhgREQkJQYYERFJiQFGRERSMhtger0eiYmJUKvViIqKQlFRkdHyffv2ISwsDGq1Gl9//bXVCiV6VLBPEVmG2QDLzs5GTU0NMjIyEBcXh5SUFMOy2tpaLFq0CGvXrsWGDRuQkZGBsrIyqxZMJDv2KSLLMBtgx48fh6+vLwDAy8sLhYWFhmWXLl1Cjx494OzsDHt7ewwYMAD5+fnWq5boEcA+RWQZSnMv0Gq1UKlUhmmFQoG6ujoolUpotVp06tTJsKxjx47QarVmN+rq2qnZ5UKYXcVjrPm2u49t/Fu0rI1bq0361Dz+wX8rtnH7Y/YITKVSQafTGab1ej2USuUDl+l0OqPOR0RNsU8RWYbZAPP29kZubi4AoKCgAJ6enoZlHh4eKCoqQnl5OWpqapCfn4/+/ftbr1qiRwD7FJFl2AjR/GCSXq+HRqPBhQsXIITAwoULcfbsWVRWVkKtVmPfvn1Yvnw5hBAICwvDxIkTH1btRFJinyKyDLMBRkRE1B7xh8xERCQlBhgREUmJAUZERFKyeoBVVVVh2rRpiIyMRGxsLG7fvt3kNQsWLEBoaCiioqIQFRWFiooKa5fVhGy39zFX77p16zB69GhDm16+fLmNKjV26tQpREVFNZnf3tq3MVM1t1Ubs09Zh6x9CpCvX1msTwkrW7t2rVi6dKkQQogdO3aI+fPnN3lNeHi4uHXrlrVLadbu3bvFrFmzhBBCnDx5Urz11luGZTU1NSIwMFCUl5eL6upqERoaKkpLS9uqVCFE8/UKIURcXJw4ffp0W5Rm0urVq0VwcLAYP3680fz22L73mapZiLZrY/Yp65CxTwkhX7+yZJ+y+hFY49vmvPTSSzh8+LDRcr1ej6KiIiQmJiI8PBxbtmyxdkkPJNvtfZqrFwDOnDmD1atXIyIiAp999llblNhEjx49kJaW1mR+e2zf+0zVDLRdG7NPWYeMfQqQr19Zsk+ZvZXUr7F582asX7/eaN4TTzxhuJNAx44dmwxlVFZWYtKkSXjttddQX1+P6Oho9O3bF88884wlSzPLGrf3sabm6gWA0aNHIzIyEiqVClOnTkVOTg78/PzaqlwAwIgRI1BcXNxkfnts3/tM1Qw8nDZmn3p4ZOxTgHz9ypJ9yqJHYOPHj8eOHTuM/uvUqZPh1jg6nQ6dO3c2eo+joyOio6Ph6OgIlUqFwYMH49y5c5Ysq0Vku71Pc/UKIRATEwMXFxfY29tj2LBhOHv2bFuValZ7bF9zHlYbs089PI9SnwLaZxs3pzVtbPUhRG9vbxw4cAAAkJubiwEDBhgt//nnnxEZGYn6+nrU1tbixIkT6NOnj7XLemCdMt3ep7l6tVotgoODodPpIITA3/72N/Tt27etSjWrPbavOW3ZxuxT1vEo9SmgfbZxc1rTxhYdQnyQiIgIzJo1CxEREbCzs8Of//xnAA1Xm/To0QMBAQEICQnBhAkTYGdnh7Fjx6JXr17WLquJoKAg5OXlITw83HB7n6ysLMPtfRISEvD6668bbu/j5ub20Gv8NfW+9957iI6Ohr29PXx8fDBs2LA2rfdB2nP7mtIe2ph9qm3qlaFPAfL1q9/SxryVFBERSYk/ZCYiIikxwIiISEoMMCIikhIDjIiIpMQAIyIiKTHAmhEREYGdO3cazausrMSgQYMeeAPV5ixdurRd3MaltRISEpCZmfmb1xMbG4sbN25YoCKSUXFxMfz9/ZvM7927NwBg7969+PTTTwEA/v7+KC4uRmZmJhISElq9zYqKCrzzzjutfr813f+Mv8WNGzcQGxtroYrkwgBrRlhYGLKysozm7dmzB4MGDYKLi8uvWtexY8dQX19vyfKktGbNmnb5WxRqHwICAvDuu+9adJ13797F3//+d4uusz1xc3PDmjVr2rqMNmH1HzLLbOTIkVi8eDHKy8vRpUsXAMD27dsRExMDnU6H5ORkXLx4EfX19YiNjUVwcDCqq6uRlJSE48ePw87ODm+//TZqampQWFiIOXPmYNmyZbC3t0diYiLKy8vh5OSE2bNno1+/fkhISEB5eTmKiooQHx9vtKf68ccfIy8vD7a2tggMDMTUqVORlpaG//3f/8WlS5dw584dqNVqvPHGG9Bqtfjwww9x48YNlJaWwsfHBx999BEAYMmSJcjOzoZCoYBarUZMTAyKioqg0WhQXl6ODh06YO7cuXjuueeatMf+/fuxceNG1NbW4t///d8xatQok9u6ceMGZs6cicrKStja2mLOnDnw8vKCv78/vvzyS2i1WiQmJqKurg4ODg5YtGgRnnrqqYfxZ6V2LDMzE0ePHkVKSorR/KKiIkycOBF3797F8OHDERcXBxsbG/z1r3/F+vXrodfr0adPH8ybNw8ODg4YPHgw+vbti7KyMri6uqK0tBTvvPMOli9fblinVqvF+++/j5s3bwIA3nnnHQQEBCAqKgrPPPMM8vPzUV1djQ8//BBDhw7FhQsXMH/+fFRWVuL27duYMmUKIiIiUF5ejtmzZ+Py5cuwt7dHQkICfHx8kJubi6VLl6Kurg7du3fH/Pnz0bVr1yafedmyZTh37hwcHByQlJSEZ555xuS2Dh8+jNTUVACAs7Mz/vznP6OyshLR0dHYt28fsrKy8Pnnn0OhUKB79+5ITU2Fg4ODFf9ibaz1N8V/PCQkJIhNmzYJIYS4fv26GD58uKirqxOpqali/fr1QgghKioqxOjRo8WVK1fEmjVrxLvvvivq6+tFaWmpGDVqlKiurhaTJk0SR44cEUIIERYWJnbv3i2EaHhsw/Dhw0V1dbWYNWuW4XEOjRUXF4tRo0YJIYSorKwU7777rqiqqhJLly4VwcHBQqvVinv37onAwEBRWFgosrKyxIoVK4QQQlRXV4vAwEBx+vRpsWvXLhEeHi6qq6uFVqsVY8aMEaWlpUKtVoszZ84IIYS4ePGiePnll5vUMGvWLPHGG2+I2tpacf36dTFkyBBRWlpqcltpaWlizZo1QgghDhw4ID7//HMhhBB+fn7i6tWrIiEhQezatUsIIURmZqb45ptvfvsfi9q9q1evij59+ogxY8YY/efp6SmEEGLr1q2GPnD/38rWrVvF0KFDxa1bt0R1dbUIDw8X33//vbhw4YKIiIgQVVVVQgghlixZIpYvXy6EEMLT09PQ365evSr8/Pya1JKZmSk0Go0QQoizZ8+KlJQUIYQQkyZNEgkJCYb5Q4YMEdXV1WLBggXi0KFDQgghrly5Iry8vIQQQmg0GsN7z507JyZMmCBu3bolxowZI8rLy4UQQmzatEl8+OGHTWrw8/Mz9J/9+/eLsWPHCiGEyW1NmjRJnDp1SgjR8FiSgwcPGn0+f39/cfPmTSGEECkpKeLs2bMt/+NIiEdgZoSGhuLTTz9FeHg4srKyMGbMGCgUChw6dAhVVVXYunUrgIZzYxcvXsSxY8cwYcIE2NrawtXVtck5NJ1OhytXruDll18G0PDYBmdnZ8OD2/r169ekBjc3Nzg4OCA8PBx+fn6YOXOmYa8qODgYHTt2BNAwnn7kyBG8/vrr+PHHH/HFF1/g8uXLKC8vR2VlJY4dO4aRI0fC3t4e9vb22LZtG3Q6HQoLC/HBBx8YtldZWYk7d+402VscN24clEol3Nzc4OXlhVOnTiE4OPiB2/Lx8cG0adPw97//HcOGDcOkSZOM1jVs2DAkJyfj4MGD8Pf3bxd39aaHo1u3bti2bZvRvPvnwEzx9/c3DNuPHDkSR48exfXr11FUVIQJEyYAAGpra41GDv74xz82u87+/fvjP//zP3Hjxg0MHz7c6DzZ/XU+++yzcHV1xfnz55GQkICDBw/is88+w4ULF1BZWQmg4fTAkiVLDJ8jIyMDOTk5uHbtGqKjowE03BjY2dn5gXWMHz8eQEOfiI+Px71790xuKyAgAFOnTkVgYCACAgIwZMgQo3Nofn5+iIiIQGBgIEaMGIFnn3222TaQHQPMjIEDB6KsrAzXrl3D9u3bsWzZMgAN/yBTU1MNN0m9efMmnJ2dsXXrVtjY2BjeX1RUhN/97neGafGAO3cJIQznxzp06NBkuVKpxObNm3H06FHk5uYiPDwcGzZsANDwyIf79Ho9FAoFNmzYgN27d2PChAl48cUXceHCBQghoFQqjWorLi42PCuo8RfK9evXDUOmjf3jtuzs7Exua8CAAdi5cyf279+PXbt24ZtvvsG6desM73/llVfQv39/5OTk4IsvvsD+/fuxYMECE38Fetzdvys88H93ia+vr8fIkSMxZ84cAA07h43PMz+oLzX21FNP4dtvv8XBgweRk5ODtWvXYteuXQCa/ltXKpWYMWMGOnfuDD8/P4waNQo7duww1Na4X126dAn19fXw9vbGqlWrAADV1dVGd4ZvrPG27vdTU9uaPHky/Pz8kJOTg9TUVPz4448ICQkxvH/OnDk4d+4cDhw4gPj4eEydOhVjx45tth1kxos4WuDVV1/FypUr4ezsjB49egAABg8ejE2bNgEASktLMWbMGFy7dg0DBw7Erl27IITArVu3MGnSJNTU1EChUKC+vh4qlQrdu3fHnj17ADTc9frmzZvN3mz17NmzmDRpEgYOHIhZs2bBw8MDP/30EwAgOzsbNTU1uHv3LnJycjB06FDk5eVBrVZjzJgxqK6uxrlz56DX6zFw4EDs2bMHtbW1+OWXX/DGG2/g5s2beOqppwwBlpeXh4kTJz6wjp07d0IIgZKSEhQWFuL55583ua3Fixdj+/btGDduHBITE5s8FmHGjBk4ffo0wsPD8e6777b7R1NQ2zpw4ADu3buH6upq7Nq1Cy+++CIGDRqE77//Hrdu3YIQAhqNpsmz04CGgKmrq2syf+PGjUhLS8PIkSMxb9483L592/C8rPtBdvr0ady7dw+enp7Iy8vD9OnTERgYaLhrfX19Pf70pz8ZRlouXbqE2NhY9OvXDwUFBYZ+umLFCixevPiBn+3+hWLff/89PDw84OTkZHJb48ePh06nw+TJkzF58mSjflNXV4eXX34ZXbt2xZtvvomxY8c+0hevADwCa5HQ0FD4+/sbLoQAgKlTp0Kj0SA4OBj19fWIj49Hjx49EBkZiQULFmDMmDEAgLlz50KlUsHX1xfz5s3Dxx9/jNTUVGg0GqSlpcHOzg5paWmwt7c3uf3nnnsOXl5eCA4OhqOjI7y9vfHSSy/hzJkzcHBwQGRkJLRaLd5880307NkTMTEx0Gg0WL16NVQqFfr374/i4mKMHz8ehYWFCA0NhV6vR3R0NP7whz8Y6vn8889hZ2eHv/zlL0Z7lPc5OTkhNDQUdXV1SE5OhouLi8ltRUVFIS4uDpmZmVAoFPj444+N1vXWW29h9uzZWL58Oezs7KDRaCzzx6JH0tNPP40pU6bg3r17CA4OxtChQwE09MOYmBjo9Xo8++yzmDJlSpP3PvHEE/iXf/kXREVFGUYugIYd0/fffx8hISFQKBSIj483PFvt6tWrGDduHADgL3/5CxQKBaZNm4bIyEg4ODjgmWeewe9//3sUFxdj+vTpmDNnDsaMGQOlUonFixejW7duWLhwIWbMmAG9Xg83NzfDxRf/6Oeff8bYsWPRsWNHw8Urprb1/vvvIyEhAUqlEk5OTkajFkqlEtOnT8e//du/wcHBAU888USTi2EeNbwbvcTuP5Z72rRpbVwJ0aMjKioKU6dOxaBBg9q6FDKDQ4hERCQlHoEREZGUeARGRERSYoAREZGUGGBERCQlBhgREUmJAUZERFL6f0S5t5leVTSMAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 504x288 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "num_qubits = int(np.log2(len(vector)))\n",
    "fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(7, 4))\n",
    "\n",
    "ax1.bar(np.arange(0, 2 ** num_qubits), classical_solution, color=\"blue\")\n",
    "ax1.set_xlim(-0.5, 2 ** num_qubits - 0.5)\n",
    "ax1.set_xlabel(\"Vector space basis\")\n",
    "ax1.set_title(\"Classical solution\")\n",
    "\n",
    "ax2.bar(np.arange(0, 2 ** num_qubits), quantum_solution, color=\"green\")\n",
    "ax2.set_xlim(-0.5, 2 ** num_qubits - 0.5)\n",
    "ax2.set_xlabel(\"Hilbert space basis\")\n",
    "ax2.set_title(\"Quantum solution\")\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Therefore with that we have reached the end of this VQLS demonstrator. To recap, (i) we learnt the conversion of PDEs to linear system of equations, (ii) along with the intution behind the algorithm, we have learnt the rigorous theory of the algorithm, (iii) global cost function was studied step by step and leant how to calculate the expectation values on the quantum hardware, (iv) finally we looked into local cost function and it's theory. Although the VQLS algorithm was succesfully demonstrated, we still have a long way to go interms of efficiency of the algorithm. \n",
    "<br>\n",
    "\n",
    "Following are the some of the challenges we listed:"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Challanges:\n",
    "\n",
    "- Find an efficient matrix decomposition of $A$ \n",
    "- Initialize the state $\\vert b \\rangle $ efficiently $\\vec{b} \\neq \\vec{0}$\n",
    "- Compare local cost function vs. global cost function (barren-plateaus vs. local minima). \n",
    "- Read out the entire vector $\\vert \\psi_{\\mathrm{optimal}} \\rangle = A \\vert x_{\\mathrm{optimal}} \\rangle$ with only a finite amount of shots (not an exponential scaling)\n"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Copyright 2023 Dr. Bharadwaj Chowdary Mummaneni, Fraunhofer IAO**\n",
    "\n",
    "Licensed under the Apache License, Version 2.0 (the \"License\");\n",
    "you may not use this file except in compliance with the License.\n",
    "You may obtain a copy of the License at\n",
    "\n",
    "    http://www.apache.org/licenses/LICENSE-2.0\n",
    "\n",
    "Unless required by applicable law or agreed to in writing, software\n",
    "distributed under the License is distributed on an \"AS IS\" BASIS,\n",
    "WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.\n",
    "See the License for the specific language governing permissions and\n",
    "limitations under the License."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "base",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.12 (main, Apr  4 2022, 05:22:27) [MSC v.1916 64 bit (AMD64)]"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "137a1295176ce66a7ded57fa2bf39d327238a1a5a04200a5a558013d7d4f9d1d"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}