vqe_dcs_qiskit_optimized.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "adfa7d5d",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Installation of the requirements\n",
    "#!python -m pip install -r requirements.txt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "74dcd404",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "SparsePauliOp(['IIII', 'IIIZ', 'IIZZ', 'IIZI', 'IZII', 'IZIZ', 'ZXIX', 'IXZX', 'ZXZX', 'IXIX', 'ZZZI', 'ZZZZ', 'IZZZ', 'ZZIZ', 'ZIZI'],\n",
      "              coeffs=[-0.81054798+0.j,  0.17218393+0.j, -0.22575349+0.j,  0.12091263+0.j,\n",
      "  0.17218393+0.j,  0.16892754+0.j,  0.0452328 +0.j, -0.0452328 +0.j,\n",
      " -0.0452328 +0.j,  0.0452328 +0.j, -0.22575349+0.j,  0.16614543+0.j,\n",
      "  0.16614543+0.j,  0.17464343+0.j,  0.12091263+0.j])\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/renata/.local/lib/python3.8/site-packages/pyscf/lib/misc.py:46: H5pyDeprecationWarning: Using default_file_mode other than 'r' is deprecated. Pass the mode to h5py.File() instead.\n",
      "  h5py.get_config().default_file_mode = 'a'\n"
     ]
    }
   ],
   "source": [
    "'''\n",
    "(C) Renata Wong 2023\n",
    "\n",
    "Qiskit code for testing fidelity of derandomised classical shadow on the ground state energy of molecules. \n",
    "This notebook implements an optimization: since the derandomized Hamiltonian may contan very few terms,\n",
    "instead of generating a quantum circuit for each and measuring once, we generate a single circuit and specify \n",
    "a shot number that matches the number of occurrences of a derandomized operator. This speeds up the computation \n",
    "significantly. \n",
    "\n",
    "Procedure:\n",
    "1. Derandomize the molecule-in-question's Hamiltonian.\n",
    "2. Choose a variational ansatz with initial parameters selected at random.\n",
    "3. Apply the derandomized Hamiltonian as basis change operators to the ansatz.\n",
    "4. Measure the ansatz in the Pauli Z basis and store the results as a shadow.\n",
    "5. Obtain the expectation value of the molecular Hamiltonian from the shadow.\n",
    "6. Optimize for minimum Hamiltonian expectation value. \n",
    "7. Feed the calculated angles/parameters back to the ansatz.\n",
    "8. Repeat steps 3-7 till the optimization is completed. \n",
    "9. Output the minimized expectation value of the molecular Hamiltonian and the mean-square-root-error. \n",
    "\n",
    "Note: Below we perform calculations on the molecular Hamiltonian of H_2.\n",
    "To perform calculations on other molecules, you will need to specify their geometry, charge and spin \n",
    "to replace the values in the driver. \n",
    "'''\n",
    "\n",
    "from qiskit import QuantumCircuit, execute\n",
    "from qiskit_aer import QasmSimulator\n",
    "\n",
    "from qiskit_nature.units import DistanceUnit\n",
    "from qiskit_nature.second_q.drivers import PySCFDriver\n",
    "from qiskit_nature.second_q.mappers import BravyiKitaevMapper\n",
    "\n",
    "from qiskit.algorithms.minimum_eigensolvers import NumPyMinimumEigensolver\n",
    "from qiskit_nature.second_q.algorithms import GroundStateEigensolver\n",
    "\n",
    "from qiskit.circuit.library import EfficientSU2\n",
    "from qiskit.algorithms.optimizers import SLSQP, COBYLA, SPSA\n",
    "\n",
    "import time\n",
    "import numpy as np\n",
    "from functools import partial\n",
    "import matplotlib.pyplot as plt\n",
    "from collections import Counter\n",
    "\n",
    "from predicting_quantum_properties.prediction_shadow import estimate_exp\n",
    "from modified_derandomization import modified_derandomized_classical_shadow\n",
    "\n",
    "# handle deprecation issues\n",
    "import qiskit_nature\n",
    "qiskit_nature.settings.use_pauli_sum_op = False\n",
    "import h5py\n",
    "H5PY_DEFAULT_READONLY=1\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "# classically obtained ground state energy\n",
    "EXPECTED_EIGENVALUE = -1.86\n",
    "\n",
    "# specifying the geometry of the molecule in question\n",
    "driver = PySCFDriver(\n",
    "    atom=\"H 0 0 0; H 0 0 0.735\",\n",
    "    basis=\"sto3g\",\n",
    "    charge=0,\n",
    "    spin=0,\n",
    "    unit=DistanceUnit.ANGSTROM,\n",
    ")\n",
    "\n",
    "\n",
    "problem = driver.run()\n",
    "hamiltonian = problem.hamiltonian\n",
    "\n",
    "# electronic Hamiltonian of the system\n",
    "second_q_op = hamiltonian.second_q_op()\n",
    "\n",
    "# The Bravyi-Kitaev repserentation of the Fermionic Hamiltonian\n",
    "mapper = BravyiKitaevMapper()\n",
    "bkencoded_hamiltonian = mapper.map(second_q_op)\n",
    "\n",
    "print(bkencoded_hamiltonian)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "0b205cb1",
   "metadata": {},
   "outputs": [],
   "source": [
    "def process_hamiltonian(hamiltonian, derandomize = False):\n",
    "\n",
    "    hamiltonian_observables = []\n",
    "    hamiltonian_coefficients = []\n",
    "\n",
    "    for observable in hamiltonian.paulis:\n",
    "        op_list = []\n",
    "        for op_index, pauli_op in enumerate(observable):\n",
    "            pauli_op = str(pauli_op)\n",
    "            if pauli_op == 'I' or pauli_op == 'X' or pauli_op == 'Y' or pauli_op == 'Z':\n",
    "                op_list.append((pauli_op, op_index))\n",
    "        hamiltonian_observables.append(op_list)\n",
    "    \n",
    "    hamiltonian_coefficients = hamiltonian.coeffs.real\n",
    "    \n",
    "    system_size = len(hamiltonian_observables[0])\n",
    "\n",
    "    # removing all occurrences of Pauli-I, for all-Pauli-I there is an empty list left\n",
    "    # these observables are needed for estimate_exp()\n",
    "    observables_xyze = []\n",
    "\n",
    "    for observable in hamiltonian_observables:\n",
    "        XYZE = []\n",
    "        for pauli in observable:\n",
    "            if pauli[0] != 'I':\n",
    "                XYZE.append(pauli)\n",
    "        observables_xyze.append(XYZE)   \n",
    "    \n",
    "    # derandomisation procedure requires that coefficients are non-negative\n",
    "    if derandomize == True:\n",
    "        absolute_coefficients = [abs(coeffcient) for coeffcient in hamiltonian_coefficients]\n",
    "\n",
    "        # removing the empty list as well\n",
    "        # these observables are needed for derandomisation procedure\n",
    "        observables_xyz = []\n",
    "        for idx, observable in enumerate(observables_xyze):\n",
    "            if observable:\n",
    "                observables_xyz.append(observable)\n",
    "            else:\n",
    "                absolute_coefficients.pop(idx)\n",
    "        \n",
    "        return system_size, observables_xyze, observables_xyz, hamiltonian_coefficients, absolute_coefficients\n",
    "    \n",
    "    return system_size, observables_xyze, hamiltonian_coefficients\n",
    "\n",
    "\n",
    "\n",
    "# process the Hamiltonian to obtain properly formatted data\n",
    "hamiltonian_data = process_hamiltonian(bkencoded_hamiltonian, derandomize = True)\n",
    "system_size, observables_xyze, observables_xyz, hamiltonian_coefficients, absolute_coefficients = hamiltonian_data\n",
    "\n",
    "#print('HAMILTONIAN\\n', observables_xyze)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "265d3af8",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAbsAAADWCAYAAABIdCwuAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAvXElEQVR4nO3deVxUZf//8dfMALILiEqCC8hSopDiiguYlubXr+JtWmreaX7TwF+5ttyhlWHaz6y4v7/b7G4x7rtU7iRLM80WAffKShM3VFxSMfcFxQWG3x8TKAIzA8zMmXP6PB+PeSRnzlzzhs65PnOuuc45urKysjKEEEIIDdMrHUAIIYSwNyl2QgghNE+KnRBCCM2TYieEEELzpNgJIYTQPCl2QgghNE+KnRBCCM2TYieEEELzpNgJIYTQPCl2QgghNE+KnRBCCM2TYieEEELzpNgJIYTQPCl2QgghNE+KnRBCCM2TYieEEELzpNgJIYTQPCl2QgghNE+KnRBCCM2TYieEEELzpNgJIYTQPCl2QgghNE+KnRBCCM2TYieEEELzpNgJIYTQPCl2QgghNM9F6QDObt86uHxKmff2aQJR99XttUrlrk9m4TiyXdeObNfqJ8XOgsun4MIxpVPUnlpzC8dQ6/ah1txCeTKMKYQQQvOk2AkhhNA8GcYUqnf1Oly+DjrA1wPcXZVOJIRwNlLshOqUlcHhM7B5Pxw8Beeu3HpOBwT6QGQQdI+AZv6KxRRCOBEpdjYwbWEie45swWBwRa83EOQfysg+qSTEDlM6mllqzP37RcjcCofOVP98GXD6sumxaT9EB8OwzuDn6dCYmqDG7QPUm1vYlxQ7GxnVdyaj+s6gtLSEFZv/wdwlIwkPbk9wYLjS0cxSU+6tByDrRygxWv+aXcehYBWM7g5tgu2XTavUtH3cTq25hf3IBBUbMxhceLDLE5QaSzh4YrvScazm7Lk35UPm97UrdOWKb8L7uZAnU9brzNm3j5qoNbewPSl2Nnaz5AarNi8EICQwUuE01nPm3IfPmI7ozEkfZXrUxFgGH22Cs0W2zfZn4czbhzlqzS1sT4qdjSz57lWSZvox8AUPPlw7g6nD3iesWQwAx88cICU9jpslNwD4JOd1Mta+qGTcCuZyz1k8kq27V1Ws+1JGEtv2fe3QfCWlsGSL6bu4+rpeYvq+r8wWjdmIM2WpjmzXQiucutgZjUbmz59PREQE7u7uxMbGkpubS1RUFOPHj1c6XiUj+6TyedoFsl4+Q+e7B7DjQHbFc8GB4fRoN5TMdXMpPHeInO2ZjOyTqmDaW8zlTh6cTsbamRRfL2LDzuV4uTekY9QDDs338xE4dcl27e3/HQpO2669uvj9InzyAzz/CUxZAi8thzW/QtE1ZXNVR7ZroRVOXezGjRtHWloaEyZMYM2aNQwfPpwRI0ZQUFBAXFyc0vGq5ePpz9Rh7/P93i/ZnLeiYvnwxGfYumcVcxaPIHlQOm4uDRRMWVV1uf29mzCkxyQWrHiaJd/N5slBbzk816Z8dbRprT0n4PXVsGU/XLtpWnaxGNbuNC0/c1m5bObIdi3UzmmL3dKlS8nIyGDlypVMnz6d3r17k5qaSrdu3SgpKaFDhw5KR6yRr2cAQ3tOZdFXL2A0mmZUuBhcaRfWi6Li87QN7aFwwupVl7tfpzEcO51PUven8fUMcGieq9fhyFnbt7u3UJnhw6JrsGg9lBqrH5a9dA0+WO+8Q5uyXQs1c9piN2fOHPr3709CQkKl5eHh4bi6uhITYxp/P3z4MAkJCURGRtKuXTs2bNigRNwqhvScxLlLhXzz078BOHxyF7sOb6J9eF9Wf/+ewulqdmdugGaNwhWZsn3svH3avXqj8onojrL1INwsrfn7x7IyKLyg/DCrObJdC7XSlZU53+fIY8eO0bx5cz744AMef/zxSs+NGDGCvXv38ssvvwDQr18/Bg8eTEpKCps3b2bYsGEcOnQINzc3s++h0+msyjL/yWxiWyfW6fcoZzQamfZOAsmD0gkJjGTSgnjmjf8Wf5+mZl+342AO09/pXaf3tEXucvMyxzCgy/9Y9cm9Ppnv1KbXWO4fv6jSMnMzLs2ZvLjyz5++2ptje3Lq1lgdDU3NJjiqFzp9zZ8xy8qMbPviNTZ/Yt/vvmS7Vm67FrZlbQlzyiO7Y8dMJ0QFBQVVWl5cXExubm7FEOaZM2fYuHEj48aNAyA+Pp5mzZqRnZ2NM/liy0IiguOIDInD092HMf3SeHvlZKVjOT2d3qDKtmticGlgttABUFaG3mD+g5qzkO1aqIlTHtkdOHCAiIgI3nrrLSZPnlyxfNasWbz88sssWLCAlJQUfv75Zx566CEKCgoq1hk+fDh9+/a12WzNbZnK3T/LLwQ6PlK31yqVuz6Z77TjKHxo5ah0+RHfnUdwNZn2IDR38Fc1mVvh+4OWT6MY0RW6tLZvFtmua8eW27VQhlNeLiwsLIyYmBjmzJlDQEAAwcHBZGVlsXr1agCnnYkpbCvYThdxNujhrob2aduc7hGm7+3MaeAC7Vs6Jo8QfyZOOYyp1+tZtmwZ0dHRJCcnM3bsWAIDA5k4cSIGg6FickqLFi34/fffuX79esVrDx06RMuW0ltoQSNv8HW3fbvNA8DF8aOYNG8EPWu4iEf5N8jDOoObU34EFULdnLLYAURGRpKdnc2VK1c4evQoaWlp7Ny5kzZt2uDh4QFAYGAg3bt354MPPgBg8+bNHD9+nN695YtkLdDpoKsdJsvZo01rDekIA+8Fzzu+lmvkDY/3go6hisQSQvNU9Rly27ZtdO3atdKyd955hzFjxpCeno6bmxtLly61OBNTqEd8BKzbY7psmC34uEMHBQ/89TroGw2Jd8P0TNOyp++H0Mam4i6EsA/VFLuioiLy8/NJSUmptDwsLIz169crlMrkg9V/Y9fhTUS36k5I4ygys+cyeei7xLZO4JOc19m8awVN/VvyzMMZ3Cy5zrPv9iW4UTjPj/zYaXMDbNi5nIUrJrFkxm8UXy9SJLefJwyMhc9/tk17zjJMePswalgT5XKYU9P20S60J/My/8qpC0dxMbiR+mgmbi7uim7XNWWNbhXP1Ld7cejkTt6Zsp3gwHAKzx1iXuZf0aEjsGEIz434CIPewIxFAykqvkD6xI0Ozy/sz2mHMe/k7e1NaWkpTz31lNJRKjl0Mo8r1y7xZsp6Ll09y7UbVxiW8AyxrRM4X3SK7QezSZ+4kdC7YtiU9zkeDbxJHZWpdGyzuctt+DWLxn7NARTN3SsK7mlmfp3Jiy3PxIyPgJjmtsulZea2j4MntuPi4sabKevp12ks3/28WNHtw1xWg96FWWM+p2e7hyrW93b3Y/bYVbyZsp6ggFB+2Gua+Db78VU1vYXQANUUO2eVd2gjHSNNF5HtEHE/+tvO38r/bRuxYYl/PNeXPUe2KBGxWuZyA3y/ZzUdIvqi0ym/iej1MLYn3H1X3dvoHAYPdbRdJq0zt30ENgzGaDSNKxcVX8DXs5EiGcuZy6rT6aqc5O7j6Y+Xh2k6rovBFb1OgdlKwuGcYEBH3S5fPceqLe/w6Ya3KCq+QELscPy8TeNSV65dwNPdFwAv94YUXbugYNLKzOUG+Oanf/HcIx/x1Y+LzLTiOG4u8EQirNsNX+00XV/SGu6ukNTBdN6afCdmPXPbh69XINdvFvP46/dg0Bn4f0//4LRZzTlz8QQ/5X/DqD4zHJBSKE2KXT35eAbwWL9XiI8exNbdqzh98dYZr17uDTn9xxmwV69dwtvdT6GUVZnL/cuBdbRp2Q1XF+ea6GPQw/1tIbYFbNgHPxSY7lFXHe8GplmXPSJN3/uJ2jG3ffyU/zUNvRqz6Jk9rP81i2W58xl9v3L3sTOXtSY3Sq7z+n8eY+qw9zAYpBv8M1B+jErl2ob2YGeBaYLMjoM5FcM7AJHNO/FrQS4AP+//lntadq22DSWYy334ZB5bdq3kb+/158jvu/jwK+f65NvEF4Z2grShMLkfDL1tePLhLqaro7w8xDTFXwpd3ZjbPsrKyiruFNDQK5Ar1y4qkrGcuaw1Sc8az6D4ibRs2sbe8YSTkGJXT6FBbXExuDJtYSIuBlfc3bwqnvP3bkK7sF5MXtCDgye2Ex+dpFzQO5jLPaTH07z+5DrmPvEVLZtGM7b/bAWT1szNBVoFQs+oW8u6hSt30riWmNs+OkY+wNFTe5i2MJF/rX2RQd1SzLRkf+ayAqR9NJyf9n/NvMzH2Jy3gt2Ht7AxbznLN6QzbWEiG3d+plBy4Uhy/G4D4wbMrfj3+l+zyMx+jeDACGJbJ/BI7+d4pPdzFc8XXy/itaWPEtW8kxJRKzGXu1z5NGxnyi0cw9z28eJfsyqtq/T2YS7rzNGfVFl/5eyqd8mdsWggAb71mAUlnJpTXgjamcgFc2tHqQvmlp92UNdbAClBycyyXdeOXAha/eTIzgIfBU/4rc97K5Vbyb+XsJ5s1+p4X2E7UuwsiLpP6QR1o9bcwjHUun2oNbdQnkxQEUIIoXlS7IQQQmieFDshhBCaJ8VOCCGE5kmxE0IIoXlS7IQQQmieFDshhBCaJ8VOCCGE5kmxE0IIoXlS7IQQQmieFDshhBCaJ8VOCCGE5kmxE0IIoXly1wML9q2Dy6eUeW+fJnW/yrtSueuTWQhL1Lhdq7UP0RopdhZcPqXcTS7rQ625hTBHjdu1GjNrkQxjCiGE0DwpdkIIITRPhjGFcKBrN2FvIfx29tayjA3QzA9aBEJEUzDIR1AhbE6KnQ1MW5jIniNbMBhc0esNBPmHMrJPKgmxw5SOZpZac6vRhavwdR5sOwQ3Sio/t/2o6QHQ0APiI6D3PeAme2edqHG7VmNmtZHdyUZG9Z3JqL4zKC0tYcXmfzB3yUjCg9sTHBiudDSz1JpbTX4ogOXbTEd1llwshjW/moriyG4Q2tj++bRIjdu1GjOriQyY2JjB4MKDXZ6g1FjCwRPblY5jNbXmdnard8CSLdYVutudvgz/+BZ2/mafXH8Watyu1ZhZDaTY2djNkhus2rwQgJDASIXTWE+tuZ3Z+n2mocu6KjVCxkYoUOgcLS1Q43atxsxqIMXORpZ89ypJM/0Y+IIHH66dwdRh7xPWLAaA42cOkJIex82SGwB8kvM6GWtfVDJuBXO55yweydbdqyrWfSkjiW37vlYqqqr8fhFW/mx+nfRRpoc5pUZYsrXq93zCPDXuj7Iv2pdTFzuj0cj8+fOJiIjA3d2d2NhYcnNziYqKYvz48UrHq2Rkn1Q+T7tA1stn6Hz3AHYcyK54LjgwnB7thpK5bi6F5w6Rsz2TkX1SFUx7i7ncyYPTyVg7k+LrRWzYuRwv94Z0jHpAwbTqsfwnKDHapq0zl2HdHtu0VR83S03fJX7xi2l4tuAUlJUpnap6atwfZV+0L6cuduPGjSMtLY0JEyawZs0ahg8fzogRIygoKCAuLk7peNXy8fRn6rD3+X7vl2zOW1GxfHjiM2zds4o5i0eQPCgdN5cGCqasqrrc/t5NGNJjEgtWPM2S72bz5KC3FE6pDr9fgn2Ftm1z837TUZ5Sth+BFz+FjzfDd7tNw7P/+w3MX2Mqxs5Kjfuj7Iv24bTFbunSpWRkZLBy5UqmT59O7969SU1NpVu3bpSUlNChQwelI9bI1zOAoT2nsuirFzAaTT2Ui8GVdmG9KCo+T9vQHgonrF51uft1GsOx0/kkdX8aX88AhROqw0+HbN/mpWLIP2n7dq2x65jpu8PiaibZnLgA/+8buHzN4bGspsb9UfZF23PaYjdnzhz69+9PQkJCpeXh4eG4uroSE2May37xxReJjIxEr9eTlZWlRNRqDek5iXOXCvnmp38DcPjkLnYd3kT78L6s/v49hdPV7M7cAM0ahcv051o4ctbyOs7UrjllZbDiF9CZef5iMWzMd2isWlPj/ij7om055Xl2x44dIy8vjylTplR57ujRo0RHR9OggWnYoX///owZM4bHH3/c0TErvJGcU2WZl7svy185B5i+e/z78id5asgCQgIjmbQgnvjowfj7NHVw0sos5RZ1U3jBTu2et0+75hw5C6cuWV5v8354MMb+eayhxv1R9kX7c8oju2PHTJcIDwoKqrS8uLiY3NzcSkOY8fHxhIWF1fo9dDqdVY/c3Jx6/S4AX2xZSERwHJEhcXi6+zCmXxpvr5xs8XW5uTlW57RH7rqoT+b6PMop8d53Pk6fq1wdymdd3vmw9vlyK778yuG/y/0DHrLq//vla6DT6+2eR6n9UY37Yn1zq+VhLac8sgsMDAQgPz+fAQMGVCyfN28ehYWFTjs5pSaDu0+s9HP3tkl0b5ukTJg6ePaRDKUjqIrxjyntamnXnBvFF61a7+b1K847NfMOat4fZV+sO6csdmFhYcTExDBnzhwCAgIIDg4mKyuL1atXA9ik2JVZuWNuy1TuXlQJCYmULaxbB6JU7vpkro/Ji03/tfb/qz2lr4XDZ279XJ6tynqjzD9/p3GjBrHyDcf+fiWl8NJncOV6zevogJ5tvFjggL+9GrdrtfYhWuOUw5h6vZ5ly5YRHR1NcnIyY8eOJTAwkIkTJ2IwGCompwjhjJrbaaKcvdo1x8UAve+p+XkdoNdDwt0OiyREnTjlkR1AZGQk2dnZlZaNHj2aNm3a4OHhoVAqISyLaQEbbDw70dUAdzezbZvWuq+N6a4NG/NNxe324wSDHsb2hLv8lMkmhLWc8siuJtu2basyhDlz5kxCQkLYsmULEyZMICQkhIMHDyqUUAgIbwJNfW3bZlwr8HSzbZvW0uvgoU4wtT90um0u2IAYeDEJokOUySVEbTjtkd2dioqKyM/PJyUlpdLytLQ00tLSFEpl8sHqv7Hr8CaiW3UnpHEUmdlzmTz0XaJbxTP17V4cOrmTd6ZsJzgwnOLrRTz7bl+CG4Xz/MiPnTJ3bOsEBs9sSHiz9gC89NhyfD0DmLFoIEXFF0ifuFHR3M5Op4PBHeDdHNu018AF+rWzTVv10aKR6bZDPxSYfn7ACTJVx9x2/Y/Pn+JQ4U7uahTGlIfew6A3KLpd16bvuHjlDC9+OAiDwRUv94bMePQ/GI2lTtOfODvVHNl5e3tTWlrKU089pXSUSg6dzOPKtUu8mbKeS1fPcu3GFYYlPENs6wQMehdmjfmcnu1uTd/2aOBN6qhMBRObmMsNEBrUjjeSc3gjOafiag2zH19lrklxmzbB0Ln2Z8RUKykO/L1s05bWmduu9/32IyUlN3gjOYeWTaP5/o8LKyu1Xde27/D28OetlI28mZxLZHAcW3evcpr+RA1UU+ycVd6hjXSMNF2QtUPE/ej1horndDqd4ieO18RcboCjp/Yw5e2evL/6eaeY4ahGwzpDZFDNz09ebHkm5n1toGtr2+bSMnPbdeHZAkLvMk1ua93sXnYd2axIxnK17TsMegN6vanLLi0rJTgwwnFhNUCKXT1dvnqOf339EtMWJrLku1e5fFUdVzywlDvjuf28mbyeoqvn2bL7C4VSqpurAZ5IrNsRnkFvGgr973tNw6LCOua265DGUfxakAvA9gPruFJ8QaGUJnXpO/Ye/YGUv3dk+4F13BUQ6oCU2qGa7+yclY9nAI/1e4X46EFs3b2K0xcVOqGmlizlLh+6jG+bxIHjvxAfPUiJmKrnajB9zxXbAlb8bN2ltyKawl86ygzHujC3XYcH30uroLZMf6c3rYLa4qfwqEtd+o67W3Tm7UnbWJb7Bl/9uIihvapeUlFUT47s6qltaA92FqwHYMfBHIzGUoUTWcdc7uIbVyj94+ddhzfRrJGMo9VXdDD8bSBM7AM9IqBlI/BwMxVD7wam4c4+beD5gTCxrxS6urK0P46+/0XmP5mNr2cjutzzX0pErFDbvuPmbVfQ8XL3xc1VTsGqDTmyq6fQoLa4GFyZtjCRNi274e7mRanx1m2l0z4aTt7hjRw/s5+HE58lvu1gBdPeYi738dP7eWPZ43i4eRMUEMpfH5ilcFpt0OkgIsj0EPZhbrs2Go0888/70OsNtA/vwz0tujhtVqjadwT43sW7Xz6DXqfHxyOA50Z8pGB69ZFiZwPjBsyt+Pf6X7PIzH6N4MAIYlsnMHP0J5XWLb5exGtLHyWqeSdHx6zCXO6Fk3+usv6MRQMJ8L3LkRGFqDVz23V1dxdQcruuTd8B8GZybqWfnak/cXa6MplqZ5aS17XzC4GOj9TttUrlrk/m+iif1Vjd3QKEbSn5t1bjdq3WPkRr5MjOAp8m6nxvpXIr+fcS2qfG7VqtfYjWSLGzIOo+pRPUjVpzC2GOGrdrNWbWIpmNKYQQQvOk2AkhhNA8KXZCCCE0T4qdEEIIzZNiJ4QQQvOk2AkhhNA8KXZCCCE0T4qdEEIIzZNiJ4QQQvOk2AkhhNA8KXZCCCE0T4qdEEIIzZNiJ4QQQvPkrgcW7FsHl08p894+Tep+xXSlctcnsxBapMY+RI2ZLZFiZ8HlU8rdeLE+1JpbCK1R476oxsyWyDCmEEIIzZNiJ4QQQvOk2AkhrFJqhBPnb/18tgjKypTLI0RtyHd2Qoga3SyFX47A1gNw9CyUGG89l7YCvBpAVBD0iITQxqDTKZdVCHOk2NnAtIWJ7DmyBYPBFb3eQJB/KCP7pJIQO0zpaGapNbdwjN3H4T/fw8Ximte5ch1+PmJ6RAbBI10gwNtxGbVCrfuimnLLMKaNjOo7ky9eLWL5y2fp12kMc5eM5PiZA0rHskituYX9GMtg+TZ4N8d8obtT/kn4v1/CLo3N4nMUte6Lasktxc7GDAYXHuzyBKXGEg6e2K50HKupNbewrbIyyPoB1u+r2+uvl8AH62HXcdvm+jNR677o7Lml2NnYzZIbrNq8EICQwEiF01hPrbmFbf14CDZb+FCePsr0qImxDD7aBBeu2jbbn4Va90Vnz+3Uxc5oNDJ//nwiIiJwd3cnNjaW3NxcoqKiGD9+vNLxKlny3askzfRj4AsefLh2BlOHvU9YsxgAjp85QEp6HDdLbgDwSc7rZKx9Ucm4FczlnrN4JFt3r6pY96WMJLbt+1qpqMLOLhXDZ9ts09a1m/DJ97Zpyxau3YTLxaYZpc5K+hD7cupiN27cONLS0pgwYQJr1qxh+PDhjBgxgoKCAuLi4pSOV8nIPql8nnaBrJfP0PnuAew4kF3xXHBgOD3aDSVz3VwKzx0iZ3smI/ukKpj2FnO5kwenk7F2JsXXi9iwczle7g3pGPWAgmmFPW3aD8U3bdfe7hNw/Lzl9expzwlY8C08/wnMXA6pWfDZT7X7LtJRpA+xL6ctdkuXLiUjI4OVK1cyffp0evfuTWpqKt26daOkpIQOHTooHbFaPp7+TB32Pt/v/ZLNeSsqlg9PfIate1YxZ/EIkgel4+bSQMGUVVWX29+7CUN6TGLBiqdZ8t1snhz0lsIphb2UGmHLftu3u8kObVordy/8MxsO3HaNx2s3TcvfXGM6T9AZSR9iH05b7ObMmUP//v1JSEiotDw8PBxXV1diYmI4f/48AwcOJDIyktjYWB544AEOHFB+FpCvZwBDe05l0VcvYDSaxk1cDK60C+tFUfF52ob2UDhh9arL3a/TGI6dziep+9P4egYonFDYy+8X4dI127e7/6Tt27TGsXOmIzio/sT3S8Xw8WbHZqoN6UNszymL3bFjx8jLy2PYsKrnahw9epTo6GgaNGiATqdj8uTJ5Ofns2PHDgYOHMjYsWMVSFzVkJ6TOHepkG9++jcAh0/uYtfhTbQP78vq799TOF3N7swN0KxROMGB4QqmEvb22zn7tHv6MhTfsE/b5mzMB3Pnt5cBh04rP8xqjvQhtuWUJ5UfO2Y6UScoKKjS8uLiYnJzc3nwwQcB8PPzo2/fvhXPx8fHM2/ePKveQ2flpR7mP5lNbOtEs+u8kZxTZZmXuy/LXzH1IEajkb8vf5KnhiwgJDCSSQviiY8ejL9PU7Pt5ubm0GlEb6ty2iN3XdQnc31M+tj08d3a/6+iss5JM+j2UFqlZeZmXJp7fvLiyj8Hh97D+RN765Gu9sa8VUDDxqEW1xs4YhLb1/6vXbOosQ+xJrM1ueuitpnLrLxmnVMe2QUGBgKQn59fafm8efMoLCyscXJKeno6SUlJ9o5Xa19sWUhEcByRIXF4uvswpl8ab6+crHQsIW5jvw8JOju2XRO93mDVejor11Oa9CH1pyuztiw6kNFopH379hQWFjJ//nyCg4PJyspi9erVHD16lK1bt9KlS5dKr5k1axZr1qxh3bp1eHp62izLtkzl7uvkFwIdH6nba5XKXZ/M9VF+NGHpaERUb1M+LPvRunXL/8Z3HsHVZPZQ8HavW666+nA9/HrM8oWq/09fCDd/cFRvauxD1JjZEqc8stPr9Sxbtozo6GiSk5MZO3YsgYGBTJw4EYPBQExMTKX1Z8+ezapVq/jqq69sWuiE+LMIsdO8AT9Pxxc6gO6R5gudTgdNfKF1E8dlEspyyu/sACIjI8nOzq60bPTo0bRp0wYPD4+KZbNmzWL16tV88803+Pn5OTilENrQzB/cXU1T821JqWIS0dR0J4aN+VWf0+nARQ+PxstdGv5MnLbYVWfbtm107dq14uddu3bx8ssv07p1axITEyuWb9++3fHhhFAxVwN0CoMNdbwmZk26KTQBT6eDoR2hsQ9k76l86bK774KB90KwvzLZhDJUU+yKiorIz88nJSWlYll0dLTVM3GEEOb1ioTN+213Sa3mAcoOE+p0kHA39IyEqUtNy15KAn8v5TIJ5aim2Hl7e1NaWqp0jGp9sPpv7Dq8iehW3QlpHEVm9lwmD32XJv4tmJf5V3ToCGwYwnMjPsKgNzBj0UCKii+QPnGjU+aObZ3AT/nfkLluLsYyIxP++w2aN47i2Xf7EtwonOdHfqxobmEfjX2hfzv4ckf92zLoYURX5xgm1N82M8HZCl1N+2B0q3imvt2LQyd38s6U7RXnqH2S8zqbd62gqX9Lnnk4g5sl1x2+X9Y2M8CGnctZuGISS2b8RvH1IkX6EqecoKImh07mceXaJd5MWc+lq2e5duMKwxKeIbZ1At7ufsweu4o3U9YTFBDKD3tXAzD78VUWWrU/c7mv3yzmy63/5LXx3/BGcg6RIXF4NPAmdVSm0rGFnd3XxjTMZ87kxZZnYg6JM30PKGpmbh806F2YNeZzerZ7qGL980Wn2H4wm/SJGwm9K4ZNeZ87fL+sbeZyG37NorFfcwDF+hIpdvWUd2gjHSNNFzbtEHF/pfN7fDz98fJoCJgu9aPXOc85PeZy7z6yBZ1OzwvvP8hrS0dTfOOKUjGFgxn08HgvaNOsbq/XYSp0PZzvDi9Ox9w+qNPpqpwwnv/bNmLDEv9Yvy97jmxxWNZytc0M8P2e1XSI6ItOp2y5kWJXT5evnuNfX7/EtIWJLPnuVS5frXrlgDMXT/BT/jcVG4kzMJf7/OXfOXe5kDn/s4bolvF8ueWfCiYVjubmAv+TAEkdwKUWn88CfeCp+03fkwnLrOk7bnfl2gU83X0B8HJvSNG1Cw5IWVltMwN889O/6NPhUQekM08139k5Kx/PAB7r9wrx0YPYunsVpy9WPhPzRsl1Xv/PY0wd9h4Gg/P8uc3l9nJvSNtWPTDoDdwbfh/LcucrmFQoQa+HxHsgtoVp0sqWg1BUw4WiQwKgRwR0aGUqlMI6lvqOO3m5N+T0H2d6X712CW93PwekrKy2mX85sI42Lbvh6uLmoIQ1kyO7emob2oOdBesB2HEwB6Ox8iSa9KzxDIqfSMumbZSIVyNzuaOad+LoqT0AHDyxnaAAy9cYFNrk7wX/dS+88hdI/W94rAcM6wTDO8OE3qaro0x/ELqGS6GrLUt9x50im3fi14JcAH7e/y33tOxqdn17qG3mwyfz2LJrJX97rz9Hft/Fh1/NcETMasnmWU+hQW1xMbgybWEibVp2w93Ni1JjCQC7D29hY95yfj9/hOUb0hnSYxI92g1ROLGJudx+3o2JCUtg6tu9aODqyd9GLVE4rVCaXmeardnYV+kk2mFuHwRI+2g4eYc3cvzMfh5OfJb4toNpF9aLyQt60MSvBX/pOdnpMw/p8TRDejwNwOQFPRjbf7bDM5eTYmcD4wbMrfj3+l+zyMx+jeDACGJbJ7By9uUq689YNJAAXwtT3hzAXO6hvaYwtNeUiueLrxfx2tJHiWreSYmoQmiSuX1w5uhPqqz/SO/neKT3cxU/K7Ff1jZzufJTrZTqS5zyQtDORK0XRJULQQtxi5Lbhxr7EDVmtkSO7CzwUfAKEPV5b6VyK/n3EsIZqbEPUWNmS6TYWRB1n9IJ6katuYXQGjXui2rMbInMxhRCCKF5UuyEEEJonhQ7IYQQmifFTgghhOZJsRNCCKF5UuyEEEJonhQ7IYQQmifFTgghhOZJsRNCCKF5UuyEEEJonhQ7IYQQmifFTgghhOZJsRNCCKF5ctcDC/atg8unlHlvnyZ1v/q4Urnrk1kILVJjH6LGzJZIsbPg8inlbmJYH2rNLYTWqHFfVGNmS2QYUwghhOZJsRNCCKF5MowphNCkU5dg13H47eytZf/4Fpr5QatAaBsCbtID/mnI/2obmLYwkT1HtmAwuKLXGwjyD2Vkn1QSYocpHc0steYWwpzDZ2DNDth3supzB343PdbvAw836BYOD7QFd1fH57ydWvdFNeWWYmcjo/rOZFTfGZSWlrBi8z+Yu2Qk4cHtCQ4MVzqaWWrNLcSdSo3w5Q7I3g1lVqxffAPW7YZfjsCobhDe1O4RzVLrvqiW3PKdnY0ZDC482OUJSo0lHDyxXek4VlNrbiHAVOg+3GAqXtYUutudvwILv4Odv9klWq2pdV909txS7GzsZskNVm1eCEBIYKTCaayn1txCAGT9CHn1mCpfWgYZG+HoWcvr2pta90Vnzy3FzkaWfPcqSTP9GPiCBx+uncHUYe8T1iwGgONnDpCSHsfNkhsAfJLzOhlrX1QybgVzuecsHsnW3asq1n0pI4lt+75WKqoQ1dp9HLYcML9O+ijTw5xSIyzZAjdLbZetNqQPsS+nLnZGo5H58+cTERGBu7s7sbGx5ObmEhUVxfjx45WOV8nIPql8nnaBrJfP0PnuAew4kF3xXHBgOD3aDSVz3VwKzx0iZ3smI/ukKpj2FnO5kwenk7F2JsXXi9iwczle7g3pGPWAgmmrKr4BG/bd+vmHArhRolwe4VjGMli+zXbtnbwIm/Jt115tSB9iX05d7MaNG0daWhoTJkxgzZo1DB8+nBEjRlBQUEBcXJzS8arl4+nP1GHv8/3eL9mct6Ji+fDEZ9i6ZxVzFo8geVA6bi4NFExZVXW5/b2bMKTHJBaseJol383myUFvKZyysp8Pw4vL4dPbOrslW+Dlz2D/74rFEg60rxDOFNm2zU37TUVUKdKH2IfTFrulS5eSkZHBypUrmT59Or179yY1NZVu3bpRUlJChw4dlI5YI1/PAIb2nMqir17AaDQC4GJwpV1YL4qKz9M2tIfCCatXXe5+ncZw7HQ+Sd2fxtczQOGEt+wthI82VT/kVHwT/pkNJ847PpdwrJ8O277N05crn5unBOlDbM9pi92cOXPo378/CQkJlZaHh4fj6upKTIxpTDgpKYmYmBjat29P586d+fbbb5WIW8WQnpM4d6mQb376NwCHT+5i1+FNtA/vy+rv31M4Xc3uzA3QrFG4000jXrOj5ufKyqC0FL7b7bg8Qhn2mlDy2zn7tFsb0ofYllOeZ3fs2DHy8vKYMmVKleeOHj1KdHQ0DRqYDuEzMjLw8/MD4JdffiExMZFz585hMBgclveN5Jwqy7zcfVn+immPMRqN/H35kzw1ZAEhgZFMWhBPfPRg/H2UPbHHUm5ndfoyHLHQyZUB24/Aw13kKhlaVWqE05fs03bhBfu0WxPpQ+zPKY/sjh0zzSEOCgqqtLy4uJjc3NxKQ5jlhQ7g4sWL6HQ6ysosD7jrdDqrHrm5OfX+fb7YspCI4DgiQ+LwdPdhTL803l452eLrcnNzrM5pj9x1UZ/M1j5i47pblaW0DAKaBNs9jzyUebh7elU5p6581uWdD2ufL/feBx/aLKca+xCl+o+6ZLaWU37mDQwMBCA/P58BAwZULJ83bx6FhYVVJqdMnDiRNWvWcPHiRT799FNcXJzr1xrcfWKln7u3TaJ72yRlwtTBs49kKB2hkuKiM1atZzSWcP3qBfuGEYox/jENv6ysrFadnjVKblyzaXv1JX1I/enKrDkMcjCj0Uj79u0pLCxk/vz5BAcHk5WVxerVqzl69Chbt26lS5cuVV6Xm5vLlClTWL9+Pd7e3jbJsi1Tufs6+YVAx0fq9lqlctcns7XKymD+GtMElJo2Xh0Q0xzG9rJvFqGsV1eahrUtKT96m7zYunb/0hF6RdU91+3U2IeoMbMlTjmMqdfrWbZsGdHR0SQnJzN27FgCAwOZOHEiBoOhYnLKnRISEtDr9WzatMnBiYUj6XTQv535QqfTQZ9oR6YSSmhup8l99mpXKMe5xvtuExkZSXZ2dqVlo0ePpk2bNnh4eABQVFTE2bNnadmyJWCaoHLw4EHuueceh+cVjtWuuWnyybIfTOdElY9ilZWBqwEe6wktGimbUdhf+5bw8xHbtunvBS1l29Ecpy121dm2bRtdu3at+PnKlSs8/PDDFBUV4eLigru7Ox9//DEtWrRQMKVwlG7hEB0MWw+azovS6aB1E+gUBp5uSqcTjtAmGPw94fxV27XZPQL0TjnmJepDNcWuqKiI/Px8UlJSKpY1bdqUrVu3KpjK5IPVf2PX4U1Et+pOSOMoMrPnMnnou7QKiubFDwdhMLji5d6QGY/+B6OxlGff7Utwo3CeH/mxU+b2cm/Iwj9mep06f4QhPSfxl56TmbFoIEXFF0ifuFHR3Lfz9TDdj0z8ORn0kBRnuuOBLQR6Q08bfVdnjZr2wehW8Ux9uxeHTu7knSnbCQ4Mp6T0ZpVlxdeLHN6f1CYzwOCZDQlv1h6Alx5bjq9ngCJ9iWo+v3h7e1NaWspTTz2ldJRKDp3M48q1S7yZsp5LV89y7cYVhiU8Q2zrBLw9/HkrZSNvJucSGRzH1t2r8GjgTeqoTKVjm80dHnwvbyTn8EZyDqF3xdDlnoEAzH58lYVWhXC82BYQ18r8OpMXW56cotfBiG7QwEGHAOb2QYPehVljPqdnu4cq1q9umaP7k9pmBggNalfRn5RfQUWJvkQ1xc5Z5R3aSMdI04VNO0Tcj15/62R2g96A/o/xkNKyUoIDIxTJWB1zucsV37jC+csnFb/ygRCWPNIVIoMsr1cTHTCym2kY3FHM7YM6na7KCePVLXO02mYGOHpqD1Pe7sn7q5+36hxoe1HNMKazunz1HKu2vMOnG96iqPgCCbHD8fO+tcfsPfoD//tZCm4u7gzrNU3BpJVZyg3w4941dIzqr1BCIaznaoAnEuGzbbDZwu1+7uTjDiO6mr7/cyRr9kFnU5fMGc/tx8fDn79/+iRbdn9BfPQgB6WtTIpdPfl4BvBYv1eIjx7E1t2rOH2x8skpd7fozNuTtrEs9w2++nERQ3tVvQSaEizlBtiU9xnDE59VIJ0QtedqgOFdTMOaq3dYvqScqwE6hcJ/3QteCtxAwJp90NnUJXP50GV82yQOHP9FsWInw5j11Da0BzsL1gOw42AORuOty/CX32gRTNeLc3P1cHi+mpjLDVBSepOjp/bQulmsEvGEqLOou2BKf5jWH/q1gzbNoIkvNPKGEH/oHAbDO8OsIabiqEShA8v7oDOqbebiG1co/WOdXYc30axRa7tnrIkc2dVTaFBbXAyuTFuYSJuW3XB386LUaLp76MET23n3y2fQ6/T4eATw3IiPFE57i7ncAL8cWMe9re9TMKEQ9dO8kenhrCztg2kfDSfv8EaOn9nPw4nPEt92cLXLnDlzE/+WvLHscTzcvAkKCOWvD8xyaN7bSbGzgXED5lb8e/2vWWRmv0ZwYASxrRN4Mzm30rrF14t4bemjRDXv5OiYVZjL3SmqH52i+lVaf8aigQT43uXomEJolrl9cOboT6qsf+cyJfqT2mZeOPnnKsuU6Euc8tqYzkSt14jT8rUxhVATNfYhasxsiRzZWeCj4OSo+ry3UrmV/HsJ4YzU2IeoMbMlcmQnhBBC82Q2phBCCM2TYieEEELzpNgJIYTQPCl2QgghNE+KnRBCCM2TYieEEELzpNgJIYTQPCl2QgghNE+KnRBCCM2TYieEEELzpNgJIYTQPCl2QgghNE+KnRBCCM2TYieEEELzpNgJIYTQPCl2QgghNE+KnRBCCM2TYieEEELz/j+SiJIgyOB6dAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 568.197x264.88 with 1 Axes>"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "'''\n",
    "VARIATIONAL ANSATZ\n",
    "Note that for molecules other than H_2 you may need to specify a different number of reps.\n",
    "'''\n",
    "\n",
    "reps = 1   \n",
    "ansatz = EfficientSU2(system_size, su2_gates=['rx', 'ry'], entanglement='circular', reps=reps, skip_final_rotation_layer=False)  \n",
    "    \n",
    "ansatz.decompose().draw('mpl')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "d5f6e337",
   "metadata": {},
   "outputs": [],
   "source": [
    "'''\n",
    "COST FUNCTION\n",
    "'''\n",
    "\n",
    "def basis_change_circuit(pauli_op):\n",
    "\n",
    "    #   Generating circuit with just the basis change operators\n",
    "    #\n",
    "    #       pauli_op: n-qubit Pauli operator\n",
    "\n",
    "    basis_change = QuantumCircuit(ansatz.num_qubits, ansatz.num_qubits)\n",
    "    \n",
    "    for idx, op in enumerate(pauli_op):\n",
    "        if op == 'X':\n",
    "            basis_change.h(idx)\n",
    "        elif op == 'Y':\n",
    "            basis_change.h(idx)\n",
    "            basis_change.p(-np.pi/2, idx)\n",
    "\n",
    "    return basis_change\n",
    "\n",
    "\n",
    "\n",
    "def min_cost():\n",
    "    \n",
    "    return min(cost_history)\n",
    "\n",
    "\n",
    "def log_cost(cost):\n",
    "\n",
    "    global cost_history\n",
    "    cost_history.append(cost)\n",
    "    \n",
    "    \n",
    "\n",
    "\n",
    "def objective_function(operators, params):\n",
    "    \n",
    "    backend = QasmSimulator(method='statevector', shots=1)\n",
    "    \n",
    "    pauli_op_dict = Counter(tuple(x) for x in operators)    \n",
    "    \n",
    "    shadow = []\n",
    "    for pauli_op in pauli_op_dict:\n",
    "        \n",
    "        qc = ansatz.bind_parameters(params)\n",
    "        qc = qc.compose(basis_change_circuit(pauli_op))                                   \n",
    "        qc.measure(reversed(range(system_size)), range(system_size))\n",
    "        result = execute(qc, backend, shots=pauli_op_dict[pauli_op]).result()\n",
    "        counts = result.get_counts()    \n",
    "\n",
    "        \n",
    "        for count in counts:\n",
    "            for _ in range(counts[count]): # number of repeated measurement values\n",
    "                \n",
    "                output_str = list(count)\n",
    "                output = [int(i) for i in output_str]\n",
    "                eigenvals = [x+1 if x == 0 else x-2 for x in output]\n",
    "                snapshot = [(op, eigenval) for op, eigenval in zip(pauli_op, eigenvals)]\n",
    "                \n",
    "                shadow.append(snapshot)\n",
    "    \n",
    "    \n",
    "    cost = 0.0\n",
    "    for term, weight in zip(observables_xyze, hamiltonian_coefficients):\n",
    "        \n",
    "        sum_product, match_count = estimate_exp(shadow, term)\n",
    "        if match_count != 0:\n",
    "            exp_val = sum_product / match_count\n",
    "            cost += weight * exp_val\n",
    "            \n",
    "    log_cost(cost)\n",
    "    \n",
    "    return cost\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "cb65b94c",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "NUMBER OF OPERATORS | DERANDOMISED OPERATORS | AVERAGE RMSE ERROR\n",
      "\n",
      "EXPERIMENT 0: GROUND STATE ENERGY FOUND = -0.8557807804833842\n",
      "EXPERIMENT 1: GROUND STATE ENERGY FOUND = -0.8557807804833842\n",
      "EXPERIMENT 2: GROUND STATE ENERGY FOUND = -0.8557807804833842\n",
      "50 | Counter({('X', 'X', 'X', 'X'): 50}) | 1.004219219516616\n",
      "EXPERIMENT 0: GROUND STATE ENERGY FOUND = -1.852959385123308\n",
      "EXPERIMENT 1: GROUND STATE ENERGY FOUND = -1.8593559426914374\n",
      "EXPERIMENT 2: GROUND STATE ENERGY FOUND = -1.863011118444654\n",
      "100 | Counter({('X', 'X', 'X', 'X'): 99, ('Z', 'Z', 'Z', 'Z'): 1}) | 0.004436661731014558\n",
      "EXPERIMENT 0: GROUND STATE ENERGY FOUND = -1.8514333028608607\n",
      "EXPERIMENT 1: GROUND STATE ENERGY FOUND = -1.8404121130262883\n",
      "EXPERIMENT 2: GROUND STATE ENERGY FOUND = -1.8505148703746463\n",
      "200 | Counter({('X', 'X', 'X', 'X'): 197, ('Z', 'Z', 'Z', 'Z'): 3}) | 0.013503595817144557\n",
      "EXPERIMENT 0: GROUND STATE ENERGY FOUND = -1.9819272519438194\n",
      "EXPERIMENT 1: GROUND STATE ENERGY FOUND = -1.9826926123489979\n",
      "EXPERIMENT 2: GROUND STATE ENERGY FOUND = -1.980855747376569\n",
      "600 | Counter({('X', 'X', 'X', 'X'): 590, ('Z', 'Z', 'Z', 'Z'): 9, ('X', 'Z', 'X', 'Z'): 1}) | 0.12182753323933261\n",
      "EXPERIMENT 0: GROUND STATE ENERGY FOUND = -1.9757002983546128\n",
      "EXPERIMENT 1: GROUND STATE ENERGY FOUND = -1.9787342056680677\n",
      "EXPERIMENT 2: GROUND STATE ENERGY FOUND = -1.9791938885943487\n",
      "1000 | Counter({('X', 'X', 'X', 'X'): 981, ('Z', 'Z', 'Z', 'Z'): 16, ('X', 'Z', 'X', 'Z'): 3}) | 0.11788632054259642\n",
      "EXPERIMENT 0: GROUND STATE ENERGY FOUND = -1.9763711727510254\n",
      "EXPERIMENT 1: GROUND STATE ENERGY FOUND = -1.974590579991244\n",
      "EXPERIMENT 2: GROUND STATE ENERGY FOUND = -1.9775773807495867\n",
      "1600 | Counter({('X', 'X', 'X', 'X'): 1570, ('Z', 'Z', 'Z', 'Z'): 25, ('X', 'Z', 'X', 'Z'): 5}) | 0.11618618870030964\n",
      "Execution time =  01:01:01\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x7f444f3b0220>"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEGCAYAAABo25JHAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAoOElEQVR4nO3de5xVdb3/8deb4S5zVC6aCQkIIiB4Q6TUxFS0i6h5AY551Cyrk6Z28lb+PGbHk5qnTLOUSlGPQmam5DHxBh6OJQJq3hFEjCFTNPGCchn4/P5Ya4Y9w1w246y9Z2a9n4/Hfux123t9Zs3s/Zn1+a71/SoiMDOz/OpU7gDMzKy8nAjMzHLOicDMLOecCMzMcs6JwMws5zqXO4At1bdv3xg4cGC5wzAza1cWLlz4ZkT0a2hdu0sEAwcOZMGCBeUOw8ysXZH0amPrXBoyM8s5JwIzs5xzIjAzy7l210Zgljfr16+nqqqKNWvWlDsUawe6d+9O//796dKlS9GvcSIwa+OqqqqorKxk4MCBSCp3ONaGRQRvvfUWVVVVDBo0qOjXuTRk1satWbOGPn36OAlYsyTRp0+fLT57zCwRSLpB0huSnm1kvSRdLWmJpKcl7ZVVLGbtnZOAFaslfytZnhFMAw5vYv1ngaHp4zTgFxnGAv/3f/C978GGDZnuxsysvcksEUTE/wL/aGKTI4GbI/EYsI2kHbKKh3nz4D//E1avzmwXZmbtUTnbCHYElhfMV6XLNiPpNEkLJC1YuXJly/ZWWZk8v/dey15vlmOXXnopI0eOZPTo0eyxxx7MmzeP8ePHM2zYMHbffXf2228/Fi1aBMA999zDnnvuye67786IESO4/vrr67zXUUcdxbhx4zbbx1VXXcXNN9+c6c/xqU99CoBly5Zx2223ZbKPc845h1133ZXRo0dz9NFHs2rVqtp1P/zhDxkyZAjDhg1j1qxZtcvvu+8+hg0bxpAhQ7jssstql7/yyivsu+++DBkyhEmTJrFu3ToAfvazn3HDDTe0XtARkdkDGAg828i6e4D9C+YfAsY095577713tMhtt0VAxAsvtOz1ZmXy/PPPl3X/f/rTn2LcuHGxZs2aiIhYuXJlrFixIg488MCYP39+RERcf/31ccQRR8S6detihx12iOXLl0dExJo1a+LFF1+sfa+33347+vfvH7vuumu8/PLLtcvXr18fo0aNivXr15fkZ5o9e3Z8/vOf36LXFBvbrFmzarc999xz49xzz42IiOeeey5Gjx4da9asiaVLl8bgwYOjuro6qqurY/DgwfHyyy/H2rVrY/To0fHcc89FRMRxxx0X06dPj4iIr33ta/Hzn/88IiJWr14de+yxR6MxNPQ3AyyIRr5Xy3lGsAIYUDDfP12WjZozgvffz2wXZpk76ywYP751H2ed1eQuX3vtNfr27Uu3bt0A6Nu3Lx//+MfrbPPpT3+aJUuW8N5771FdXU2fPn0A6NatG8OGDavd7s477+SII45g8uTJzJgxo3b5ww8/zF577UXnzskV7ePHj+fss89mzJgxDB8+nPnz5/PFL36RoUOHcuGFF9a+7qijjmLvvfdm5MiRTJ06FYBXX32VoUOH8uabb7Jx40YOOOAA7r//fgB69eoFwPnnn8/cuXPZY489+MlPfsKGDRs455xz2GeffRg9enTtWcycOXM44IADmDhxIiNGjGjyONWYMGFC7c8xbtw4qqqqALj77ruZPHky3bp1Y9CgQQwZMoTHH3+cxx9/nCFDhjB48GC6du3K5MmTufvuu4kIHn74YY499lgATjrpJO666y4AevbsycCBA3n88ceLiqk55byPYCZwuqQZwL7AOxHxWmZ7S/8AXBoy2zITJkzgkksuYZddduGQQw5h0qRJHHjggXW2+cMf/sCoUaPo3bs3EydOZKedduLggw/mC1/4AlOmTKFTp+R/zunTp3PRRRex/fbbc8wxx/Dd734XgEcffZS99967znt27dqVBQsW8NOf/pQjjzyShQsX0rt3b3beeWfOPvts+vTpww033EDv3r358MMP2WeffTjmmGPYaaedOO+88/jGN77B2LFjGTFiBBMmTKjz3pdddhlXXnkl99xzDwBTp05l6623Zv78+axdu5b99tuv9jVPPPEEzz77bO11+QcccADvNfA9cuWVV3LIIYfUWXbDDTcwadIkAFasWFGnJNa/f39WrEj+9x0wYECd5fPmzeOtt95im222qU0qhdsDjBkzhrlz5zJ27NjGf3lFyiwRSJoOjAf6SqoC/h3oAhAR1wH3Ap8DlgAfAKdkFQvgNgLrGK66quS77NWrFwsXLmTu3LnMnj2bSZMm1daxTzjhBHr06MHAgQO55pprAPjVr37FM888w4MPPsiVV17JAw88wLRp03j99ddZvHgx+++/P5Lo0qULzz77LLvtthuvvfYaw4cPr7PfiRMnAjBq1ChGjhzJDjsk15IMHjyY5cuX06dPH66++mp+//vfA7B8+XIWL15Mnz59+MpXvsJvf/tbrrvuOp566qlmf8b777+fp59+mjvuuAOAd955h8WLF9O1a1fGjh1b5+asuXPnFnXcLr30Ujp37swJJ5xQ1PZbarvttuPFF19slffKLBFExJRm1gfwzaz2vxmXhsxarKKigvHjxzN+/HhGjRrFTTfdBMCtt97KmDFjNtt+1KhRjBo1ihNPPJFBgwYxbdo0br/9dt5+++3aL9V3332X6dOnc+mll9KjR4/NboKqKUV16tSpdrpmvrq6mjlz5vDggw/y5z//mZ49ezJ+/Pja9/jggw9qSzLvv/8+lTWf/0ZEBNdccw2HHXZYneVz5sxhq622qrOsmDOCadOmcc899/DQQw/VXte/4447snz5putjqqqq2HHH5PqYhpb36dOHVatWUV1dTefOnetsD8mNhj169Gjy5ypWfu4sdmnIrEUWLVrE4sWLa+efeuopdtpppwa3ff/995kzZ06D206fPp377ruPZcuWsWzZMhYuXFjbTjB8+HCWLFmyRXG98847bLvttvTs2ZMXX3yRxx57rHbdeeedxwknnMAll1zCV7/61c1eW1lZWefL/LDDDuMXv/gF69evB+Cll15idSOXms+dO5ennnpqs0dNErjvvvu44oormDlzJj179qx93cSJE5kxYwZr167llVdeYfHixYwdO5Z99tmHxYsX88orr7Bu3TpmzJjBxIkTkcRBBx1Ue5Zy0003ceSRR9a+30svvcRuu+22RcesMfnpa8ilIbMWef/99znjjDNYtWoVnTt3ZsiQIUydOrW2EbNQRHDFFVfwta99jR49erDVVlsxbdo0li1bxquvvlqnRj5o0CC23npr5s2bx2c/+1lOPPHELYrr8MMP57rrrmP48OEMGzas9r0feeQR5s+fz6OPPkpFRQW/+93vuPHGGznllE3V59GjR1NRUcHuu+/OySefzJlnnsmyZcvYa6+9iAj69etX2zC7pU4//XTWrl3LoYceCiQNxtdddx0jR47k+OOPZ8SIEXTu3Jlrr72WiooKILkc9LDDDmPDhg18+ctfZuTIkQBcfvnlTJ48mQsvvJA999yTU089tXY/jz76KBdffHGLYqxPSYWm/RgzZky0aISyjRuhogIuugi+//3WD8wsIy+88MJm9fOO6Oijj+aKK65g6NCh5Q6lzXvyySf58Y9/zC233NLg+ob+ZiQtjIjN63jkqTTUqVNSHnIbgVmbdNlll/Haa9ldONiRvPnmm/zgBz9otffLT2kIkkTg0pC1QxHR4TueGzZsWJ17DqxxNWWnhrSkypOfMwJI2gmcCKyd6d69O2+99VaLPuCWL5GOR9C9e/ctel2+zggqK10asnanf//+VFVV0eJ+tixXakYo2xL5SgQuDVk71KVLly0abcpsS7k0ZGaWc04EZmY5l79E4DYCM7M68pUI3EZgZraZfCWCmjOCjRvLHYmZWZuRv0QA8MEH5Y3DzKwNyVcicA+kZmabyVcicA+kZmabcSIwM8u5fCYCX0JqZlYrX4nAbQRmZpvJVyJwacjMbDP5TAQuDZmZ1cpXInBpyMxsM04EZmY5l69EUFEBPXs6EZiZFchXIgD3QGpmVk/+EoF7IDUzqyN/icCD05iZ1ZHPRODSkJlZrfwlApeGzMzqyF8icGnIzKyOTBOBpMMlLZK0RNL5Daz/hKTZkp6U9LSkz2UZD+BEYGZWT2aJQFIFcC3wWWAEMEXSiHqbXQjcHhF7ApOBn2cVTy23EZiZ1ZHlGcFYYElELI2IdcAM4Mh62wTwT+n01sDfMown0atXkggiMt+VmVl7kGUi2BFYXjBflS4rdDHwJUlVwL3AGRnGk6isTAav97jFZmZA+RuLpwDTIqI/8DngFkmbxSTpNEkLJC1YuXLlR9ujeyA1M6sjy0SwAhhQMN8/XVboVOB2gIj4M9Ad6Fv/jSJiakSMiYgx/fr1+2hReUwCM7M6skwE84GhkgZJ6krSGDyz3jZ/BQ4GkDScJBF8xH/5m+EeSM3M6sgsEURENXA6MAt4geTqoOckXSJpYrrZvwFflfQXYDpwckTGrbg+IzAzq6Nzlm8eEfeSNAIXLruoYPp5YL8sY9iM2wjMzOood2Nx6bk0ZGZWR/4SgUtDZmZ15DcRuDRkZgbkMRG4NGRmVkf+EkHnztC9uxOBmVmqyUQgqZOk40sVTMm4B1Izs1pNJoKI2AicW6JYSsc9kJqZ1SqmNPSgpO9IGiCpd80j88iy5FHKzMxqFXND2aT0+ZsFywIY3PrhlIhLQ2ZmtZpNBBExqBSBlFRlJfzjH+WOwsysTWg2EUjqAnwD+HS6aA5wfUSszzCubFVWwquvljsKM7M2oZjS0C+ALmwaRvLEdNlXsgoqc24jMDOrVUwi2Ccidi+YfzjtLbT9chuBmVmtYq4a2iBp55oZSYOBDdmFVAI1l4963GIzs6LOCL4DzJa0FBCwE3BKplFlrVcv2LAB1qyBHj3KHY2ZWVk1mQgkVQC7A0OBYeniRRGxNuvAMlXYA6kTgZnlXHN3Fm8ApkTE2oh4On207yQA7oHUzKxAMaWhRyX9DPgNsLpmYUQ8kVlUWfOYBGZmtYpJBHukz5cULAvgM60eTam4K2ozs1rFtBHMjIiflCie0vAZgZlZraLaCEoUS+m4jcDMrFY+2whcGjIzq5XPNgKXhszMahXT++hBpQikpFwaMjOr1WwXE5K2l/RrSX9M50dIOjX70DLUpQt06+YzAjMziutraBowC/h4Ov8ScFZG8ZSOeyA1MwOKSwR9I+J2YCNARFTT3judA/dAamaWKiYRrJbUh6SBGEnjgHcyjaoUPIC9mRlQ3FVD3wZmAjtLehToBxybaVSl4NKQmRlQ3FVDT0g6kKT3UZH0Ptp+h6msUVkJ77T/Exszs4+qmDOCmnaB5zKOpbQqK2HFinJHYWZWdsW0EbSYpMMlLZK0RNL5jWxzvKTnJT0n6bYs46nDjcVmZkCRZwQtkXZYdy1wKFAFzJc0MyKeL9hmKHABsF9EvC1pu6zi2YzbCMzMgCbOCCR9qWB6v3rrTi/ivccCSyJiaUSsA2YAR9bb5qvAtRHxNkBEvFFs4B9ZzRmBxy02s5xrqjT07YLpa+qt+3IR770jsLxgvipdVmgXYBdJj0p6TNLhDb2RpNMkLZC0YOXKlUXsugiVlVBdDevWtc77mZm1U00lAjUy3dB8S3UmGQ95PEl317+UtE39jSJiakSMiYgx/fr1a509uwdSMzOg6UQQjUw3NN+QFcCAgvn+6bJCVSQD36yPiFdIuq8YWsR7f3TugdTMDGi6sXhXSU+T/Pe/czpNOj+4iPeeDwyVNIgkAUwG/rneNneRnAncKKkvSaloafHhfwTugdTMDGg6EQz/KG8cEdVpo/IsoAK4ISKek3QJsCAiZqbrJkh6nqT/onMi4q2Pst+i+YzAzAxoIhFExKuF82l/Q58G/hoRC4t584i4F7i33rKLCqaDpFH625Sa2wjMzICmLx+9R9Ju6fQOwLMkVwvdIums0oSXIZ8RmJkBTTcWD4qIZ9PpU4AHIuIIYF+Ku3y0bXMbgZkZ0HQiKOxY7mDSEk9EvEc6NkG75tKQmRnQdGPxcklnkFziuRdwH4CkHkCXEsSWLZeGzMyAps8ITgVGAicDkyJiVbp8HHBjtmGVQLduydjFLg2ZWc41ddXQG8DXG1g+G5idZVAl4x5IzcwaTwSSZjb1woiY2PrhlJh7IDUza7KN4JMkncZNB+bRev0LtR0+IzAzazIRfIxkLIEpJF1D/A8wPSI6zkhlHsDezKzxxuKI2BAR90XESSQNxEuAOUWORdA+uDRkZtb0CGWSugGfJzkrGAhcDfw++7BKpLIS/v73ckdhZlZWTTUW3wzsRnIj2fcL7jLuOFwaMjNr8ozgS8Bq4EzgW1JtW7FI+ov7p4xjy54bi83MmryPoKmbzToGtxGYmTV5Z3GDJG0j6XtZBFNylZXJmMUet9jMcqypbqgHSJqadkf9FUlbSfovYDGwXelCzJB7IDUza/KM4Gbgb8A1JH0OLQA+DoyKiDNLEFv23AOpmVmTjcW9I+LidHqWpOOAEyKi/XdBXcM9kJqZNXsfwbZs6lriLWBrpZcPRcQ/Mo4tey4NmZk1mQi2BhZSt4+hJ9LnAAZnFVTJ+IzAzKzJy0cHljCO8nAbgZnZll8+2qH4jMDMzIkAcBuBmeVavhOBS0NmZsUlAkn7Szolne4naVC2YZVI9+5QUeFEYGa51mwikPTvwHnABemiLsB/ZxlUyUjugdTMcq+YM4KjgYkkPZESEX8DKrMMqqTcA6mZ5VwxiWBdRATJvQNI2irbkErMPZCaWc4Vkwhul3Q9sI2krwIPAr/MNqwS8hmBmeVck11MAETElZIOBd4FhgEXRcQDmUdWKm4jMLOcazYRAKRf/B3ny79Qr16wcmW5ozAzK5tirhp6T9K79R7LJf1eUpP9DUk6XNIiSUsknd/EdsdICkljWvJDfCQuDZlZzhVzRnAVUAXcRtIB3WRgZ5IO6G4Axjf0IkkVwLXAoenr50uaGRHP19uukmRc5Hkt+gk+KpeGzCznimksnhgR10fEexHxbkRMBQ6LiN8A2zbxurHAkohYGhHrgBnAkQ1s9wPgcmDNlgbfKnxGYGY5V0wi+EDS8ZI6pY/j2fSlHU28bkdgecF8VbqslqS9gAER8T9NBSDpNEkLJC1Y2dr1/F69YM0aqK5u3fc1M2snikkEJwAnAm8Ar6fTX5LUAzi9pTuW1An4MfBvzW0bEVMjYkxEjOnXr19Ld9kw90BqZjlXzOWjS4EjGln9f028dAUwoGC+f7qsRiWwGzAnHfTsY8BMSRMjYkFzcbWawh5It22q0mVm1jE1mwgkdQdOJRnAvnvN8oj4cjMvnQ8MTTuoW0HSyPzPBa9/B+hbsJ85wHdKmgTAPZCaWe4VUxq6heS/9cOAR0j+s2/2WzMiqklKR7OAF4DbI+I5SZdImtjykFuZS0NmlnPFXD46JCKOk3RkRNwk6TZgbjFvHhH3AvfWW3ZRI9uOL+Y9W50HpzGznCvmjGB9+rxK0m4kg9pvl11IJeYzAjPLuWLOCKZK2ha4EJgJ9AL+X6ZRlZLbCMws55pMBOklnu9GxNvA/wJNdinRLvmMwMxyrsnSUERsBM4tUSzl4TYCM8u5YtoIHpT0HUkDJPWueWQeWan06AGdOvmMwMxyq5g2gknp8zcLlgUdpUwkeZQyM8u1Yu4sHlSKQMrKPZCaWY4VMx5BT0kXSpqazg+V9IXsQysh90BqZjlWTBvBjcA64FPp/ArgPzKLqBxcGjKzHCsmEewcEVeQ3lgWER+QDFDTcfiMwMxyrJhEsC7tcjoAJO0MrM00qlJzG4GZ5VgxVw1dDNwHDJB0K7AfcHKGMZWeS0NmlmPFXDV0v6SFwDiSktCZEfFm5pGVkktDZpZjxYxH8AeSgetnRsTq7EMqAycCM8uxYtoIrgQOAJ6XdIekY9PBajqOykr48EPYsKHckZiZlVyziSAiHomIfyW5k/h64HiS8Ys7jpoeSN1gbGY5VMwZAelVQ8cAXwf2AW7KMqiScw+kZpZjxbQR3A6MJbly6GfAI2mvpB2HeyA1sxwr5vLRXwNTImIDgKT9JU2JiG8287r2w4PTmFmOFXP56CxJe0qaQtI+8ApwZ+aRlZJLQ2aWY40mAkm7AFPSx5vAbwBFxEEliq10nAjMLMeaOiN4EZgLfCEilgBIOrskUZWa2wjMLMeaumroi8BrwGxJv5R0MB2ts7kabiMwsxxrNBFExF0RMRnYFZgNnAVsJ+kXkiaUKL7ScGnIzHKsmBvKVkfEbRFxBNAfeBI4L/PISqlnz2TISpeGzCyHirqhrEZEvB0RUyPi4KwCKotOnWCrrXxGYGa5tEWJoENzx3NmllNOBDWcCMwsp5wIaniUMjPLKSeCGh6lzMxyKtNEIOlwSYskLZF0fgPrvy3peUlPS3pI0k5ZxtMkl4bMLKcySwSSKoBrgc8CI4ApkkbU2+xJYExEjAbuAK7IKp5muTRkZjmV5RnBWGBJRCyNiHXADODIwg0iYnZEfJDOPkZyn0J5uDRkZjmVZSLYEVheMF+VLmvMqcAfG1oh6TRJCyQtWLlyZSuGWMClITPLqTbRWCzpS8AY4EcNrU9vYhsTEWP69euXTRCVlbB6NWzsWGPumJk1J8tEsAIYUDDfP11Wh6RDgO8BEyNibYbxNK2mv6HVq8sWgplZOWSZCOYDQyUNktQVmAzMLNxA0p7A9SRJ4I0MY2meeyA1s5zKLBFERDVwOjALeAG4PSKek3SJpInpZj8CegG/lfSUpJmNvF323AOpmeVUMWMWt1hE3AvcW2/ZRQXTh2S5/y3iwWnMLKfaRGNxm+DSkJnllBNBDZeGzCynnAhqOBGYWU45EdRwG4GZ5ZQTQQ23EZhZTjkR1HAiMLOcciKoUTNusUtDZpYzTgSF3AOpmeVQpjeUtTvugbRt+fDD5LmiIjljq6gAqbwxmXVATgSFnAjK6913YfZsuP/+5LFkScPbFSaGwuktfe7Ir+2UnuxLmx6F802ta26+XK9tK3E0t64dciIo5FHKSqu6GhYs2PTF/9hjsGFD0lZz0EFw8snJF9vGjcnyYp+3ZNvGXrN+fbb7iSj30bcsZZWAfvxjOOWUVg/XiaBQr17w+uvljqJje+WVTV/8Dz8Mq1Ylf+RjxsD558OECTBuHHTtWu5IsxWRbYKLqPuo2Wf96S2d/yivbStxtOfXDh1KFpwIClVWNl6OsJZ5552k3PPAA3XLPQMGwLHHJl/8n/kM9OlT3jhLTdpU6jErMyeCQi4NfXSNlXt69UrKPd/6VvLlv8su7bqmataROBEUcmNxyxSWex56KDkLyGO5x6ydciIo1KtXckYQ4f9Wm1JT7qn58n/55WT5gAFw3HH5LfeYtVNOBIUqK5MksHr1pi4nLCn3zJ+ffOk/8MDm5Z4zz3S5x6wdcyIoVNgDad4TwdKlmxp4Xe4x69CcCAoVdjz3sY+VN5ZSa6zc84lPuNxj1sE5ERTK0+A0heWe+++HefPqlnvOOgsOPdTlHrMccCIo1NEHp1m6dFOd3+UeM0s5ERTqaGcELveYWRGcCAq198FpXO4xsxZwIijUHs8IGiv37LOPyz1mVhQngkLtoY3A5R4za2VOBIW22ip5bktnBMWUeyZMSHoldLnHzFrAiaBQ587Qo0f5E0FNuaemq+bCcs8FFyR1fpd7zKyVOBHUV44eSN95J/nCr6n1u9xjZiXkRFBfKXogra6Gxx/f1IWDyz1mVkZOBPX16pVNImiu3DNhAuy7r8s9ZlZyTgT1teSMICL5Yn/9dfj735PnmukVK2Du3CQRQFLuOf74pM5/8MHQu3fr/wxmZlsg00Qg6XDgp0AF8KuIuKze+m7AzcDewFvApIhYlmVMzaqshDffTL7c33138y/2hqZffx3Wrt38vSoqYLvtki4czj7b5R4za5MySwSSKoBrgUOBKmC+pJkR8XzBZqcCb0fEEEmTgcuBSVnFVJRevZLyTY8eDX+5d+qUfLlvv33yGD48ef7YxzYtq5nu0yfZ3sysDcvyjGAssCQilgJImgEcCRQmgiOBi9PpO4CfSVJERIZxNe2005IkUPNlX/OlXvjl7gHHzawDyTIR7AgsL5ivAvZtbJuIqJb0DtAHeLNwI0mnAacBfOITn8gq3sQhhyQPM7OcaBd1i4iYGhFjImJMv379yh2OmVmHkmUiWAEMKJjvny5rcBtJnYGtSRqNzcysRLJMBPOBoZIGSeoKTAZm1ttmJnBSOn0s8HBZ2wfMzHIoszaCtOZ/OjCL5PLRGyLiOUmXAAsiYibwa+AWSUuAf5AkCzMzK6FM7yOIiHuBe+stu6hgeg1wXJYxmJlZ09pFY7GZmWXHicDMLOecCMzMck7t7SIdSSuBV+st7ku9m9DamLYcn2NrubYcn2NrmbYcG3y0+HaKiAZvxGp3iaAhkhZExJhyx9GYthyfY2u5thyfY2uZthwbZBefS0NmZjnnRGBmlnMdJRFMLXcAzWjL8Tm2lmvL8Tm2lmnLsUFG8XWINgIzM2u5jnJGYGZmLeREYGaWc+0+EUg6XNIiSUsknV+G/Q+QNFvS85Kek3Rmury3pAckLU6ft02XS9LVabxPS9qrBDFWSHpS0j3p/CBJ89IYfpP2Doukbun8knT9wBLEto2kOyS9KOkFSZ9sK8dO0tnp7/RZSdMldS/XsZN0g6Q3JD1bsGyLj5Okk9LtF0s6qaF9tWJ8P0p/r09L+r2kbQrWXZDGt0jSYQXLW/3z3FBsBev+TVJI6pvOl/TYNRabpDPSY/ecpCsKlmdz3CKi3T5IejV9GRgMdAX+AowocQw7AHul05XAS8AI4Arg/HT5+cDl6fTngD8CAsYB80oQ47eB24B70vnbgcnp9HXAN9LpfwWuS6cnA78pQWw3AV9Jp7sC27SFY0cyet4rQI+CY3ZyuY4d8GlgL+DZgmVbdJyA3sDS9HnbdHrbDOObAHROpy8viG9E+lntBgxKP8MVWX2eG4otXT6ApHfkV4G+5Th2jRy3g4AHgW7p/HZZH7dMP+RZP4BPArMK5i8ALihzTHcDhwKLgB3SZTsAi9Lp64EpBdvXbpdRPP2Bh4DPAPekf+BvFnxAa49h+qH4ZDrdOd1OGca2NcmXreotL/uxY9Mwqr3TY3EPcFg5jx0wsN4XxhYdJ2AKcH3B8jrbtXZ89dYdDdyaTtf5nNYcuyw/zw3FRjJO+u7AMjYlgpIfuwZ+r7cDhzSwXWbHrb2XhhoaF3nHMsVCWg7YE5gHbB8Rr6Wr/g5sn06XOuargHOBjel8H2BVRFQ3sP86Y0gDNWNIZ2UQsBK4MS1d/UrSVrSBYxcRK4Argb8Cr5Eci4W0nWMHW36cyvl5+TLJf9o0EUfJ4pN0JLAiIv5Sb1XZYwN2AQ5IS4yPSNon69jaeyJoMyT1An4HnBUR7xauiyRNl/w6XUlfAN6IiIWl3neROpOcFv8iIvYEVpOUOGqV8dhtCxxJkqw+DmwFHF7qOIpVruNUDEnfA6qBW8sdC4CknsB3gYua27ZMOpOciY4DzgFul6Qsd9jeE0Ex4yJnTlIXkiRwa0TcmS5+XdIO6fodgDfS5aWMeT9goqRlwAyS8tBPgW2UjBFdf/+lHkO6CqiKiHnp/B0kiaEtHLtDgFciYmVErAfuJDmebeXYwZYfp5J/XiSdDHwBOCFNVm0hvp1JEvxf0s9Gf+AJSR9rA7FB8rm4MxKPk5zN980ytvaeCIoZFzlTaab+NfBCRPy4YFXheMwnkbQd1Cz/l/TqhHHAOwWn960qIi6IiP4RMZDk2DwcEScAs0nGiG4otpKNIR0RfweWSxqWLjoYeJ42cOxISkLjJPVMf8c1sbWJY9fAPos5TrOACZK2Tc94JqTLMiHpcJKy5MSI+KBe3JOVXGk1CBgKPE6JPs8R8UxEbBcRA9PPRhXJBR9/p20cu7tIGoyRtAtJA/CbZHncWqOxo5wPklb+l0hazb9Xhv3vT3JK/jTwVPr4HEl9+CFgMckVAL3T7QVcm8b7DDCmRHGOZ9NVQ4PTP6AlwG/ZdHVC93R+Sbp+cAni2gNYkB6/u0iuyGgTxw74PvAi8CxwC8nVGmU5dsB0kraK9SRfXKe25DiR1OqXpI9TMo5vCUntuuZzcV3B9t9L41sEfLZgeat/nhuKrd76ZWxqLC7psWvkuHUF/jv9u3sC+EzWx81dTJiZ5Vx7Lw2ZmdlH5ERgZpZzTgRmZjnnRGBmlnNOBGZmOedEYK0m7cXxvwrmvyPp4lZ672mSjm1+y4+8n+OU9II6O+t9tUeSjpI0otxxWOtyIrDWtBb4Yk2Xvm1FwZ3AxTgV+GpEHJRVPK1lC3+u1nIUSS+Y1oE4EVhrqiYZU/Xs+ivq/0cv6f30eXzasdbdkpZKukzSCZIel/SMpJ0L3uYQSQskvZT2o1Qz1sKPJM1X0n/81wred66kmSR3BNePZ0r6/s9KujxddhHJDYK/lvSjetsXFaekfpJ+l8YzX9J+6fKxkv6spHO9P9XcTS1pZPoeT6XxD5U0UHX79a89s5I0R9JVkhYAZ0raO41roaRZ2tTlxBxJP0mP1wuS9pF0p5K+9P+j4L2/VLD/6yVV1Px+JF0q6S+SHpO0vaRPAROBH6Xb7yzpW0rG4nha0oxi/kisDcryzkw/8vUA3gf+ieROza2B7wAXp+umAccWbps+jwdWkXT1242kj5Tvp+vOBK4qeP19JP+8DCW5C7M7cBpwYbpNN5K7lAel77saGNRAnB8n6UKiH0kHXw8DR6Xr5tDAHctbEOdtwP7p9CdIuh4hPS413VcfAvwunb6GpB8eSO4o7cHm3RIXHsc5wM/T6S7An4B+6fwk4IaC7S4viO9vBbFXkdyVPBz4A9Al3e7nwL+k0wEckU5fUXCM6/8e/8amu6u3KfffoB8te5Tj1NI6sIh4V9LNwLeAD4t82fxI+wyS9DJwf7r8GdI+V1K3R8RGYLGkpcCuJH2+jC4429iaJFGsAx6PiFca2N8+wJyIWJnu81aSAULuaoU4DwFGaFNnkf+kpGfarYGbJA0l+ZLtkq7/M/A9Sf1JOhpbrOY7mvxN+jwM2A14IH1NBUl3BTVq+pt5BniuIPalJJ2U7Q/sDcxPX9+DTR3XrSMZgwGS7rcPbSSWp4FbJd1F88fP2ignAsvCVSR9pNxYsKyatBQpqRPJf7811hZMbyyY30jdv9H6/aEESd8wZ0REnQ7AJI0nOSNoTcXE2QkYFxFr6sXzM2B2RBytZNyKOQARcZukecDngXvT0tZL1C3bdq8XR83PJZIv+E82E29hrIXxCrgpIi5o4LXrI/03H9hA498VnydJokeQJLRRsWm8Bmsn3EZgrS4i/kEyytKpBYuXkfz3CUmduQtb7jhJndJ6/GCSjrdmAd9Q0hU4knZRMrhNUx4HDpTUN62JTwEeaUE8DbkfOKNmRtIe6eTWbOoa+OSC9YOBpRFxNUnvoaOB14HtJPWR1I2kG+eGLAL6Sfpk+l5dJI3cglgfAo6VtF36+t6SdmrmNe+RDMlak9AHRMRs4Lz0Z+y1Bfu3NsKJwLLyXyR9qNf4JcmX719IhtZryX/rfyX5Ev8j8PX0v+5fkTQGP5E2sF5PM2e6aYnkfJIupf8CLIyIu5t6zRb4FjAmbTx9Hvh6uvwK4IeSnqwX3/HAs5KeIinz3BzJ+AeXkPysD5D0gNrQz7GOpMvry9Pj+hTwqWIDjYjngQuB+yU9ne5rh2ZeNgM4J/05hgL/LekZ4Eng6ohYVez+re1w76NmZjnnMwIzs5xzIjAzyzknAjOznHMiMDPLOScCM7OccyIwM8s5JwIzs5z7/8mPtxrGE99PAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "'''\n",
    "RUNNING EXPERIMENTS\n",
    "''' \n",
    "\n",
    "start_time = time.time()\n",
    "rmse_errors = []\n",
    "\n",
    "print('NUMBER OF OPERATORS | DERANDOMISED OPERATORS | AVERAGE RMSE ERROR\\n')\n",
    "measurement_range = [50, 100, 200, 600, 1000, 1600]\n",
    "\n",
    "for num_operators in measurement_range:   \n",
    "    \n",
    "    derandomized_hamiltonian = modified_derandomized_classical_shadow(observables_xyz, num_operators, \n",
    "                                                                      system_size, weight=absolute_coefficients)\n",
    "\n",
    "    tuples = (tuple(pauli) for pauli in derandomized_hamiltonian)\n",
    "    counts = Counter(tuples)\n",
    "\n",
    "    optimizer = SPSA(maxiter=2000) \n",
    "    cost_function = partial(objective_function, derandomized_hamiltonian)\n",
    "    \n",
    "    expectation_values = []\n",
    "    num_experiments = 3\n",
    "\n",
    "\n",
    "    for iteration in range(num_experiments):\n",
    "        \n",
    "        cost_history = []\n",
    "        \n",
    "        params = np.random.rand(ansatz.num_parameters)\n",
    "        result = optimizer.minimize(fun=cost_function, x0=params)\n",
    "        minimal_cost = min_cost()\n",
    "        expectation_values.append(minimal_cost)\n",
    "                \n",
    "        print(\"EXPERIMENT {}: GROUND STATE ENERGY FOUND = {}\".format(iteration, minimal_cost))\n",
    "\n",
    "\n",
    "    rmse_derandomised_cs = np.sqrt(np.sum([(expectation_values[i] - EXPECTED_EIGENVALUE)**2 \n",
    "                                           for i in range(num_experiments)])/num_experiments)\n",
    "\n",
    "    rmse_errors.append(rmse_derandomised_cs)\n",
    "    \n",
    "    print('{} | {} | {}'.format(num_operators, counts, rmse_derandomised_cs))\n",
    "\n",
    "\n",
    "elapsed_time = time.time() - start_time\n",
    "print(\"Execution time = \", time.strftime(\"%H:%M:%S\", time.gmtime(elapsed_time)))\n",
    "\n",
    "          \n",
    "points = measurement_range\n",
    "num_points = len(measurement_range)\n",
    "\n",
    "plt.plot([i for i in points], [rmse_errors[i] for i in range(num_points)], 'r', label='SPSA(maxiter=2000)')\n",
    "plt.xlabel('Number of measurements')\n",
    "plt.ylabel('Average RMSE error')\n",
    "plt.legend(loc=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d21f81e5",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}