04.quantum_chemistry.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "b9e89d96",
   "metadata": {},
   "source": [
    "# Exercise 4 - Quantum Chemistry. Electronic Structure Calculations \n",
    "\n",
    "\n",
    "## Background\n",
    "\n",
    "Quantum computing holds promise to revolutionize high-performance computing and possesses enormous potential for transforming various fields beyond the current capabilities of classical computing. Chemistry is one of the most promising applications of quantum computing. In particular, simulating the electronic structure of molecules to better understand chemical reactions and design new pharmaceutical drugs and materials, is one of the main near-term exciting applications for quantum computers.\n",
    "\n",
    "The Variational Quantum Eigensolver (VQE) is the main workhorse algorithm for quantum chemistry using near-term quantum computers. It is a hybrid quantum-classical algorithm that is used to find the ground state energy (lowest energy) of a molecular system. The algorithm was first published in 2014 by Alberto Peruzzo and coworkers [1] and used by IBM Quantum to simulate the ground state energy of small molecules in 2017.[2] Obtaining the ground state energies allows us to understand and study chemical reactions and gain additional insignt into various molecular properties, for instance the exact forces on nuclei than can serve to run molecular dynamics simulations to explore what happens in chemical systems with time.[3]\n",
    "\n",
    "VQE algorithm starts with a parameterized quantum circuit called an ansatz and a set of tunable parameters that get iteratively updated with a classical optimizers until the lowest energy is obtained. The VQE's advantage over classical algorithms comes from the fact that a quantum processing unit can represent and store the problem's wavefunction more efficiently than a classical computer. \n",
    "\n",
    "The quantum Equation Of Motion (qEOM) algorithm is an extension of the VQE to find the excitation energies (difference between ground and excited states energies, also called electronic transition energies) of the molecule. [4]\n",
    "\n",
    "In this exercise you will simulate a water molecule using VQE and qEOM to calculate its ground and excited state energies with a variety of tools available in Qiskit Nature.\n",
    "\n",
    "\n",
    "### References\n",
    "\n",
    "1. Peruzzo, Alberto, et al. \"A variational eigenvalue solver on a photonic quantum processor.\" Nature communications 5.1 (2014): 1-7.\n",
    "2. Kandala, Abhinav, et al. \"Hardware-efficient variational quantum eigensolver for small molecules and quantum magnets.\" Nature 549.7671 (2017): 242-246.\n",
    "3. Sokolov, Igor O., et al. \"Microcanonical and finite-temperature ab initio molecular dynamics simulations on quantum computers.\" Physical Review Research 3.1 (2021): 013125.\n",
    "4. Ollitrault, Pauline J., et al. \"Quantum equation of motion for computing molecular excitation energies on a noisy quantum processor.\" Physical Review Research 2.4 (2020): 043140.\n",
    "\n",
    "In addition, a review on VQE: Cerezo, M., et al. \"Variational quantum algorithms.\" Nat Rev Phys 3, 625–644 (2021). \n",
    "\n",
    "## Introduction\n",
    "\n",
    "One of the most important tasks when running a VQE algorithm on a noisy quantum computer is to reduce the computational cost by minimizing the quantum resources (number of qubits, quantum gates and measurements). This in turn, minimizes the computational errors. The goal is to construct the most compact and optimal quantum circuit capable of representing the ground state accurately. \n",
    "\n",
    "<div class=\"alert alert-block alert-success\">\n",
    "\n",
    "<b>Goal</b> \n",
    "\n",
    "We will first learn how to define an electronic structure problem with Qiskit Nature. We will then find the most efficient ansatz for approximating the ground state of a water molecule. Finally, we will use the ground state results to perform excited states calculations. \n",
    "    \n",
    "<b>Plan</b> \n",
    "        \n",
    "**Part 1: VQE for water molecule :** familiarize yourself with VQE and select the best combination of ansatz/classical optimizer by running statevector simulations.\n",
    "\n",
    "**Part 2: qEOM for water molecule :** extend the ground state calculation to find the electronic transition energies of this water molecule. Learn how to compare the quantum algorithm results with those obtained classically.\n",
    "\n",
    "</div>\n",
    "\n",
    "\n",
    "<div class=\"alert alert-block alert-danger\">\n",
    "\n",
    "Below is an introduction to the theory behind VQE simulations. You don't have to understand the whole thing before moving on. Don't be scared!\n",
    "\n",
    "</div>\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c0961d75",
   "metadata": {},
   "source": [
    "## Theory\n",
    "\n",
    "### Quantum and Classical\n",
    "\n",
    "\n",
    "Here is the general workflow representing how the molecular simulations using VQE are performed on quantum computers.\n",
    "\n",
    "<img src=\"resources/workflow.png\" width=600/>\n",
    "\n",
    "The core idea behind a hybrid quantum-classical approach is to divide the task between **CPU (classical processing unit)** and **QPU (quantum processing unit)** based on what they can do best. The CPU takes care of listing the terms that need to be measured to compute the energy and also optimizing the circuit parameters. The QPU implements a quantum circuit representing the quantum state of a system and measures the energy.\n",
    "\n",
    "For example, **CPU** can be used to find a set of functions (also called orbitals) which help re-expressing and simplifying the problem. In practice these orbitals are found through the [Hartree–Fock (HF) method](https://en.wikipedia.org/wiki/Hartree%E2%80%93Fock_method#:~:text=In%20computational%20physics%20and%20chemistry,system%20in%20a%20stationary%20state.). Further information on this method is provided below. \n",
    "\n",
    "After a HF calculation, the Hamiltonian is expressed in the basis of the resulting orbitals. It is then mapped to a qubit Hamiltonian using fermion-to-qubit transformations (see Hamiltonian section below). One can further analyze the properties of the system to reduce the number of qubits or shorten the ansatz circuit:\n",
    "\n",
    "- For Z2 symmetries and two-qubit reduction, see [Bravyi *et al*, 2017](https://arxiv.org/abs/1701.08213v1).\n",
    "- For entanglement forging, see [Eddins *et al.*, 2021](https://arxiv.org/abs/2104.10220v1).\n",
    "- For the adaptive ansatz see, [Grimsley *et al.*,2018](https://arxiv.org/abs/1812.11173v2), [Rattew *et al.*,2019](https://arxiv.org/abs/1910.09694), [Tang *et al.*,2019](https://arxiv.org/abs/1911.10205). You may use the ideas found in those works to find ways to shorten the quantum circuits.\n",
    "\n",
    "**QPU** is employed to implements quantum circuits (see Ansatzes section below), parameterized by angles $\\vec\\theta$. The goal is to find the set of these angles which best approximate the exact ground state wavefunction (or in other words, that minimizes the energy). Expectation values of different operators can then be measured on the trial wavefucntion prepared in the QPU.\n",
    "\n",
    "Below we go slightly more into mathematical details of each component of the VQE and qEOM algorithms. It might be also helpful if you watch our [video episode about VQE](https://www.youtube.com/watch?v=Z-A6G0WVI9w).\n",
    "\n",
    "\n",
    "### Hartree-Fock\n",
    "\n",
    "Hartree-Fock theory is based on the assumption that the many-electron wavefunction of the system can be described by a single Slater determinant made of one-electron molecular orbitals.  \n",
    "\n",
    "Electrons are indistinguishable fermions. This means that any permutation of two electrons must change the sign of the wavefunction (this is a really important constraint called the Pauli exclusion principle, it is the\n",
    "reason for the specific structure of the periodic table).\n",
    "\n",
    "We can always write the electronic wavefunction as the following expansion:\n",
    "\n",
    "\\begin{equation}\n",
    "\\psi_{\\text{elec}}(\\mathbf{r}_1, \\mathbf{r}_2, ...,\\mathbf{r}_n) = \\sum_{m_1, m_2,...,m_N} C_{m_1, m_2, ..., m_N} | \\phi_{m_1}(\\mathbf{r}_1) \\phi_{m_2}(\\mathbf{r}_2) ... \\phi_{m_n}(\\mathbf{r}_n) |\n",
    "\\end{equation}\n",
    "\n",
    "where $C_{m_1, m_2, ..., m_N}$ are complex coefficients, $| \\phi_{m_1}(\\mathbf{r}_1) \\phi_{m_2}(\\mathbf{r}_2) ... \\phi_{m_n}(\\mathbf{r}_n) |$ are Slater determinants and the $\\phi_{m_i}(\\mathbf{r_i})$ are one-electron functions. These one-electron functions are referred to as <b>molecular orbitals (MOs)</b>. MOs are computed with the famous <b>Hartree-Fock (HF)</b> method. Usually they are expanded in terms of <b> atomic orbitals (AOs) </b>:\n",
    "\n",
    "\\begin{equation}\n",
    "\\phi(\\mathbf{r}, s) = \\sum_n D_{mn}\\chi_n(\\mathbf{r})\\otimes s\n",
    "\\end{equation}\n",
    "\n",
    "($s$ is the spin variable $s \\in \\{ \\alpha, \\beta \\}$).\n",
    "The MOs can host two electrons of opposite spin (either spin up ($\\alpha$) or spin down ($\\beta$)) and can either be filled (occupied) or empty (virtual). This spin component of a molecular orbital is also referred to as a spin orbital (SO). When mapping chemistry problems on quantum computers using Jordan-Wigner mapping, one qubit represents one spin orbital.\n",
    "\n",
    "The Hartree-Fock method optimizes the coefficients $D_{mn}$ such that the electronic wavefunction expressed as a <b>single Slater determinant</b> ($\\psi_{\\text{elec}}(\\mathbf{r}_1, \\mathbf{r}_2, ...,\\mathbf{r}_n) = | \\phi_{m_1}(\\mathbf{r}_1) \\phi_{m_2}(\\mathbf{r}_2) ... \\phi_{m_n}(\\mathbf{r}_n) |$) minimizes the energy. It often provides a good starting point for more elaborate theoretical methods which are better approximations to the solution of the electronic Schrödinger equation. The Hartree- Fock method determines the set of MOs which minimize the energy and give us this \"best single determinant.\"\n",
    "\n",
    "The AOs are functions centered on the atoms. They are themselves expanded in terms of the functions composing the <b> basis set </b>. The basis functions used in quantum chemsitry are taken as sums of Gaussians centered on the atoms times a polynomial function, $P_l$, in the electron coordinates relative to that atom:\n",
    "\n",
    "\\begin{equation}\n",
    "\\chi_n(\\mathbf{r}) = \\sum_l N_l^n \\exp\\big( - \\alpha_l^n (|\\mathbf{r} - \\mathbf{R}_I^n|^2)\\big) P_l(\\mathbf{r} - \\mathbf{R}_I^n)\n",
    "\\end{equation}\n",
    "\n",
    "The choice of the $\\alpha_l^n$ and of the order of the polynomial determine the <b> basis set </b> (that take different names. Famous ones are: 6-31G, cc-pVQZ, TZ2P...). \n",
    "\n",
    "The basis set (together with the number of particles) determines the total number of AOs. For the same molecule different basis sets lead to different number of AOs and it turn of MOs. \n",
    "\n",
    "The type of basis functions used also influences the accuracy.\n",
    "The better a single basis function is able to reproduce the unknown function, the fewer basis functions necessary for achieving a given level of accuracy. The smaller the basis, the poorer the representation of the MO. \n",
    "\n",
    "\n",
    "### Molecular Hamiltonian \n",
    "\n",
    "The Hamiltonian of a system is an operator corresponding to the total energy of that system. The molecular Hamiltonian defined as:\n",
    "$$\n",
    "\\hat{H} =\\sum_{r s} h_{r s} \\hat{a}_{r}^{\\dagger} \\hat{a}_{s}\n",
    "+\\frac{1}{2} \\sum_{p q r s} g_{p q r s} \\hat{a}_{p}^{\\dagger} \\hat{a}_{q}^{\\dagger} \\hat{a}_{r} \\hat{a}_{s}+E_{N N}\n",
    "$$\n",
    "where\n",
    "$$\n",
    "h_{p q}=\\int \\phi_{p}^{*}(r)\\left(-\\frac{1}{2} \\nabla^{2}-\\sum_{I} \\frac{Z_{I}}{R_{I}-r}\\right) \\phi_{q}(r)\n",
    "$$\n",
    "$$\n",
    "g_{p q r s}=\\int \\frac{\\phi_{p}^{*}\\left(r_{1}\\right) \\phi_{q}^{*}\\left(r_{2}\\right) \\phi_{r}\\left(r_{2}\\right) \\phi_{s}\\left(r_{1}\\right)}{\\left|r_{1}-r_{2}\\right|} \n",
    "$$\n",
    "\n",
    "The $h_{r s}$  and $g_{p q r s}$ are the one-/two-body integrals (Hamiltonian terms expressed in the orbitals resulting from the HF calculation) and $E_{N N}$ the nuclear repulsion energy. \n",
    "The one-body integrals represent the kinetic energy of the electrons and their interaction with nuclei. \n",
    "The two-body integrals represent the electron-electron interaction.\n",
    "The $\\hat{a}_{r}^{\\dagger} / \\hat{a}_{r}$ operators create / annihilate an electron in spin orbital $r$. They need to be mapped to Pauli operators, so that we can measure their expectation values on a quantum computer.\n",
    "Note that VQE minimizes the electronic energy so you have to retrieve and add the nuclear repulsion energy $E_{NN}$ to compute the total energy. \n",
    "    \n",
    "### Hamiltonian Mapping \n",
    "\n",
    "For every non-zero matrix element in the $ h_{r s}$ and $g_{p q r s}$ tensors, we can construct corresponding Pauli string (tensor product of Pauli operators) with the following fermion-to-qubit transformation. \n",
    "For instance, in Jordan-Wigner mapping for an orbital $r = 3$, we obtain the following Pauli string:\n",
    "$$\n",
    "\\hat a_{3}^{\\dagger}= \\hat \\sigma_z \\otimes \\hat \\sigma_z \\otimes\\left(\\frac{ \\hat \\sigma_x-i \\hat \\sigma_y}{2}\\right) \\otimes 1 \\otimes \\cdots \\otimes 1\n",
    "$$\n",
    "where $\\hat \\sigma_x, \\hat \\sigma_y, \\hat \\sigma_z$ are the well-known Pauli operators. The tensor products of $\\hat \\sigma_z$ operators are placed to enforce the fermionic anti-commutation relations.\n",
    "A representation of the Jordan-Wigner mapping between the 14 spin-orbitals of a water molecule and some 14 qubits is given below:\n",
    "\n",
    "<img src=\"resources/mapping.png\" width=600/>\n",
    "\n",
    "\n",
    "Then, one simply replaces the one-/two-body excitations (e.g. $\\hat{a}_{r}^{\\dagger} \\hat{a}_{s}$, $\\hat{a}_{p}^{\\dagger} \\hat{a}_{q}^{\\dagger} \\hat{a}_{r} \\hat{a}_{s}$) in the Hamiltonian by corresponding Pauli strings (i.e. $\\hat{P}_i$, see picture above). The resulting operator set is ready to be measured on the QPU.\n",
    "For additional details see [Seeley *et al.*, 2012](https://arxiv.org/abs/1208.5986v1).\n",
    "\n",
    "\n",
    "\n",
    "### Active Space Selection\n",
    "\n",
    "\n",
    "For some molecular systems, not all orbitals need to be considered to capture most of the essential chemistry. In fact, we can use our chemical intuition and select a subset of orbitals that play the most important role in determining the energy of our system. This concept is known as an <b>active space</b> selection and it allows us to reduce the size and cost of our quantum chemistry problem. \n",
    "\n",
    "As illustrated in the figure below, instead of simulating all orbitals on a quantum computer (left side) we construct an active space (right side) with some orbitals that are simulated on a quantum computer using VQE (colored in blue), while the rest are approximated using classical methods (colored in yellow).\n",
    "\n",
    "</div>\n",
    "<div style=\"float: right; width: 50%; padding-right=20px;\">\n",
    "    <img src=\"resources/active_space.png\" width=\"300\"/>\n",
    "</div>\n",
    "\n",
    "\n",
    "In general, selecting the best active space is a challenging task. However, the choice of the active space is crucial to obtain reliable wavefunctions and corresponding energies.\n",
    "\n",
    "\n",
    "### Ansätze\n",
    "\n",
    "There are mainly 2 types of ansätze you can use for chemical problems. \n",
    "\n",
    "- **q-UCC ansätze** are physically inspired, and roughly map the electron excitations to quantum circuits. The q-UCCSD ansatz (`UCCSD`in Qiskit) possess all possible single and double electron excitations. The paired double q-pUCCD (`PUCCD`) and singlet q-UCCD0 (`SUCCD`) just consider a subset of such excitations (meaning significantly shorter circuits) and have proved to provide good results for dissociation profiles. For instance, q-pUCCD doesn't have single excitations and the double excitations are paired as in the image below.\n",
    "- **Heuristic ansätze (`TwoLocal`)** were invented to shorten the circuit depth but still be able to represent the ground state. \n",
    "As in the figure below, the R gates represent the parametrized single qubit rotations and $U_{CNOT}$ the entanglers (two-qubit gates). The idea is that after repeating certain $D$-times the same block (with independent parameters) one can reach the ground state. \n",
    "\n",
    "For additional details refer to [Sokolov *et al.* (q-UCC ansatzes)](https://arxiv.org/abs/1911.10864v2) and [Barkoutsos *et al.* (Heuristic ansatzes)](https://arxiv.org/pdf/1805.04340.pdf).\n",
    "\n",
    "<img src=\"resources/ansatz.png\" width=700/>\n",
    "\n",
    "\n",
    "\n",
    "### Putting it all together: Variational Quantum Eigensolver (VQE)\n",
    "\n",
    "Given a Hermitian operator $\\hat H$ with an unknown minimum eigenvalue $E_{min}$, associated with the eigenstate $|\\psi_{min}\\rangle$, VQE provides an estimate $E_{\\theta}$, bounded by $E_{min}$:\n",
    "\n",
    "\\begin{align*}\n",
    "    E_{min} \\le E_{\\theta} \\equiv \\langle \\psi(\\theta) |\\hat H|\\psi(\\theta) \\rangle\n",
    "\\end{align*}  \n",
    "\n",
    "where $|\\psi(\\theta)\\rangle$ is the trial state associated with $E_{\\theta}$. By applying a parameterized circuit, represented by $U(\\theta)$, to some arbitrary starting state $|\\psi\\rangle$, the algorithm obtains an estimate $U(\\theta)|\\psi\\rangle \\equiv |\\psi(\\theta)\\rangle$ on $|\\psi_{min}\\rangle$. The estimate is iteratively optimized by a classical optimizer by changing the parameter $\\theta$ and minimizing the expectation value of $\\langle \\psi(\\theta) |\\hat H|\\psi(\\theta) \\rangle$. \n",
    "\n",
    "Besides finding the ground state energies, VQE can be used in molecular dynamics simulations, see [Sokolov *et al.*, 2021](https://arxiv.org/abs/2008.08144v1), and excited states calculations, see [Ollitrault *et al.*, 2019](https://arxiv.org/abs/1910.12890) to name a few.\n",
    "\n",
    "<div class=\"alert alert-block alert-danger\">\n",
    " \n",
    "<b> References for additional details</b> \n",
    "\n",
    "For the Qiskit Nature tutorials that implement this algorithm see [here](https://qiskit.org/documentation/nature/tutorials/01_electronic_structure.html).\n",
    "For additional information, please refer to Qiskit Github [first page of github repository](https://github.com/Qiskit/qiskit-nature) >>> the [test folder](https://github.com/Qiskit/qiskit-nature/tree/main/test) with the base code for the use of each functionality.\n",
    "    \n",
    "\n",
    "</div>\n",
    "\n",
    "\n",
    "### Quantum equation of motion (qEOM)\n",
    "\n",
    "The excitation energies can be solved using the quantum equation of motion (qEOM) algorithm [Phys. Rev. Research 2, 043140]. The EOM method finds the excitation energies (differences in energy between the ground state and all $n$th excited states) by solving the following pseudo-eigenvalue problem.\n",
    "\n",
    "$$\n",
    "\\begin{pmatrix}\n",
    "    \\text{M} & \\text{Q}\\\\ \n",
    "    \\text{Q*} & \\text{M*}\n",
    "\\end{pmatrix}\n",
    "\\begin{pmatrix}\n",
    "    \\text{X}_n\\\\ \n",
    "    \\text{Y}_n\n",
    "\\end{pmatrix}\n",
    "= E_{0n}\n",
    "\\begin{pmatrix}\n",
    "    \\text{V} & \\text{W}\\\\ \n",
    "    -\\text{W*} & -\\text{V*}\n",
    "\\end{pmatrix}\n",
    "\\begin{pmatrix}\n",
    "    \\text{X}_n\\\\ \n",
    "    \\text{Y}_n\n",
    "\\end{pmatrix}\n",
    "$$\n",
    "\n",
    "with \n",
    "\n",
    "$$\n",
    "M_{\\mu_{\\alpha}\\nu_{\\beta}} = \\langle0| [(\\hat{\\text{E}}_{\\mu_{\\alpha}}^{(\\alpha)})^{\\dagger},\\hat{\\text{H}}, \\hat{\\text{E}}_{\\nu_{\\beta}}^{(\\beta)}]|0\\rangle\n",
    "$$\n",
    "$$\n",
    "Q_{\\mu_{\\alpha}\\nu_{\\beta}} = -\\langle0| [(\\hat{\\text{E}}_{\\mu_{\\alpha}}^{(\\alpha)})^{\\dagger}, \\hat{\\text{H}}, (\\hat{\\text{E}}_{\\nu_{\\beta}}^{(\\beta)})^{\\dagger}]|0\\rangle\n",
    "$$\n",
    "$$\n",
    "V_{\\mu_{\\alpha}\\nu_{\\beta}} = \\langle0| [(\\hat{\\text{E}}_{\\mu_{\\alpha}}^{(\\alpha)})^{\\dagger}, \\hat{\\text{E}}_{\\nu_{\\beta}}^{(\\beta)}]|0\\rangle\n",
    "$$\n",
    "$$\n",
    "W_{\\mu_{\\alpha}\\nu_{\\beta}} = -\\langle0| [(\\hat{\\text{E}}_{\\mu_\\alpha}^{(\\alpha)})^{\\dagger}, (\\hat{\\text{E}}_{\\nu_{\\beta}}^{(\\beta)})^{\\dagger}]|0\\rangle\n",
    "$$\n",
    "\n",
    "Note, each matrix element must be measured on the quantum computer with the corresponding ground state. In order to use the qEOM as a solver in Qiskit, we have to run a ground state calculation first and use our results as an input to the qEOM solver.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6bf1ba48",
   "metadata": {},
   "source": [
    "## Part 1: VQE for Water Molecule \n",
    "\n",
    "\n",
    "In part 1, you will simulate H$_2$O molecule using the cc-pVDZ basis with the PySCF driver and Jordan-Wigner mapping.\n",
    "\n",
    "#### 1. Chemistry Drivers\n",
    "\n",
    "The interfaces to the classical chemistry codes that are available in Qiskit are called drivers.\n",
    "Qiskit Nature supports a few of them, namely `PSI4Driver`, `PyQuanteDriver` and `PySCFDriver`. \n",
    "\n",
    "The first step is to define a water molecule by specifying its atomic coordinates, total spin (multiplicity) and charge information. The second step is to perform a Hartree-Fock calculation and obtain all the necessary information about our molecule that we will then use to perform VQE. We implement the second step by running the PySCF driver."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "id": "a184f130",
   "metadata": {},
   "outputs": [],
   "source": [
    "from qiskit_nature.drivers import Molecule\n",
    "from qiskit_nature.drivers.second_quantization import ElectronicStructureMoleculeDriver, ElectronicStructureDriverType\n",
    "\n",
    "molecule = Molecule(\n",
    "    # coordinates are given in Angstrom\n",
    "    geometry=[\n",
    "        [\"O\", [0.0, 0.0, 0.0]],\n",
    "        [\"H\", [0.758602, 0.0, 0.504284]],\n",
    "        [\"H\", [0.758602, 0.0, -0.504284]]\n",
    "    ],\n",
    "    multiplicity=1,  # = 2*spin + 1\n",
    "    charge=0,\n",
    ")\n",
    "\n",
    "driver = ElectronicStructureMoleculeDriver(\n",
    "    molecule=molecule,\n",
    "    basis=\"sto3g\",\n",
    "    driver_type=ElectronicStructureDriverType.PYSCF,\n",
    ")\n",
    "\n",
    "properties = driver.run()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "93efcca8",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-block alert-danger\">\n",
    "    \n",
    "<b> Tutorial questions 1</b> \n",
    "    \n",
    "Look into the properties of the molecule (see also: https://qiskit.org/documentation/nature/tutorials/08_property_framework.html) and answer the questions below.\n",
    "\n",
    "    \n",
    "1. What is the total number of (alpha and beta) electrons in your system? (Coding)\n",
    "2. What is the number of molecular and spin orbitals? (Coding)\n",
    "3. How many qubits would you need to simulate this molecule with Jordan-Wigner mapping? (Discussion in Slack)\n",
    "4. What is the value of the nuclear repulsion energy? (Coding)\n",
    "\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "id": "2423820b-3526-469e-9cb6-fa35debef946",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "9.819476509782405\n",
      "(5, 5) 14\n"
     ]
    }
   ],
   "source": [
    "from qiskit_nature.properties.second_quantization.electronic import ElectronicEnergy\n",
    "electronic_energy = properties.get_property(ElectronicEnergy)\n",
    "print(electronic_energy.nuclear_repulsion_energy)\n",
    "from qiskit_nature.problems.second_quantization import ElectronicStructureProblem\n",
    "problem = ElectronicStructureProblem(driver)\n",
    "second_q_ops = problem.second_q_ops()\n",
    "main_op = second_q_ops[0]\n",
    "particle_number = problem.grouped_property_transformed.get_property(\"ParticleNumber\")\n",
    "print(particle_number.num_particles, particle_number.num_spin_orbitals)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "id": "64353d56",
   "metadata": {},
   "outputs": [],
   "source": [
    "# WRITE YOUR CODE BETWEEN THESE LINES - START\n",
    "\n",
    "\n",
    "num_alpha_electrons = 5\n",
    "\n",
    "num_beta_electrons = 5\n",
    "\n",
    "num_spin_orbitals = 14\n",
    "\n",
    "nuclear_rep_energy = 9.819476509782405\n",
    "\n",
    "\n",
    "# WRITE YOUR CODE BETWEEN THESE LINES - END"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 51,
   "id": "cee64ade-5a50-4c9d-bfeb-ef83c0fd8c54",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Submitting your answer. Please wait...\n",
      "Congratulations 🎉! Your answer is correct and has been submitted.\n"
     ]
    }
   ],
   "source": [
    "from qc_grader.challenges.spring_2022 import grade_ex4a\n",
    "\n",
    "grade_ex4a(\n",
    "    num_alpha_electrons=num_alpha_electrons,\n",
    "    num_beta_electrons=num_beta_electrons,\n",
    "    num_spin_orbitals=num_spin_orbitals,\n",
    "    nuclear_rep_energy=nuclear_rep_energy\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "814a0c47",
   "metadata": {},
   "source": [
    "#### 2. Electronic Structure Problem and Active Space Transformer\n",
    "\n",
    "\n",
    "The next step is to select our active space which will describe what molecular orbitals will be simulated with VQE and treat the rest with classical methods. The more orbitals we include in our active space, the more computationally expensive our quantum circuit will be. As such, the goal is to construct the most accurate activate space with as few orbitals as possible. In Qiskit Nature, we use  `ActiveSpaceTransformer` to specify how many electrons and molecular orbitals we would like to have in our active space, as well as which of those particular orbitals we would like to consider."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 52,
   "id": "12a41c9d",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ParticleNumber:\n",
      "\t14 SOs\n",
      "\t5 alpha electrons\n",
      "\t\torbital occupation: [1. 1. 1. 1. 1. 0. 0.]\n",
      "\t5 beta electrons\n",
      "\t\torbital occupation: [1. 1. 1. 1. 1. 0. 0.]\n"
     ]
    }
   ],
   "source": [
    "from qiskit_nature.transformers.second_quantization.electronic import ActiveSpaceTransformer\n",
    "\n",
    "# Check the occupation of the spin orbitals\n",
    "PN_property = properties.get_property(\"ParticleNumber\")\n",
    "print(PN_property)\n",
    "\n",
    "# Define the active space around the Fermi level \n",
    "# (selected automatically around the HOMO and LUMO, ordered by energy)\n",
    "transformer = ActiveSpaceTransformer(\n",
    "    num_electrons=2, #how many electrons we have in our active space\n",
    "    num_molecular_orbitals=2, #how many orbitals we have in our active space\n",
    ")\n",
    "\n",
    "\n",
    "# We can hand-pick the MOs to be included in the AS\n",
    "# (in case they are not exactly around the Fermi level)\n",
    "\n",
    "# transformer = ActiveSpaceTransformer(\n",
    "#     num_electrons=2, #how many electrons we have in our active space\n",
    "#     num_molecular_orbitals=2, #how many orbitals we have in our active space\n",
    "#     active_orbitals=[4,7], #We are specifically choosing MO number 4 (occupied with two electrons) and 7 (empty)\n",
    "# )"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "904fabad",
   "metadata": {},
   "source": [
    "You can then create an `ElectronicStructureProblem` that produces a list of fermionic operators before mapping them to qubits (Pauli strings)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "id": "1629ce0d",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Fermionic Operator\n",
      "register length=4, number terms=14\n",
      "  (0.038537319755372834+0j) * ( +_0 -_1 +_2 -_3 )\n",
      "+ (-0.03853731975537283+0j) * ( +_0 -_1 -_2 +_3 )\n",
      "+ (-0.03853731975537283+0j) * ( -_0 +_1 +_2 -_3 )\n",
      "+ (0.03853731975537284+0j) * ( -_0 +_1 -_2 +_3 ...\n"
     ]
    }
   ],
   "source": [
    "from qiskit_nature.problems.second_quantization.electronic import ElectronicStructureProblem\n",
    "\n",
    "problem = ElectronicStructureProblem(driver, [transformer])\n",
    "\n",
    "second_q_ops = problem.second_q_ops() # this calls driver.run() internally\n",
    "\n",
    "hamiltonian = second_q_ops[0]\n",
    "\n",
    "print(hamiltonian)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6b2b03ab",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-block alert-danger\">\n",
    "    \n",
    "<b> Tutorial questions 2</b> \n",
    "    \n",
    "Create the following operator as a FermionicOperator with Qiskit and check if it is Hermitian.\n",
    "    \n",
    "</div>\n",
    "\n",
    "\\begin{align*}\n",
    "O = 0.5i a^{\\dagger}_1 a_0 + 2.0 a^{\\dagger}_2 a_0 \n",
    "\\end{align*}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 88,
   "id": "3150e803",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Fermionic Operator\n",
      "register length=3, number terms=2\n",
      "  0.5j * ( -_0 +_1 )\n",
      "+ (2+0j) * ( -_0 +_2 )\n"
     ]
    }
   ],
   "source": [
    "from qiskit_nature.operators.second_quantization import FermionicOp\n",
    "\n",
    "# WRITE YOUR CODE BETWEEN THESE LINES - START\n",
    "\n",
    "\n",
    "# list of valid tuples to initialize the FermionicOp\n",
    "list_operator = [('-_0 +_1', (0.+0.5j)),('-_0 +_2', (2.0+0j))]\n",
    "\n",
    "# positive integer that represents the length of registers\n",
    "num_register_length = 3\n",
    "print(FermionicOp(list_operator, display_format='sparse'))\n",
    "\n",
    "# WRITE YOUR CODE BETWEEN THESE LINES - END"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 89,
   "id": "6987d9bf-0d6f-4c73-93b9-0b7dbe73f657",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Submitting your answer. Please wait...\n",
      "Congratulations 🎉! Your answer is correct and has been submitted.\n"
     ]
    }
   ],
   "source": [
    "from qc_grader.challenges.spring_2022 import grade_ex4b\n",
    "\n",
    "grade_ex4b(list_operator, num_register_length)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "eaf627a0",
   "metadata": {},
   "source": [
    "#### 3. Fermion-Qubit Mapping \n",
    "\n",
    "Because electrons are fermions, the electronic systems are described by Hamiltonians consisting of fermionic operators expressed in second quantization. Since quantum computers are made up of qubits, we have to transform such fermionic Hamiltonians into the qubit operators. There are different types of mappers to transform fermionic operators to qubits operators. You can try different mapping options, but we will stick to `JordanWignerMapper`, which provides a simple correspondence: one qubit represents one spin orbital."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 90,
   "id": "73830efb",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.009634329938843209 * YYYY\n",
      "+ 0.009634329938843209 * XXYY\n",
      "+ 0.009634329938843209 * YYXX\n",
      "+ 0.009634329938843209 * XXXX\n",
      "- 0.8720633572288636 * IIII\n",
      "- 0.17865261550743655 * ZIII\n",
      "+ 0.11923749314031712 * IZII\n",
      "+ 0.1475343334568582 * ZZII\n",
      "- 0.1786526155074366 * IIZI\n",
      "+ 0.16144783491422293 * ZIZI\n",
      "+ 0.1571686633957014 * IZZI\n",
      "+ 0.11923749314031729 * IIIZ\n",
      "+ 0.1571686633957014 * ZIIZ\n",
      "+ 0.2200397733437611 * IZIZ\n",
      "+ 0.1475343334568582 * IIZZ\n"
     ]
    }
   ],
   "source": [
    "from qiskit_nature.mappers.second_quantization import ParityMapper, BravyiKitaevMapper, JordanWignerMapper\n",
    "from qiskit_nature.converters.second_quantization import QubitConverter\n",
    "\n",
    "# Setup the mapper and qubit converter\n",
    "mapper_type = 'JordanWignerMapper'\n",
    "\n",
    "if mapper_type == 'ParityMapper':\n",
    "    mapper = ParityMapper()\n",
    "elif mapper_type == 'JordanWignerMapper':\n",
    "    mapper = JordanWignerMapper()\n",
    "elif mapper_type == 'BravyiKitaevMapper':\n",
    "    mapper = BravyiKitaevMapper()\n",
    "\n",
    "converter = QubitConverter(mapper)\n",
    "\n",
    "qubit_op = converter.convert(hamiltonian)\n",
    "print(qubit_op)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "eb3295c1",
   "metadata": {},
   "source": [
    "#### 4. Initial State\n",
    "As we described in the Theory section, a good initial state in chemistry is the HF state (i.e. $|\\Psi_{HF} \\rangle = |0101 \\rangle$). We can initialize it as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 91,
   "id": "f9e3af40",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "     ┌───┐\n",
      "q_0: ┤ X ├\n",
      "     └───┘\n",
      "q_1: ─────\n",
      "     ┌───┐\n",
      "q_2: ┤ X ├\n",
      "     └───┘\n",
      "q_3: ─────\n",
      "          \n"
     ]
    }
   ],
   "source": [
    "from qiskit_nature.circuit.library import HartreeFock\n",
    "\n",
    "particle_number = problem.grouped_property_transformed.get_property(\"ParticleNumber\")\n",
    "num_spin_orbitals = particle_number.num_spin_orbitals\n",
    "num_particles = particle_number.num_particles\n",
    "\n",
    "\n",
    "init_state = HartreeFock(num_spin_orbitals, num_particles, converter)\n",
    "print(init_state)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "103eb34f",
   "metadata": {},
   "source": [
    "#### 5. Ansatz\n",
    "One of the most important choices in the quantum circuit design is ansatz selection which is used to approximate the ground state. Here is the example of a qiskit circuit library that contains many possibilities for choosing different ansatze."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 92,
   "id": "d4015bc5",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "     ┌───┐┌───────────────────────────────────┐»\n",
      "q_0: ┤ X ├┤0                                  ├»\n",
      "     └───┘│                                   │»\n",
      "q_1: ─────┤1                                  ├»\n",
      "     ┌───┐│  exp(-it (IIXY + IIYX))(1.0*t[0]) │»\n",
      "q_2: ┤ X ├┤2                                  ├»\n",
      "     └───┘│                                   │»\n",
      "q_3: ─────┤3                                  ├»\n",
      "          └───────────────────────────────────┘»\n",
      "«     ┌───────────────────────────────────┐»\n",
      "«q_0: ┤0                                  ├»\n",
      "«     │                                   │»\n",
      "«q_1: ┤1                                  ├»\n",
      "«     │  exp(-it (XYII + YXII))(1.0*t[1]) │»\n",
      "«q_2: ┤2                                  ├»\n",
      "«     │                                   │»\n",
      "«q_3: ┤3                                  ├»\n",
      "«     └───────────────────────────────────┘»\n",
      "«     ┌─────────────────────────────────────────────────────────────────────────────┐\n",
      "«q_0: ┤0                                                                            ├\n",
      "«     │                                                                             │\n",
      "«q_1: ┤1                                                                            ├\n",
      "«     │  exp(-it (YYXY + XXXY + XYYY + YXYY + XYXX + YXXX + YYYX + XXYX))(1.0*t[2]) │\n",
      "«q_2: ┤2                                                                            ├\n",
      "«     │                                                                             │\n",
      "«q_3: ┤3                                                                            ├\n",
      "«     └─────────────────────────────────────────────────────────────────────────────┘\n"
     ]
    }
   ],
   "source": [
    "from qiskit.circuit.library import TwoLocal\n",
    "from qiskit_nature.circuit.library import UCCSD, PUCCD, SUCCD\n",
    "\n",
    "# Choose the ansatz\n",
    "ansatz_type = \"UCCSD\"\n",
    "\n",
    "\n",
    "# Put arguments for twolocal\n",
    "if ansatz_type == \"TwoLocal\":\n",
    "    # Single qubit rotations that are placed on all qubits with independent parameters\n",
    "    rotation_blocks = ['ry', 'rz']\n",
    "    # Entangling gates\n",
    "    entanglement_blocks = 'cz'\n",
    "    # How the qubits are entangled \n",
    "    entanglement = 'full'\n",
    "    # Repetitions of rotation_blocks + entanglement_blocks with independent parameters\n",
    "    repetitions = 1\n",
    "    # Skip the final rotation_blocks layer\n",
    "    skip_final_rotation_layer = True\n",
    "    ansatz = TwoLocal(qubit_op.num_qubits, rotation_blocks, entanglement_blocks, reps=repetitions, \n",
    "                      entanglement=entanglement, skip_final_rotation_layer=skip_final_rotation_layer)\n",
    "    # Add the initial state\n",
    "    ansatz.compose(init_state, front=True, inplace=True)\n",
    "elif ansatz_type == \"UCCSD\":\n",
    "    ansatz = UCCSD(converter,num_particles,num_spin_orbitals,initial_state = init_state)\n",
    "elif ansatz_type == \"PUCCD\":\n",
    "    ansatz = PUCCD(converter,num_particles,num_spin_orbitals,initial_state = init_state)\n",
    "elif ansatz_type == \"SUCCD\":\n",
    "    ansatz = SUCCD(converter,num_particles,num_spin_orbitals,initial_state = init_state)\n",
    "elif ansatz_type == \"Custom\":\n",
    "    # Example of how to write your own circuit\n",
    "    from qiskit.circuit import Parameter, QuantumCircuit, QuantumRegister\n",
    "    # Define the variational parameter\n",
    "    theta = Parameter('a')\n",
    "    n = qubit_op.num_qubits\n",
    "    # Make an empty quantum circuit\n",
    "    qc = QuantumCircuit(qubit_op.num_qubits)\n",
    "    qubit_label = 0\n",
    "    # Place a Hadamard gate\n",
    "    qc.h(qubit_label)\n",
    "    # Place a CNOT ladder\n",
    "    for i in range(n-1):\n",
    "        qc.cx(i, i+1)\n",
    "    # Visual separator\n",
    "    qc.barrier()\n",
    "    # rz rotations on all qubits\n",
    "    qc.rz(theta, range(n))\n",
    "    ansatz = qc\n",
    "    ansatz.compose(init_state, front=True, inplace=True)\n",
    "\n",
    "print(ansatz.decompose())"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f383eb00",
   "metadata": {},
   "source": [
    "#### 6. Backend\n",
    "\n",
    "To run a quantum circuit, we have to specify a quantum backend which is either a simulator or a real quantum computer. We will focus on the `statevector_simulator` in this challenge."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 93,
   "id": "10c4e4aa",
   "metadata": {},
   "outputs": [],
   "source": [
    "from qiskit import Aer\n",
    "backend = Aer.get_backend('statevector_simulator')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "691f2f5a",
   "metadata": {},
   "source": [
    "#### 7. Optimizer\n",
    "\n",
    "The optimizer guides the evolution of the parameters of the ansatz. A good optimizer choice might help to reduce the number of energy evaluations needed. As such, it is very important to investigate the energy convergence as it would define the number of measurements needed to be performed on the QPU."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 94,
   "id": "41eca79e",
   "metadata": {},
   "outputs": [],
   "source": [
    "from qiskit.algorithms.optimizers import COBYLA, L_BFGS_B, SPSA, SLSQP\n",
    "\n",
    "optimizer_type = 'L_BFGS_B'\n",
    "\n",
    "# You may want to tune the parameters \n",
    "# of each optimizer, here the defaults are used\n",
    "if optimizer_type == 'COBYLA':\n",
    "    optimizer = COBYLA(maxiter=500)\n",
    "elif optimizer_type == 'L_BFGS_B':\n",
    "    optimizer = L_BFGS_B(maxfun=500)\n",
    "elif optimizer_type == 'SPSA':\n",
    "    optimizer = SPSA(maxiter=500)\n",
    "elif optimizer_type == 'SLSQP':\n",
    "    optimizer = SLSQP(maxiter=500)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ed285285",
   "metadata": {},
   "source": [
    "#### 8. Exact eigensolver\n",
    "\n",
    "\n",
    "Exact diagonalization is a numerical technique used to determine the eigenstates and energy eigenvalues of a quantum Hamiltonian represented as a matrix. Since the dimensions of this matrix scale exponentially with the number of molecular orbitals we can only use this method for small molecular systems. For very large systems you would run out of memory trying to store their wavefunctions.\n",
    "\n",
    "Since our  molecular problem is small, we can use exact diagonalization and solve the problem exactly. This gives us a great reference to compare with the VQE result."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 95,
   "id": "c75d1521",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Exact electronic energy -1.697007028852895\n",
      "=== GROUND STATE ENERGY ===\n",
      " \n",
      "* Electronic ground state energy (Hartree): -84.700895930366\n",
      "  - computed part:      -1.697007028853\n",
      "  - ActiveSpaceTransformer extracted energy part: -83.003888901513\n",
      "~ Nuclear repulsion energy (Hartree): 9.819476509782\n",
      "> Total ground state energy (Hartree): -74.881419420584\n",
      " \n",
      "=== MEASURED OBSERVABLES ===\n",
      " \n",
      "  0:  # Particles: 2.000 S: 0.000 S^2: 0.000 M: 0.000\n",
      " \n",
      "=== DIPOLE MOMENTS ===\n",
      " \n",
      "~ Nuclear dipole moment (a.u.): [2.86710004  0.0  0.0]\n",
      " \n",
      "  0: \n",
      "  * Electronic dipole moment (a.u.): [2.0649676  0.0  0.0]\n",
      "    - computed part:      [0.0022825  0.0  0.0]\n",
      "    - ActiveSpaceTransformer extracted energy part: [2.0626851  0.0  0.0]\n",
      "  > Dipole moment (a.u.): [0.80213244  0.0  0.0]  Total: 0.80213244\n",
      "                 (debye): [2.0388171  0.0  0.0]  Total: 2.0388171\n",
      " \n"
     ]
    }
   ],
   "source": [
    "from qiskit_nature.algorithms.ground_state_solvers.minimum_eigensolver_factories import NumPyMinimumEigensolverFactory\n",
    "from qiskit_nature.algorithms.ground_state_solvers import GroundStateEigensolver\n",
    "import numpy as np \n",
    "\n",
    "def exact_diagonalizer(problem, converter):\n",
    "    solver = NumPyMinimumEigensolverFactory()\n",
    "    calc = GroundStateEigensolver(converter, solver)\n",
    "    result = calc.solve(problem)\n",
    "    return result\n",
    "\n",
    "result_exact = exact_diagonalizer(problem, converter)\n",
    "exact_energy = np.real(result_exact.eigenenergies[0])\n",
    "print(\"Exact electronic energy\", exact_energy)\n",
    "print(result_exact)\n",
    "\n",
    "# The targeted electronic energy for H2 is -1.85336 Ha\n",
    "# Check with your VQE result."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "67217a0c",
   "metadata": {},
   "source": [
    "#### 9. VQE and initial parameters for the ansatz\n",
    "Now we can import the VQE class and run the algorithm."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 96,
   "id": "53c335b9",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{   'aux_operator_eigenvalues': None,\n",
      "    'cost_function_evals': 28,\n",
      "    'eigenstate': array([-3.58849955e-17-2.04424520e-17j,  1.63824335e-16-1.49692699e-16j,\n",
      "       -8.04215005e-17-6.89851534e-17j, -6.31589224e-17-3.92523115e-17j,\n",
      "       -3.14509145e-17-1.70155183e-16j,  9.99478510e-01-1.09141895e-15j,\n",
      "       -2.09121074e-07+2.38109037e-18j,  1.92408556e-16-1.42186385e-16j,\n",
      "       -4.45951815e-17+2.40940530e-17j, -2.32842340e-07-4.56597600e-17j,\n",
      "       -3.22910007e-02-1.36943951e-17j, -2.93148272e-17-2.04354544e-17j,\n",
      "       -2.53745237e-17+4.54981921e-17j, -7.51784981e-17-9.56591292e-17j,\n",
      "        9.42942939e-19-4.69779457e-17j, -1.65723353e-17-5.55111512e-17j]),\n",
      "    'eigenvalue': (-1.697007028852752+0j),\n",
      "    'optimal_parameters': {   ParameterVectorElement(t[0]): 2.091210740916418e-07,\n",
      "                              ParameterVectorElement(t[1]): 2.3284234040085874e-07,\n",
      "                              ParameterVectorElement(t[2]): 0.03229661498347012},\n",
      "    'optimal_point': array([2.09121074e-07, 2.32842340e-07, 3.22966150e-02]),\n",
      "    'optimal_value': -1.697007028852752,\n",
      "    'optimizer_evals': None,\n",
      "    'optimizer_time': 0.8725237846374512}\n"
     ]
    }
   ],
   "source": [
    "from qiskit.algorithms import VQE\n",
    "from IPython.display import display, clear_output\n",
    "\n",
    "# Print and save the data in lists\n",
    "def callback(eval_count, parameters, mean, std):  \n",
    "    # Overwrites the same line when printing\n",
    "    display(\"Evaluation: {}, Energy: {}, Std: {}\".format(eval_count, mean, std))\n",
    "    clear_output(wait=True)\n",
    "    counts.append(eval_count)\n",
    "    values.append(mean)\n",
    "    params.append(parameters)\n",
    "    deviation.append(std)\n",
    "\n",
    "counts = []\n",
    "values = []\n",
    "params = []\n",
    "deviation = []\n",
    "\n",
    "# Set initial parameters of the ansatz\n",
    "# We choose a fixed small displacement \n",
    "# So all participants start from similar starting point\n",
    "try:\n",
    "    initial_point = [0.01] * len(ansatz.ordered_parameters)\n",
    "except:\n",
    "    initial_point = [0.01] * ansatz.num_parameters\n",
    "\n",
    "algorithm = VQE(ansatz,\n",
    "                optimizer=optimizer,\n",
    "                quantum_instance=backend,\n",
    "                callback=callback,\n",
    "                initial_point=initial_point)\n",
    "\n",
    "result = algorithm.compute_minimum_eigenvalue(qubit_op)\n",
    "\n",
    "print(result)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1014e339",
   "metadata": {},
   "source": [
    "The goal is to reach the chemical accuracy which is $\\delta E_{chem} = 0.001$ Ha $= 1$ mHa as efficiently as possible. Efficiency means that we want a low number of: \n",
    "* parameters\n",
    "* optimization steps\n",
    "* number of CNOT gates\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 97,
   "id": "43f41da4",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>optimizer</th>\n",
       "      <th>ansatz</th>\n",
       "      <th>rotation blocks</th>\n",
       "      <th>entanglement_blocks</th>\n",
       "      <th>entanglement</th>\n",
       "      <th>repetitions</th>\n",
       "      <th># of qubits</th>\n",
       "      <th># of parameters</th>\n",
       "      <th># of CNOTs</th>\n",
       "      <th># of evaluations</th>\n",
       "      <th>error (mHa)</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>L_BFGS_B</td>\n",
       "      <td>UCCSD</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>None</td>\n",
       "      <td>4</td>\n",
       "      <td>3</td>\n",
       "      <td>56</td>\n",
       "      <td>28</td>\n",
       "      <td>1.429967e-10</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "  optimizer ansatz rotation blocks entanglement_blocks entanglement  \\\n",
       "0  L_BFGS_B  UCCSD            None                None         None   \n",
       "\n",
       "  repetitions  # of qubits  # of parameters  # of CNOTs  # of evaluations  \\\n",
       "0        None            4                3          56                28   \n",
       "\n",
       "    error (mHa)  \n",
       "0  1.429967e-10  "
      ]
     },
     "execution_count": 97,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAElCAYAAAD6NKUrAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAx5ElEQVR4nO3deZhcVZ3/8fe39yTVWSCxIYCEHcJOIpiJS8KqiCKZAWR0RAFxGR1hFHVEEH+OPDg6g+OCCg7gAgQEIqsEwTQBlSWBsCRhF0hIyAIk3dWdXqrr+/vj3qpUuqurq5PurntvPq/nqSdV95y695yqzv3WOefec8zdERER6U9VpQsgIiLRpkAhIiIlKVCIiEhJChQiIlKSAoWIiJSkQCEiIiUpUIhsAzP7lJk9VIHjvtfMnhvp48r2SYFiO2Vmr5jZsWXmnWVmWTNLh4/Xzew7vfK4mbUV5NlQkNZoZv8THrPNzF4zs5vN7KiCPCeb2RIzazGz9Wb2ZzPbY4ByXWJm3QXHXG5m/zjAe6aEZa0pp+4jxcz+w8z+2GvbC/1s+5i7P+ju+41sKfNlmBV+hvN6bT803N5ciXLJ8FGgkHKtcveUu6eA9wBnm9lHe+U5NJfH3ccDmFk98GfgYOAkYCxwADAX+GCYZ2/gN8BXgHHAHsDPgJ4yynVjQbnOA35nZk3bUtH+DHNwWQj8g5lVh8faGagFDu+1be8wb0UUfAbrgBlmtmNB8pnA8yNfqtKi9qMgjhQoZNDc/e/AX4GpZWT/F2BX4KPu/oy797h7m7vf7O6XhHkOA/7u7vd7oNXdb3H31wZZrvlAK7BXue8xs3Fm9hszW2dmr5rZt8ysKkz7lJn9xcwuN7M3gUvMbEczuz1s+Tza+1hm9r9mtiJMX2xm7y1Iu8TMbgqP12pmS81sepj8GEFgOCx8/V5gAfBcr20vufuq8Ff9yoJ9H2FmT4T7/b2Z3Whm/xmmzTKzlWb2FTNba2arzezTBe+tN7Mfhi29NWb2CzMb1eu9XzezN4Brwrd1AX8APhbmqwZOB64b5Odxc1jWVjN73MwOLUh/JWxpLTOzt83sGjNrKEg/KWyFbjCzv5rZIb3e+3UzewpoU7DYNgoUMmhmtg8wE3i4jOzHAvPdva1EnseB/cMT8mwzS21FmczMPgTUAcsG8dafELRi9gTeD3wS+HRB+lHAy0AT8D2Clk4HsDNwVvgo9BjBiX0H4Hrg94UnN+AjBK2p8cDtwE8B3L0LeAR4X5jvfcCDwEO9tvVpTZhZHTAPuDY87g3AKb2y7RTWcxfgbOBnZjYhTLsM2Dcs995hnot7vXcHYHfg3ILtvyH4vABOAJ4BVg3y8zgZ+H1B+h/MrLYg/ePhvvcKy/itsM6HA1cDnwV2BH4J3B62YHPOAD4EjHf3DLL13F2P7fABvAIcW2beWUAW2AC0AA7cCtQV5PEwbUP4+HG4/T7gsoJ8hxXs57mC7e8GbiLo0uggOOmlBijXJQS/bDcAbQRdVV8b4D1TwrLWANXh+6cWpH8WaA6ffwp4rSCtGugG9i/YdinwUInjvU3QJZcr730FaVOBTb3qMy98/iSwD/CBXtvOLPhOVobP3we8DljBvh4C/rMg7yagpiB9bfiZW/jZ7VWQNoOghZd7bxfQ0OvvIXfsF4D9CILfx4Fzcp9fmZ/HwwVpVcBq4L0Ff6OfK0g/kaBFBfBz4Lu99v0c8P6C955V6f9nSXmoRSHlWuXu4919LMGv4U3Ar3vlOSLMM97d/y3c9ibBr28A3H2JB+MXc4D6gu0Pu/tp7j6JoIvlfcCFZZTrpvB4Ywh+dX7SzD4LYGZ/tM0D3R8v8t6JBN09rxZse5XgF3XOioLnkwgCzIpe+fPM7KsWDKpvtGBAf1x4nJw3Cp63Aw0F3SILgfeY2Q7AJHd/gaCL7x/CbQdRfHxiMvC6h2fIIuUGeNO3/FXdDqTCOo0GFoddOBuAe8LtOevcvaPIcQF+C3wRmE3QqtlCGZ9HvpzungVWhvUpVo9XC9J2B76SK3O4791KvFe2gQKFDJq7byToJvhwGdnvB443szGD2P9jBC2WgwZZrleAP+bK5e4f9M2D69cVect6ghbC7gXb3knw6zy/24Ln64AMwQmpMD8QXLIKfA04DZgQBsSNBL/ay/E3ghPpZ4C/hHVoIejO+QxBsP57kfetBnYxs8Lj7FYkXzHrCYL+gQVBfpwHFwfklJpi+rfAF4C73b29MKHMz2O3gvxVBONZq4qlE3zWubQVwPcKyjze3Ue7+w1lllsGQYFi+1ZrZg0Fj7IG/MIxhI8BS8vI/huCE9k8MzvIzKrDPurcIC5m9h4z+4yZvSN8vT9BX345YyCF5dqVoKumnHLh7j0E3V3fs+AS3t2Bfwd+VyL/rQSD2qPNbCrBlT45jQSBZB1QY2YXE1zlVRZ33wQsCsvwYEHSQ+G2/q52+htBt9sXzazGzE4GjizzmFngKuDygs9/FzM7ocz3/51gbKdY66+cz2Oamc0J//bOAzrZ8nv/VzPbNWxRXQjcGG6/CvicmR0Vjk+NMbMPmVljOeWWwVGg2L7dTfBrMve4pETeybluHIIugB0I+qRLCrssZhMMMN9FODYBvIvglyYEYwwfAZ4O938PQTfGf5VRh9MLyvUYwS/x7wzwHtj8a/NLBH30LxOckK8nGCTtzxcJumzeIBhHuaYgbX5Y9ucJPqMOBt/98QDwjrAsOQ+G24oGCg8GwucQDFJvAD4B3Elw0i3H14EXgYfNrIVgXKnsezTc/SF37z2IDeV9HrcRXC31NsEVcnPcvbsg/XrgXoLv5yXgP8NjLiJoZf00fO+LBGNKMgxsy25NkWQLL6FcGHaDJJaZPQL8wt2vGTBzhZjZJcDe7v6JftJfAc5x9/tGslzSl1oUst0I+8BPI+jeSRQze7+Z7RR2PZ0JHELwa15kmylQCABm9s2CK4QKH38c+N3DWq4/9lOub27F7l4jGMM4f4iLGQX7EVw+u4HgDvd/cvfVFS2RJIa6nkREpCS1KEREpCQFChERKUmBQqQEC6bN3rvXtkvM7HcFr8ea2Y8smFQvbWYvha8nFuT5ZzNbFKavDsde3hOmjTezq83sDQsmx3vezL7Rqwy5KdzfNLP7zez0kai/CChQiGwTCybkux84kGCgfCzBXElvEt70Zmb/DvyIYF6oJoI7jK8gmBAP4HKCezMOILgz+yME9wUUOjS8W3o/gvs3fmpm3x6maolsQYPZIiWYmQP7uPuLBdsuIbz+38zOIZhVdi93Txd5/ziCKUE+7e6/7+cYzwDfcvc/DKIM/0RwB/ku7v7m1tZPpBxqUYhsm2OBe4oFidAMoIEiE+YVeJhgGpFPWzCFezluI5igsKypOkS2hQKFyLbZkWAuq1Lp6730eghfIljw54vAMjN70cw+WOqg4TQX6wmmUhEZVgoUIqX1EExFXqiWYNZZ6DWNehFvAhNLTbjo7pvc/VJ3n0YQWG4iWOCn3yBgweI+k4C3Bq6CyLZRoBAp7TWCxY4K7cHmdSjuA04oMY363wgm5/toOQcLpxW/FBgTHqc/JxPMzPpoOfsV2RYKFCKl3Qh8K5zqusrMjiVY7+LmMP23BDOi3mJm+4d5dgynRDkxXLvjYoKlRz8aTk9ea2YfNLP/AjCzi8zsXWZWF07B/mWCqTie610YM9vBgkWYfgZ8XwPZMhK04LhIaf8vfDwETCCY6vrj7v4MgLt3hsHjO8CfwjxrCAabHwnz/LeZvUGw3vN1QCuwmOBqKQimPL+G4LLZDPAU8KFeA+RPhlc/dRHM6XS+u18/XJUWKaTLY0VEpCR1PYmISEkKFCIiUpIChYiIlKRAISIiJSXyqqeJEyf6lClT8q/b2toYM6a/y9zjT/WLv6TXUfWLvsWLF69390nF0hIZKKZMmcKiRZuXRW5ubmbWrFmVK9AwU/3iL+l1VP2iz8xe7S9NXU8iIlKSAoWIiJSkQCEiIiUpUIiISEkKFCIiUpIChYiIlKRAISIiJSXyPortyZ+fXcMfXuji8a4+SxeMqOqqKs44cjfeMbahouUQkaGnQBFzF/1hKa9v6MZefrFiZcjNVF9fW8Xn3r9XxcohIsNDgSLmWjZ1c9zuNVz1+RMqVgZ3Z58L/0jLpu6BM4tI7GiMIsayWSfdlWFUjVW0HGZGqqGGdGemouUQkeGhQBFjbV0Z3Kl4oABobKihtUOBQiSJFChiLHdiHhWBDsRUfS2tHep6EkmiigQKMzvVzJaaWdbMppfI9wEze87MXjSzb4xkGeMg19UzWi0KERlGlWpRPAPMARb2l8HMqoGfAR8EpgJnmNnUkSlePOR+wUehRTFWgUIksSpyinH35RAMgpZwJPCiu78c5p0LnAwsG/YCxkRLvuup8i2KVL0Gs0WSKgK/Rfu1C7Ci4PVK4Kj+MpvZucC5AE1NTTQ3N+fT0un0Fq+T4tHVwYnZuzdVvH4tb3XyVmtmWMqR1O+vUNLrqPrF27AFCjO7D9ipSNKF7n7bUB/P3a8ErgSYPn26F642lYTVp4pZ9chr8OTTTBw7uuL1e6TjWRaufJn3v//9A7UUBy2p31+hpNdR9Yu3YQsU7n7sNu7idWC3gte7htsklO4MxigaItD11NhQQybrdGayNNRWV7o4IjKEonx57GPAPma2h5nVAR8Dbq9wmSKltSODGTRE4Lzc2FALQIsukRVJnEpdHnuKma0EZgB3mdn8cPtkM7sbwN0zwBeB+cBy4CZ3X1qJ8kZVa0eGVH3NkHf1bI3G+qBxmtaVTyKJU6mrnuYB84psXwWcWPD6buDuESxarLR2ZBgb/pKvtMaG4E9Jl8iKJE+Uu55kAK0d3fkTdKWl6hUoRJJKgSLG0p2Z/Am60nJjFLkBdhFJDgWKGGvtyESmRZErR4taFCKJo0ARY60d3aQiNkahwWyR5FGgiLF0Z3RaFBqjEEkuBYoYa4lQ11NNdRWjaqs11bhIAilQxFRnpoeuTDZ//0IUNGqVO5FEUqCIqdxYQGNExihAa1LIyKmuruawww7LPy677LJKF6mkn/70p+y9996YGevXr+8332uvvcbxxx/PAQccwNSpU3nllVcA+POf/8wRRxzBQQcdxJlnnkkmE/w/+8EPfpD/DA466CCqq6t56623hrz80fk5KoPSmg8UNdBV4cKEUg21tKpFISNg1KhRLFmyZEj3mclkqKkZnlPizJkzOemkkwacOPCTn/wkF154IccddxzpdJqqqiqy2Sxnnnkm999/P/vuuy8XX3wxv/71rzn77LO54IILuOCCCwC44447uPzyy9lhhx2GvPxqUcRUrosnKvdRQG7xIo1RSOVMmTKFb3/72xxxxBEcfPDBPPvsswC0tbVx1llnceSRR3L44Ydz223BBNbXXnstH/nIRzj66KM55phjaG9v57TTTmPq1KmccsopHHXUUSxatIirr76a8847L3+cq666ivPPP7/sch1++OFMmTKlZJ5ly5aRyWQ47rjjAEilUowePZo333yTuro69t13XwCOO+44brnllj7vv+GGGzjjjDPKLtNgKFDEVG7yvSh1PaXq1fUkI2PTpk1bdD3deOON+bSJEyfy+OOP8/nPf54f/vCHAHzve9/j6KOP5tFHH2XBggVccMEFtLW1AfD4449z880388ADD3DFFVcwYcIEli1bxne/+10WL14MwGmnncYdd9xBd3fw/+6aa67hrLPO4uKLL+b224dmrtLnn3+e8ePHM2fOHA4//HAuuOACenp6mDhxIplMhkWLFgFw8803s2LFii3e297ezj333MM//uM/DklZeovOz1EZlMKup84KlyWnsaFG91HIiCjV9TRnzhwApk2bxq233grAvffey+23354PHB0dHbz22mtA8As9113z0EMP8eUvfxmAgw46iEMOOQQIft0fffTR3HnnnRxwwAF0d3dz8MEHc/DBBwMMyaJFmUyGBx98kCeeeIJ3vvOdnH766Vx77bWcffbZzJ07l/PPP5/Ozk6OP/54qqu3nDL6jjvuYObMmcPS7QQKFLGVLggU/Q+NjazGhlp1PUnF1dfXA8GAd27Q19255ZZb2G+//bbI+8gjjzBmzJiy9nvOOedw6aWXsv/++/PpT3+6ZN4TTjiBNWvWMH36dH71q1+Vtf9dd92Vww47jD333BOAj370ozz88MOcffbZzJgxgwcffBAIgt7zzz+/xXvnzp07bN1OoK6n2GqNaNdTW1cPPVmvdFFEtnDCCSfwk5/8BPfgb/OJJ54omm/mzJncdNNNQDBm8PTTT+fTjjrqKFasWMH1118/4El5/vz5LFmypOwgAfCud72LDRs2sG7dOiC40mnq1KkArF27FoDOzk6+//3v87nPfS7/vo0bN/LAAw9w8sknl32swVKgiKkoDmbnp/HQlU8yzHqPUXzjG98omf+iiy6iu7ubQw45hAMPPJCLLrqoaL4vfOELrFu3jqlTp/Ktb32LAw88kHHjxuXTTzvtNGbOnMmECRMAyh6j+PGPf8yuu+7KypUrOeSQQzjnnHMAWLRoUf55dXU1P/zhDznmmGM4+OCDcXc+85nPAMFlsAcccACHHHIIH/7whzn66KPz+543bx7HH3982S2jreLuiXtMmzbNCy1YsMCT5tK7lvm+F97t7tGp39xHX/Xdv36nr3irbUj3G5X6Daek1zEu9ctkMr5p0yZ3d3/xxRd9ypQp3tnZmU//0Ic+5Pfdd1+f98WlfqUAi7yfc2p0fo7KoERp+o6czVONq0Uh8dTe3s7s2bPp7u7G3bniiiuoq6tjw4YNHHnkkRx66KEcc8wxlS7miIvWmUbKFkwIGJ3xCdAqdxJ/jY2N+ctQC40fP77PAPL2pFJrZp9qZkvNLGtm00vku9rM1prZMyNZvjiI0up2OSmtmy2SSJUazH4GmAMsHCDftcAHhr00MdTaEZ3V7XJyLZwWXSIrkigVCRTuvtzdnysj30Jg6Ge4SoB0JMco1PUkkkS6PDamgq6naI5RaDBbJFmG7Sepmd0H7FQk6UJ3v20YjncucC5AU1PTFrfUp9PpIbnFPkrebutg4/o3aG5+OzL1c3eqDJ557iWafcXAbyhTVOo3nJJeR9Uv3oYtULj7scO1736OdyVwJcD06dO9cDrf5ubmAaf3jZNs1umYfzcH7DWFWbP2i1T9Ghfey45Nk5k166Ah22eU6jdckl5H1S/e1PUUQ21dGdwhFbExCtDiRSJJVKnLY08xs5XADOAuM5sfbp9sZncX5LsB+Buwn5mtNLOzK1HeqMmNAURtjAKCS2RbFChEEqUiP0ndfR4wr8j2VcCJBa+HbzrEGNtidbuIGdtQS7pTl8eKJIm6nmIoFyiidh8FqOtJJIkUKGIoilOM56QaanR5rEjCKFDEUJS7ntSiEEkeBYoY2jyYHcVAEaxy567Fi0SSQoEihiLd9VRfQ3eP05nJVrooIjJEFChiqLUjgxmMrq0eOPMIG6v5nkQSR4EihnIzx1ZVWaWL0kdK8z2JJI4CRQy1dmQYG8FuJ4DG+qBcrZpqXCQxFChiKN3ZHcl7KEBTjYskkQJFDLVGcC2KnJQChUjiKFDEUGtHJpITAgL5LjF1PYkkhwJFDKU7M5G8NBY2TyuiFoVIcihQxFCwul00WxS66kkkeRQoYqi1I0NjRAeza6urGFVbra4nkQRRoIiZrkyWzkw2si0K0MSAIkmjQBEzuV/qUb08FoJLZLV4kUhyKFDETJRXt8tprNcMsiJJokARM1GeYjynsaGWtMYoRBJDgSJmWnJdT5EOFGpRiCRJRQKFmZ1qZkvNLGtm0/vJs5uZLTCzZWHeL490OaMoHZ6AozrXEwTjJxrMFkmOSrUongHmAAtL5MkAX3H3qcC7gX81s6kjUbgoi0vXk1oUIslRkbONuy8HMOt/mmx3Xw2sDp+3mtlyYBdg2UiUMapyv9SjfNVT7vLYnqxTHcGp0EVkcKJ7tilgZlOAw4FHSuQ5FzgXoKmpiebm5nxaOp3e4nWcPflSFwCPP/pXasOTcNTqt3ZlMI4y//5mRtdue6CIWv2GQ9LrqPrF27AFCjO7D9ipSNKF7n7bIPaTAm4BznP3lv7yufuVwJUA06dP91mzZuXTmpubKXwdZ39rX07d31/huKNn57dFrX5rxrzG3Oee5tB3vZtdxo/a5v1FrX7DIel1VP3ibdgChbsfu637MLNagiBxnbvfuu2lir/Wzkx+udGoSoWLF6U1TiGSCJG9PNaCAYz/A5a7+/9UujxREaxFEd0rnqBw8SLdSyGSBJW6PPYUM1sJzADuMrP54fbJZnZ3mG0m8C/A0Wa2JHycWInyRkm6I7qr2+Vo8SKRZKnUVU/zgHlFtq8CTgyfPwTokpleory6XU6ua6xV91KIJEJku56kuNaOTORbFI1a5U4kURQoYibKq9vl5AKZBrNFkkGBImZaIry6Xc7oumqqTGMUIkmhQBEj2ayHLYpoBwozI1Vfo64nkYRQoIiR9u4e3KM9z1NOY0OtBrNFEkKBIkZyv9CjPkYBmmpcJEkUKGIkNzgc9aueIAgUGswWSQYFihhpicEU4zlB15PGKESSQIEiRjZ3PUU/UKS0brZIYihQxEhuLYq4jFGo60kkGRQoYiQOq9vlaJU7keRQoIiRuA1md/Vk6cz0VLooIrKNFChipLWjGzMYUxePQAG6O1skCRQoYqSlI0OqroaqGKxDnWv1KFCIxJ8CRYzEYfqOnNyAuwa0ReKvrEBhZv9tZgcOd2GktNaO7lhc8QRa5U4kScptUSwHrjSzR8zsc2Y2bjgLJcW1dmTyq8dFXb7rSfM9icReWYHC3X/l7jOBTwJTgKfM7Hozmz2chZMtxanraWx+8SIFCpG4K3uMwsyqgf3Dx3rgSeDfzWzuYA9qZqea2VIzy5rZ9H7yNJjZo2b2ZJj3O4M9TtIEy6DGo+sppa4nkcQo6+epmV0OfBi4H7jU3R8Nk75vZs9txXGfAeYAvyyRpxM42t3TZlYLPGRmf3T3h7fieIkQh2VQc3ItHw1mi8RfuWedp4BvuXtbkbQjB3tQd18OwQI3JfI4kA5f1oYPH+yxkqS1o5uxMel6qq2uoqG2SmMUIglQ7lnnSWC/Xif2jcCr7r5xyEsVCru7FgN7Az9z90dK5D0XOBegqamJ5ubmfFo6nd7idRxlsk5nJsu61Stobl6zRVpU61dnzvMvv9anvIMV1foNpaTXUfWLOXcf8AE8DHQBiwhO3J3A48BLwPH9vOc+gi6m3o+TC/I0A9PLOP54YAFwUDnlnTZtmhdasGCBx92b6U7f/et3+jUPvdwnLar1m/2DBf6v1y3e5v1EtX5DKel1VP2iD1jk/ZxTy21RrALOdvelAGY2Ffh/wNeAW4F7iwSgYwcVsUpw9w1mtgD4AEGw2e7EaXW7nJRWuRNJhHKveto3FyQA3H0ZsL+7vzw8xQIzm2Rm48Pno4DjgGeH63hRlzvhxuU+CginGtcYhUjslRsolpnZz83s/eHjinBbPTDo6x/N7BQzWwnMAO4ys/nh9slmdneYbWdggZk9BTwG/Mnd7xzssZIiTlOM5zTW1+ryWJEEKPescybwBeC88PVfgK8SBIlB33Tn7vOAeUW2rwJODJ8/BRw+2H0nVb7rqT5eXU+6PFYk/gYMFOGVR3e7+2zgv4tkSRfZJkNs8+p2MWpRaIxCJBEG7Hpy9x4gq/mdKiueXU81pLsyZLPb9e0vIrFX7lknDTxtZn8C8jfdufu/DUuppI9c11O8BrNrcYe2rvhMPSIifZV71rk1fEiFtHZmqKupor6mutJFKVvhKncKFCLxVVagcPdfh5eovtPdt2ZuJ9lGrR2Z2EzfkZNr/egSWZF4K3fhog8DS4B7wteHmdntw1gu6SUdowkBcxrzU43rElmROCv3PopLCCb/2wDg7kuAPYelRFJUnFa3y8kFthZd+SQSa+UGim7vO/lfdqgLI/2L0xTjOWM11bhIIpQbKJaa2T8D1Wa2j5n9BPjrMJZLeonT6nY5jVrlTiQRyg0UXwIOJJg19gaghc13acsIiOOVQ5sHszVGIRJn5V711A5cGD6kAoIxini1KMbUVWOmFoVI3JW7FOq+BHM7TSl8j7sfPTzFkkLuHsuuJzMjVa9pPETirtwzz++BXwC/AnqGrzhSTFtXD1kndoPZAGMbahUoRGKu3DNPxt1/PqwlkX6l8/M8xWuMAnITA2qMQiTOyh3MvsPMvmBmO5vZDrnHsJZM8javbhe/FkWqXosXicTdYNajALigYJujm+5GREsMV7fLaWyoYX26q9LFEJFtUO5VT3sMd0Gkf7lf5HGb6wkg1VDL39e3DZxRRCKrZNeTmX2t4PmpvdIuHa5CyZY2dz3Fc4xCXU8i8TbQGMXHCp7/R6+0D2ztQc3sVDNbamZZM5s+QN5qM3vCzLbb9bJzg9lxvOqpsaFGcz2JxNxAgcL6eV7s9WA8A8wBFpaR98vA8m04VuzFcXW7nMb6GroyWTozuqpaJK4GChTez/Nir8vm7svLWdfCzHYFPkRw/8Z2q7WjGzMYUxfDQBF2l2liQJH4GujMc6iZtRC0HkaFzwlfNwxryQI/Ar4GNA6U0czOBc4FaGpqorm5OZ+WTqe3eB03y1/qpKEaFi58oGh6lOu38vVgfOW+B/5C05hyr8beUpTrN1SSXkfVL95KBgp33+p1N83sPmCnIkkXuvttZbz/JGCtuy82s1kD5Xf3K4ErAaZPn+6zZm1+S3NzM4Wv4+bOdU8yYcP6fusQ5fp1LX2Dq55ezIGHTeOgXcZt1T6iXL+hkvQ6qn7xNmx9Ge5+7DbuYibwETM7kaD1MtbMfufun9j20sVLuiMTy3soYHPXU4vuzhaJra3rCxgB7v4f7r6ru08huPrqz9tjkABo7Yzf6nY5jVq8SCT2KhIozOwUM1sJzADuMrP54fbJZnZ3JcoUZXFc3S4nFyg0MaBIfFXk7OPu84B5RbavAk4ssr0ZaB72gkVUuiPDO3cYXelibJVcgNPEgCLxFdmuJ9msJYar2+XkL4/V3dkisaVAEQNxXN0up66mivqaKnU9icSYAkXEBXc1Z2mM6RgFhGtSqEUhElsKFBGX67KJa4sCgu4ntShE4kuBIuLyEwLGdIwCCNfN1mC2SFwpUERcS4xXt8tpbKjRfRQiMaZAEXH5mWPjPkahQCESWwoUEbd5jCLOXU+1ujxWJMYUKCKuNSFdT5rrSSS+FCgiLvdLPK6TAsLm5VCz2a1ewkREKkiBIuLivLpdTmNDDe7Q3q1V7kTiSIEi4lo6uqmrrqK+ZquXBqm43PiKLpEViScFiohLd2Ri3ZqAzRMD6hJZkXhSoIi41gQEilz5WxQoRGJJgSLiWju6Yz2QDep6Eok7BYqIS3dmaKyP7z0UULDKne6lEIklBYqIS1LXk+7OFoknBYqIa+3IxL7rSYPZIvFWqTWzTzWzpWaWNbPpJfK9YmZPm9kSM1s0kmWMitaObsbGePoOgDF1NZhpjEIkrir1U/UZYA7wyzLyznb39cNcnkhyd9Kdmfwv8riqqjJS9TW66kkkpipyBnL35QBmVonDx0Z7Vw9Zj/dd2TmN9TUazBaJqaifgRy418wc+KW7X9lfRjM7FzgXoKmpiebm5nxaOp3e4nVcvN2RBWD1ay/T3Lyi33xxqF9VTyd/X7ma5ua3B/3eONRvWyW9jqpfvA1boDCz+4CdiiRd6O63lbmb97j762b2DuBPZvasuy8sljEMIlcCTJ8+3WfNmpVPa25upvB1XLy4thWaF3LEIQcy69DJ/eaLQ/2alv+V+poqZs1696DfG4f6bauk11H1i7dhCxTufuwQ7OP18N+1ZjYPOBIoGiiSqCUBEwLmNDbU8FZbV6WLISJbIbKXx5rZGDNrzD0HjicYBN9uJGF1u5zGhlrdRyESU5W6PPYUM1sJzADuMrP54fbJZnZ3mK0JeMjMngQeBe5y93sqUd5KSXfEf3W7nFS9lkMViatKXfU0D5hXZPsq4MTw+cvAoSNctEhJwup2OWMbanQfhUhMRbbrSTZ3PcX9zmwIWhSdmSxdmWyliyIig6RAEWGtnRnMIFUX/0ChiQFF4kuBIsJaO7pJ1dVQVRX/GxM11bhIfClQRFg6ARMC5qQ0g6xIbClQRFgSphjP0VTjIvGlQBFhrZ3dsZ8QMCe3+JK6nkTiR4EiwtIdmUTcQwEazBaJMwWKCFPXk4hEgQJFhLV2JidQpNSiEIktBYoIa+3oTkzXU31NNXU1VbRojEIkdhQoIqq7J0tHdzYxg9kQTG6orieR+FGgiKh0gqYYz2lsqMnXS0TiQ4EioloTNHNsTjDVuLqeROJGgSKicn35Sep6SmndbJFYUqCIqNwJdWzCup40RiESPwoUEZXErqeUAoVILClQRFS6M+x6SlCLYqzGKERiSYEiolqTetVTZwZ3r3RRRGQQKrVm9qlmttTMsmY2vUS+8WZ2s5k9a2bLzWzGSJazkvKr2yVsMDvr0N7VU+miiMggVKpF8QwwB1g4QL7/Be5x9/0J1s9ePtwFi4rWjgx11VU01FZXuihDZvPiRRqnEImTivxcdfflAGb9r9xmZuOA9wGfCt/TBXSNQPEiIZi+IzmtCShcvKibncY1VLg0IlKuKI9R7AGsA64xsyfM7FdmNqbShRop6c7krG6Xk59BVvdSiMSKDdfAopndB+xUJOlCd78tzNMMfNXdFxV5/3TgYWCmuz9iZv8LtLj7Rf0c71zgXICmpqZpc+fOzael02lSqdQ21mhkXb64gw2dznf+YdSAeeNSvxfe7uF7j3TwlWn1HDyp/CAYl/pti6TXUfWLvtmzZy929+Jjxu5esQfQDEzvJ20n4JWC1+8F7ipnv9OmTfNCCxYs8Lj5p5//xU//5V/LyhuX+j27usV3//qdfueTqwb1vrjUb1skvY6qX/QBi7yfc2pku57c/Q1ghZntF246BlhWwSKNqNYErW6X01gwRiEi8VGpy2NPMbOVwAzgLjObH26fbGZ3F2T9EnCdmT0FHAZcOuKFrZAkrW6Xk9IqdyKxVKmrnuYB84psXwWcWPB6CdDvfRZJ1trRTWOC7qEASNXVYKbBbJG4iWzX0/bM3Ul3Jq/rqarKSNXVqOtJJGYUKCKovauHrCdr+o6clBYvEokdBYoIyk0xnrT7KEBTjYvEkQJFBOW6ZpLW9QTBfE+tnep6EokTBYoIasnNHJuwwWwIgp+6nkTiRYEigtIJnGI8R11PIvGjQBFBSVzdLqexoUaXx4rEjAJFBCVxdbucRq1yJxI7ChQRlMTV7XJS9TV0dGfp7slWuigiUiYFigjKDWaPqUteoMgFPw1oi8SHAkUEpTsypOprqK7qf2GnuNIqdyLxo0ARQUlc3S4ntwa47qUQiQ8FighqDVsUSTRWM8iKxI4CRQQFEwImM1BoqnGR+FGgiKCg6yl591DA5jGKtLqeRGJDgSKCWjszibyHAgpXuVOLQiQuFCgiqLUjk+/LT5r8YLYChUhsKFBEUGtHd2IHsxtqq6mrrlKgEIkRBYqI6e7J0tGdTewYBQQD2prGQyQ+KvKz1cxOBS4BDgCOdPdFRfLsB9xYsGlP4GJ3/9FwlasrU/lpJTZs6gKSOX1HTmNDTX5xJhGJPnP3kT+o2QFAFvgl8NVigaJX/mrgdeAod391oP3vsPsBftw3r86/3rBhA5943wH8y4wpbOrq4VPXPNrnPf80bVcuvm0pm7p7Blmb4bHXpDFMTNXnX1/84akcOHkcD72wnp/8+YUt8m7YsIErznove01Kcd+yNVz14Mt99nf56Ycxefwo7nhyFb97uO9H+PNPTGOHMXX8ftEKbl68sk/6tZ8+klF11fz2b69w51Or+6Tf+NkZAFy58CXuX752i7SG2mp+fdaRAPz4/hf4xQMv0dHdQ11N0KCtMmOncQ0AvJXuoiPT6zvIZpm8wxgA1qc7+wT02uoqJjUGn9W61s4+80jV1VTlP8u1LR1kslv+zTfUVLNDqg6ANzZ2kO31f2JUbTU7jAnSVxdJH1Nfw4TRQfrKt9v7fDaN9bWMG11L1p1VGzb1SR/bUEtNtov6UaNZvbFv+rhRtTQ21NLdk2VNS0ef9Amj6xhTX0NXJsva1r7pO46pY3RdDZ3dPaxLd/ZNT9UzqraaTd09vFkkfVKqnvraatq7MrzV1tUn/R2NDdTVVJHuzLChvW9609gGujo3ka2qY+Omvi3JnceNorrKaNnUTUuRlubk8aOoMmNDe3fRq+V2nTAagLfbu2jr9QPEzNhl/CgA3mzrYlPXlulVVcbkcUH6+nQnHb3+/9dUVeX/Nte1dtLZ62+ztrqKprENtLW3kc5U9/nbq6+pzv9tvrGxg0x2y/SG2ur83+aqjZvI9vrbHFVXw47h397rGzbR+3zd+2+vpqqKqZPH5tNPOmTnAc97p07fjbfautgxVb/Y3af3yUSFWhTuvhyCL7FMxwAvlRMktsV5x+5DS0c3ty9Z1SftsN3Gs//OY2nZ1M3dT/c9UU6fsgN7vyPFW21d3Lv0jT7p795rR6bsOIY1LR0seHZtn/T37jOJXSaM4vW3N/HQi+vzX34STR7XwFvtm//D11VXMXXn4I97+eoW3m7f8mRgPd359KWrNrJx05b/2VP11fn0p3o2kO7c8j/zuFE1+fTO7h42dW/5n3XC6FoOCNPbO3vo6vWffWKqjn2aGsGhZVOGnl7/WSel6tlrUgqg+Il0bD1TdhxDT9bZ0N73RLfTuAbqMxnGT0gVPZFOHjeKyeNHsamrp+jYzi7jR9E0toF0Z6ZIS83ZZcIoJqbqadnUTbqr7/t3mzCKCWPqeLuti/Yi6e/ccTRjR9WyPt1Z9IfU7hNHk6qvYU1LR58TKcAeE8fQuqGTTF1Dn88WYM9JY6irqWLl2+10Z/um7/2OFNVVxqtvttOzsW/6Pk3BZ//yurY+QbzaLJ/OmjRvsmV6bXVVPr3HvU+ga6itzqd392TpHafH1NWwT1OKtWs3saarlrZen9/Yhtr8+9u6Mn0C0fjRdfn0lo7uPoFmxzGb0ze0d/X525uYqmfPScGPqDfbOqmuGp7RhIq0KPIHN2umvBbF1cDj7v7TEnnOBc4FaGpqmjZ37tx8WjqdJpVKDUmZo0j1i7+k11H1i77Zs2f326LA3YflAdwHPFPkcXJBnmZg+gD7qQPWA03lHnvatGleaMGCBZ5kql/8Jb2Oql/0AYu8n3PqsHU9ufuxQ7SrDxK0JtYM0f5ERGQQ4nB57BnADZUuhIjI9qoigcLMTjGzlcAM4C4zmx9un2xmdxfkGwMcB9xaiXKKiEjlrnqaB8wrsn0VcGLB6zZgxxEsmoiI9BKHricREakgBQoRESlJgUJEREqq6A13w8XM1gGFd3FPJLgXI6lUv/hLeh1Vv+jb3d0nFUtIZKDozcwWeX93HCaA6hd/Sa+j6hdv6noSEZGSFChERKSk7SVQXFnpAgwz1S/+kl5H1S/GtosxChER2XrbS4tCRES2kgKFiIiUlOhAYWYfMLPnzOxFM/tGpcszHMzsFTN72syWmFnJBaDiwMyuNrO1ZvZMwbYdzOxPZvZC+O+ESpZxW/RTv0vM7PXwO1xiZieW2keUmdluZrbAzJaZ2VIz+3K4PUnfYX91TMz32FtixyjCdbafJ5h9diXwGHCGuy+raMGGmJm9QrD4U9xv9gHAzN4HpIHfuPtB4bb/At5y98vCgD/B3b9eyXJurX7qdwmQdvcfVrJsQ8HMdgZ2dvfHzawRWAx8FPgUyfkO+6vjaSTke+wtyS2KI4EX3f1ld+8C5gInV7hMMgB3Xwi81WvzycCvw+e/JvhPGUv91C8x3H21uz8ePm8FlgO7kKzvsL86JlaSA8UuwIqC1ytJ5pfpwL1mtjhcNzyJmtx9dfj8DaCpkoUZJl80s6fCrqnYdssUMrMpwOHAIyT0O+xVR0jg9wjJDhTbi/e4+xEES8b+a9i1kVjh2r5J6y/9ObAXcBiwGvjvipZmCJhZCrgFOM/dWwrTkvIdFqlj4r7HnCQHiteB3Qpe7xpuSxR3fz38dy3BYlBHVrZEw2JN2C+c6x9eW+HyDCl3X+PuPe6eBa4i5t+hmdUSnECvc/fc6pSJ+g6L1TFp32OhJAeKx4B9zGwPM6sDPgbcXuEyDSkzGxMOpuWWjT0eeKb0u2LpduDM8PmZwG0VLMuQy51AQ6cQ4+/QzAz4P2C5u/9PQVJivsP+6pik77G3xF71BBBenvYjoBq42t2/V9kSDS0z25PNS8rWANfHvY5mdgMwi2Da5jXAt4E/ADcB7ySYPv40d4/lgHA/9ZtF0F3hwCvAZwv682PFzN4DPAg8DWTDzd8k6MNPynfYXx3PICHfY2+JDhQiIrLtktz1JCIiQ0CBQkRESlKgEBGRkhQoRESkJAUKEREpSYFCpBczS4f/TjGzfx7ifX+z1+u/DuX+RYaDAoVI/6YAgwoUZlYzQJYtAoW7/8MgyyQy4hQoRPp3GfDecG2B882s2sx+YGaPhRO/fRbAzGaZ2YNmdjuwLNz2h3CixqW5yRrN7DJgVLi/68JtudaLhft+Jlxf5PSCfTeb2c1m9qyZXRfeGYyZXRauifCUmSVuamuJjoF+/Yhsz74BfNXdTwIIT/gb3f1dZlYP/MXM7g3zHgEc5O5/D1+f5e5vmdko4DEzu8Xdv2FmX3T3w4ocaw7BXb2HEty1/ZiZLQzTDgcOBFYBfwFmmtlygmki9nd3N7PxQ1t1kc3UohAp3/HAJ81sCcGUFDsC+4RpjxYECYB/M7MngYcJJqfch9LeA9wQTiq3BngAeFfBvleGk80tIegS2wh0AP9nZnOA9m2sm0i/FChEymfAl9z9sPCxh7vnWhRt+Uxms4BjgRnufijwBNCwDcftLHjeA9S4e4ZgdtKbgZOAe7Zh/yIlKVCI9K8VaCx4PR/4fDjFNGa2bzhrb2/jgLfdvd3M9gfeXZDWnXt/Lw8Cp4fjIJOA9wGP9lewcC2Ece5+N3A+QZeVyLDQGIVI/54CesIupGuB/yXo9nk8HFBeR/ElPe8BPheOIzxH0P2UcyXwlJk97u4fL9g+D5gBPEkw++jX3P2NMNAU0wjcZmYNBC2df9+qGoqUQbPHiohISep6EhGRkhQoRESkJAUKEREpSYFCRERKUqAQEZGSFChERKQkBQoRESnp/wObEqOvolDmRgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Store results in a dictionary\n",
    "from qiskit.transpiler import PassManager\n",
    "from qiskit.transpiler.passes import Unroller\n",
    "\n",
    "# Unroller transpile your circuit into CNOTs and U gates\n",
    "pass_ = Unroller(['u', 'cx'])\n",
    "pm = PassManager(pass_)\n",
    "ansatz_tp = pm.run(ansatz)\n",
    "cnots = ansatz_tp.count_ops()['cx']\n",
    "\n",
    "energy = result.optimal_value\n",
    "\n",
    "if ansatz_type == \"TwoLocal\":\n",
    "    result_dict = {\n",
    "        'optimizer': optimizer.__class__.__name__,\n",
    "        'mapping': converter.mapper.__class__.__name__,\n",
    "        'ansatz': ansatz.__class__.__name__,\n",
    "        'rotation blocks': rotation_blocks,\n",
    "        'entanglement_blocks': entanglement_blocks,\n",
    "        'entanglement': entanglement,\n",
    "        'repetitions': repetitions,\n",
    "        'skip_final_rotation_layer': skip_final_rotation_layer,\n",
    "        'energy (Ha)': energy,\n",
    "        'error (mHa)': (energy-exact_energy)*1000,\n",
    "        '# of parameters': len(result.optimal_point),\n",
    "        'final parameters': result.optimal_point,\n",
    "        '# of evaluations': result.cost_function_evals,\n",
    "        '# of qubits': int(qubit_op.num_qubits),\n",
    "        '# of CNOTs': cnots}\n",
    "else:\n",
    "    result_dict = {\n",
    "        'optimizer': optimizer.__class__.__name__,\n",
    "        'mapping': converter.mapper.__class__.__name__,\n",
    "        'ansatz': ansatz.__class__.__name__,\n",
    "        'rotation blocks': None,\n",
    "        'entanglement_blocks': None,\n",
    "        'entanglement': None,\n",
    "        'repetitions': None,\n",
    "        'skip_final_rotation_layer': None,\n",
    "        'energy (Ha)': energy,\n",
    "        'error (mHa)': (energy-exact_energy)*1000,\n",
    "        '# of parameters': len(result.optimal_point),\n",
    "        'final parameters': result.optimal_point,\n",
    "        '# of evaluations': result.cost_function_evals,\n",
    "        '# of qubits': int(qubit_op.num_qubits),\n",
    "        '# of CNOTs': cnots}\n",
    "\n",
    "# Plot the results\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "fig, ax = plt.subplots(1, 1)\n",
    "ax.set_xlabel('Iterations')\n",
    "ax.set_ylabel('Energy')\n",
    "ax.grid()\n",
    "fig.text(0.7, 0.75, f'Energy: {result.optimal_value:.3f}')\n",
    "plt.title(f\"{result_dict['optimizer']}-{result_dict['mapping']}\\n{result_dict['ansatz']}\")\n",
    "ax.plot(counts, values)\n",
    "ax.axhline(exact_energy, linestyle='--')\n",
    "fig_title = f\"\\\n",
    "{result_dict['optimizer']}-\\\n",
    "{result_dict['mapping']}-\\\n",
    "{result_dict['ansatz']}-\\\n",
    "Energy({result_dict['energy (Ha)']:.3f}).png\"\n",
    "fig.savefig(fig_title, dpi=300)\n",
    "\n",
    "# Display and save the data\n",
    "import pandas as pd\n",
    "import os.path\n",
    "filename = 'results_h2o.csv'\n",
    "if os.path.isfile(filename):\n",
    "    result_df = pd.read_csv(filename)\n",
    "    result_df = result_df.append([result_dict])\n",
    "else:\n",
    "    result_df = pd.DataFrame.from_dict([result_dict])\n",
    "result_df.to_csv(filename)\n",
    "result_df[['optimizer','ansatz', 'rotation blocks', 'entanglement_blocks',\n",
    "    'entanglement', 'repetitions', '# of qubits', '# of parameters', '# of CNOTs', '# of evaluations', 'error (mHa)']]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "67d1222a",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-block alert-danger\">\n",
    "    \n",
    "<b>Tutorial questions 3</b> \n",
    "\n",
    "Experiment with all the parameters and answer the following:\n",
    "\n",
    "1. Can you find an efficient heuristic ansatz? \n",
    "2. Which of the `UCCSD, PUCCD or SUCCD` ansatz is best? \n",
    "3. Which optimizer perform best with `UCCSD`? \n",
    "\n",
    "Remember, to find the most efficient ansatz you have to consider the accuracy, number of parameters and number of optimization steps, number of CNOTs. This is a discussion question and will not be graded. Feel free to discuss your answers in the Slack channel. \n",
    "</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e0e1ad7c",
   "metadata": {},
   "source": [
    "## Part 2: qEOM for water molecule \n",
    "\n",
    "\n",
    "In this part, you will compute the excited states of the previously defined water molecule\n",
    "\n",
    "To use the qEOM as a solver in Qiskit, we have to define a ground state calculation first. With this the qEOM solver can be initialized:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 98,
   "id": "46886db2",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "=== GROUND STATE ENERGY ===\n",
      " \n",
      "* Electronic ground state energy (Hartree): -84.700895930363\n",
      "  - computed part:      -1.69700702885\n",
      "  - ActiveSpaceTransformer extracted energy part: -83.003888901513\n",
      "~ Nuclear repulsion energy (Hartree): 9.819476509782\n",
      "> Total ground state energy (Hartree): -74.88141942058\n",
      " \n",
      "=== EXCITED STATE ENERGIES ===\n",
      " \n",
      "  1: \n",
      "* Electronic excited state energy (Hartree): -84.276708397288\n",
      "> Total excited state energy (Hartree): -74.457231887506\n",
      "  2: \n",
      "* Electronic excited state energy (Hartree): -84.199633757776\n",
      "> Total excited state energy (Hartree): -74.380157247994\n",
      "  3: \n",
      "* Electronic excited state energy (Hartree): -83.506845228421\n",
      "> Total excited state energy (Hartree): -73.687368718639\n",
      " \n",
      "=== MEASURED OBSERVABLES ===\n",
      " \n",
      "  0:  # Particles: 2.000 S: 0.000 S^2: 0.000 M: 0.000\n",
      " \n",
      "=== DIPOLE MOMENTS ===\n",
      " \n",
      "~ Nuclear dipole moment (a.u.): [2.86710004  0.0  0.0]\n",
      " \n",
      "  0: \n",
      "  * Electronic dipole moment (a.u.): [2.06496784  0.0  0.0]\n",
      "    - computed part:      [0.00228274  0.0  0.0]\n",
      "    - ActiveSpaceTransformer extracted energy part: [2.0626851  0.0  0.0]\n",
      "  > Dipole moment (a.u.): [0.8021322  0.0  0.0]  Total: 0.8021322\n",
      "                 (debye): [2.0388165  0.0  0.0]  Total: 2.0388165\n",
      " \n"
     ]
    }
   ],
   "source": [
    "from qiskit import Aer\n",
    "from qiskit.utils import QuantumInstance\n",
    "from qiskit_nature.algorithms import GroundStateEigensolver, QEOM, VQEUCCFactory\n",
    "\n",
    "# This first part sets the ground state solver\n",
    "# see more about this part in the ground state calculation tutorial\n",
    "quantum_instance = QuantumInstance(Aer.get_backend(\"aer_simulator_statevector\"))\n",
    "solver = VQEUCCFactory(quantum_instance)\n",
    "gsc = GroundStateEigensolver(converter, solver)\n",
    "\n",
    "# The qEOM algorithm is simply instantiated with the chosen ground state solver\n",
    "qeom_excited_states_calculation = QEOM(gsc, \"sd\")\n",
    "\n",
    "qeom_results = qeom_excited_states_calculation.solve(problem)\n",
    "\n",
    "print(qeom_results)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "137ac260",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-block alert-danger\">\n",
    "    \n",
    "<b>Tutorial questions 4</b> \n",
    "\n",
    "Obtain the exact excited states energies and compare with the ones obtained from qEOM. (Coding)\n",
    "\n",
    "Post your thoughts on the observed differences in the Slack channel.\n",
    "    \n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 99,
   "id": "149a9757",
   "metadata": {},
   "outputs": [],
   "source": [
    "# WRITE YOUR CODE BETWEEN THESE LINES - START\n",
    "\n",
    "\n",
    "\n",
    "# WRITE YOUR CODE BETWEEN THESE LINES - END"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f9cf1807",
   "metadata": {},
   "source": [
    "## Additional information\n",
    "\n",
    "**Created by:** Ieva Liepuoniute and Pauline Ollitrault\n",
    "\n",
    "**Version:** qiskit-nature 0.3.2"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "72792716",
   "metadata": {},
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Qiskit v0.35.0 (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.13"
  },
  "widgets": {
   "application/vnd.jupyter.widget-state+json": {
    "state": {},
    "version_major": 2,
    "version_minor": 0
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}