{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Excited states solvers" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Introduction\n", "\n", "\n", "\n", "In this tutorial we are going to discuss the excited states calculation interface of Qiskit Nature. The goal is to compute the excited states of a molecular Hamiltonian. This Hamiltonian can be electronic or vibrational. To know more about the preparation of the Hamiltonian, check out the [Electronic structure](01_electronic_structure.ipynb) and [Vibrational structure tutorials](02_vibrational_structure.ipynb).\n", "\n", "The first step is to define the molecular system. In the following we ask for the electronic part of a hydrogen molecule." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "from qiskit_nature.units import DistanceUnit\n", "from qiskit_nature.second_q.drivers import PySCFDriver\n", "\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", "es_problem = driver.run()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We will also be sticking to the Jordan-Wigner mapping. To learn more about the various mappers available in Qiskit Nature, check out the [Qubit Mappers tutorial](06_qubit_mappers.ipynb)." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "from qiskit_nature.second_q.mappers import JordanWignerMapper\n", "\n", "mapper = JordanWignerMapper()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## The Solver\n", "\n", "After these steps we need to define a solver. The solver is the algorithm through which the excited states are computed. \n", "\n", "Let's first start with a purely classical example: the `NumPyEigensolver`. This algorithm exactly diagonalizes the Hamiltonian. Although it scales badly, it can be used on small systems to check the results of the quantum algorithms. \n", "Here, we are only interested to look at eigenstates with a given number of particles. To compute only those states a filter function can be passed to the `NumPyEigensolver`. A default filter function is already implemented in Qiskit Nature which you can use for this purpose.\n", "\n", "We also need to specify the number of eigenvalues to be computed by the `NumPyEigensolver`. For this particular system, we are interested in the ground and first three excited states, so we will set `k=4` (which defaults to 1 so be sure to set this, otherwise you will only obtain the ground state!)." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "from qiskit_algorithms import NumPyEigensolver\n", "\n", "numpy_solver = NumPyEigensolver(k=4, filter_criterion=es_problem.get_default_filter_criterion())" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The excitation energies can also be accessed with the [qEOM algorithm](https://journals.aps.org/prresearch/abstract/10.1103/PhysRevResearch.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", "Although the previous equation can be solved classically, each matrix element must be measured on the quantum computer with the corresponding ground state. \n", "To use the qEOM as a solver in Qiskit, we have to define a ground state calculation first, which will provide the required ground state information to the algorithm. With this the qEOM solver can be initialized like so:" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "from qiskit_algorithms import VQE\n", "from qiskit_algorithms.optimizers import SLSQP\n", "from qiskit.primitives import Estimator\n", "from qiskit_nature.second_q.algorithms import GroundStateEigensolver, QEOM, EvaluationRule\n", "from qiskit_nature.second_q.circuit.library import HartreeFock, UCCSD\n", "\n", "ansatz = UCCSD(\n", " es_problem.num_spatial_orbitals,\n", " es_problem.num_particles,\n", " mapper,\n", " initial_state=HartreeFock(\n", " es_problem.num_spatial_orbitals,\n", " es_problem.num_particles,\n", " mapper,\n", " ),\n", ")\n", "\n", "estimator = Estimator()\n", "# This first part sets the ground state solver\n", "# see more about this part in the ground state calculation tutorial\n", "solver = VQE(estimator, ansatz, SLSQP())\n", "solver.initial_point = [0.0] * ansatz.num_parameters\n", "gse = GroundStateEigensolver(mapper, solver)\n", "\n", "# The qEOM algorithm is simply instantiated with the chosen ground state solver and Estimator primitive\n", "qeom_excited_states_solver = QEOM(gse, estimator, \"sd\", EvaluationRule.ALL)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## The calculation and results\n", "\n", "We are now ready to compute the results. Below, we are comparing the results obtained by the exact `NumPyEigensolver` with the default filter criterion enabled with the results obtained by the qEOM algorithm." ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "scrolled": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "=== GROUND STATE ENERGY ===\n", " \n", "* Electronic ground state energy (Hartree): -1.857275030202\n", " - computed part: -1.857275030202\n", "~ Nuclear repulsion energy (Hartree): 0.719968994449\n", "> Total ground state energy (Hartree): -1.137306035753\n", " \n", "=== EXCITED STATE ENERGIES ===\n", " \n", " 1: \n", "* Electronic excited state energy (Hartree): -0.882722150245\n", "> Total excited state energy (Hartree): -0.162753155796\n", " 2: \n", "* Electronic excited state energy (Hartree): -0.224911252831\n", "> Total excited state energy (Hartree): 0.495057741618\n", " \n", "=== MEASURED OBSERVABLES ===\n", " \n", " 0: # Particles: 2.000 S: 0.000 S^2: 0.000 M: 0.000\n", " 1: # Particles: 2.000 S: 0.000 S^2: 0.000 M: 0.000\n", " 2: # Particles: 2.000 S: 0.000 S^2: 0.000 M: 0.000\n", " \n", "=== DIPOLE MOMENTS ===\n", " \n", "~ Nuclear dipole moment (a.u.): [0.0 0.0 1.3889487]\n", " \n", " 0: \n", " * Electronic dipole moment (a.u.): [0.0 0.0 1.388948701555]\n", " - computed part: [0.0 0.0 1.388948701555]\n", " > Dipole moment (a.u.): [0.0 0.0 -0.000000001555] Total: 0.000000001555\n", " (debye): [0.0 0.0 -0.000000003953] Total: 0.000000003953\n", " \n", " 1: \n", " * Electronic dipole moment (a.u.): [0.0 0.0 1.388948701555]\n", " - computed part: [0.0 0.0 1.388948701555]\n", " > Dipole moment (a.u.): [0.0 0.0 -0.000000001555] Total: 0.000000001555\n", " (debye): [0.0 0.0 -0.000000003953] Total: 0.000000003953\n", " \n", " 2: \n", " * Electronic dipole moment (a.u.): [0.0 0.0 1.388948701555]\n", " - computed part: [0.0 0.0 1.388948701555]\n", " > Dipole moment (a.u.): [0.0 0.0 -0.000000001555] Total: 0.000000001555\n", " (debye): [0.0 0.0 -0.000000003953] Total: 0.000000003953\n", " \n", "\n", "\n", "\n", "=== GROUND STATE ENERGY ===\n", " \n", "* Electronic ground state energy (Hartree): -1.857275030145\n", " - computed part: -1.857275030145\n", "~ Nuclear repulsion energy (Hartree): 0.719968994449\n", "> Total ground state energy (Hartree): -1.137306035696\n", " \n", "=== EXCITED STATE ENERGIES ===\n", " \n", " 1: \n", "* Electronic excited state energy (Hartree): -1.244586756145\n", "> Total excited state energy (Hartree): -0.524617761696\n", " 2: \n", "* Electronic excited state energy (Hartree): -0.882724356546\n", "> Total excited state energy (Hartree): -0.162755362097\n", " 3: \n", "* Electronic excited state energy (Hartree): -0.224913459141\n", "> Total excited state energy (Hartree): 0.495055535308\n", " \n", "=== MEASURED OBSERVABLES ===\n", " \n", " 0: # Particles: 2.000 S: 0.000 S^2: 0.000 M: 0.000\n", " 1: # Particles: 2.000 S: 1.000 S^2: 2.000 M: 0.000\n", " 2: # Particles: 2.000 S: 0.000 S^2: -0.000 M: 0.000\n", " 3: # Particles: 2.000 S: 0.000 S^2: 0.000 M: -0.000\n", " \n", "=== DIPOLE MOMENTS ===\n", " \n", "~ Nuclear dipole moment (a.u.): [0.0 0.0 1.3889487]\n", " \n", " 0: \n", " * Electronic dipole moment (a.u.): [0.0 0.0 1.388948961657]\n", " - computed part: [0.0 0.0 1.388948961657]\n", " > Dipole moment (a.u.): [0.0 0.0 -0.000000261657] Total: 0.000000261657\n", " (debye): [0.0 0.0 -0.000000665065] Total: 0.000000665065\n", " \n", " 1: \n", " * Electronic dipole moment (a.u.): [0.0 0.0 1.388948701647]\n", " - computed part: [0.0 0.0 1.388948701647]\n", " > Dipole moment (a.u.): [0.0 0.0 -0.000000001647] Total: 0.000000001647\n", " (debye): [0.0 0.0 -0.000000004186] Total: 0.000000004186\n", " \n", " 2: \n", " * Electronic dipole moment (a.u.): [0.0 0.0 1.388949035771]\n", " - computed part: [0.0 0.0 1.388949035771]\n", " > Dipole moment (a.u.): [0.0 0.0 -0.000000335771] Total: 0.000000335771\n", " (debye): [0.0 0.0 -0.000000853445] Total: 0.000000853445\n", " \n", " 3: \n", " * Electronic dipole moment (a.u.): [0.0 0.0 1.388948738594]\n", " - computed part: [0.0 0.0 1.388948738594]\n", " > Dipole moment (a.u.): [0.0 0.0 -0.000000038594] Total: 0.000000038594\n", " (debye): [0.0 0.0 -0.000000098097] Total: 0.000000098097\n", " \n" ] } ], "source": [ "from qiskit_nature.second_q.algorithms import ExcitedStatesEigensolver\n", "\n", "numpy_excited_states_solver = ExcitedStatesEigensolver(mapper, numpy_solver)\n", "numpy_results = numpy_excited_states_solver.solve(es_problem)\n", "\n", "qeom_results = qeom_excited_states_solver.solve(es_problem)\n", "\n", "print(numpy_results)\n", "print(\"\\n\\n\")\n", "print(qeom_results)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "One can see from these results that one state is missing from the NumPy results. The reason for this is because the spin is also used as a filter and only singlet states are shown. \n", "In the following we use a custom filter function to check our results consistently and only filter out states with the incorrect number of particles (in this case the number of particle is 2) as well as the wrong magnetization (which we enforce to be 0)." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "=== GROUND STATE ENERGY ===\n", " \n", "* Electronic ground state energy (Hartree): -1.857275030202\n", " - computed part: -1.857275030202\n", "~ Nuclear repulsion energy (Hartree): 0.719968994449\n", "> Total ground state energy (Hartree): -1.137306035753\n", " \n", "=== EXCITED STATE ENERGIES ===\n", " \n", " 1: \n", "* Electronic excited state energy (Hartree): -1.244584549813\n", "> Total excited state energy (Hartree): -0.524615555364\n", " 2: \n", "* Electronic excited state energy (Hartree): -0.882722150245\n", "> Total excited state energy (Hartree): -0.162753155796\n", " 3: \n", "* Electronic excited state energy (Hartree): -0.224911252831\n", "> Total excited state energy (Hartree): 0.495057741618\n", " \n", "=== MEASURED OBSERVABLES ===\n", " \n", " 0: # Particles: 2.000 S: 0.000 S^2: 0.000 M: 0.000\n", " 1: # Particles: 2.000 S: 1.000 S^2: 2.000 M: 0.000\n", " 2: # Particles: 2.000 S: 0.000 S^2: 0.000 M: 0.000\n", " 3: # Particles: 2.000 S: 0.000 S^2: 0.000 M: 0.000\n", " \n", "=== DIPOLE MOMENTS ===\n", " \n", "~ Nuclear dipole moment (a.u.): [0.0 0.0 1.3889487]\n", " \n", " 0: \n", " * Electronic dipole moment (a.u.): [0.0 0.0 1.388948701555]\n", " - computed part: [0.0 0.0 1.388948701555]\n", " > Dipole moment (a.u.): [0.0 0.0 -0.000000001555] Total: 0.000000001555\n", " (debye): [0.0 0.0 -0.000000003953] Total: 0.000000003953\n", " \n", " 1: \n", " * Electronic dipole moment (a.u.): [0.0 0.0 1.388948701555]\n", " - computed part: [0.0 0.0 1.388948701555]\n", " > Dipole moment (a.u.): [0.0 0.0 -0.000000001555] Total: 0.000000001555\n", " (debye): [0.0 0.0 -0.000000003953] Total: 0.000000003953\n", " \n", " 2: \n", " * Electronic dipole moment (a.u.): [0.0 0.0 1.388948701555]\n", " - computed part: [0.0 0.0 1.388948701555]\n", " > Dipole moment (a.u.): [0.0 0.0 -0.000000001555] Total: 0.000000001555\n", " (debye): [0.0 0.0 -0.000000003953] Total: 0.000000003953\n", " \n", " 3: \n", " * Electronic dipole moment (a.u.): [0.0 0.0 1.388948701555]\n", " - computed part: [0.0 0.0 1.388948701555]\n", " > Dipole moment (a.u.): [0.0 0.0 -0.000000001555] Total: 0.000000001555\n", " (debye): [0.0 0.0 -0.000000003953] Total: 0.000000003953\n", " \n" ] } ], "source": [ "import numpy as np\n", "\n", "\n", "def filter_criterion(eigenstate, eigenvalue, aux_values):\n", " return np.isclose(aux_values[\"ParticleNumber\"][0], 2.0) and np.isclose(\n", " aux_values[\"Magnetization\"][0], 0.0\n", " )\n", "\n", "\n", "new_numpy_solver = NumPyEigensolver(k=4, filter_criterion=filter_criterion)\n", "new_numpy_excited_states_solver = ExcitedStatesEigensolver(mapper, new_numpy_solver)\n", "new_numpy_results = new_numpy_excited_states_solver.solve(es_problem)\n", "\n", "print(new_numpy_results)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/html": [ "

Version Information

Qiskit SoftwareVersion
qiskit-terra0.24.1
qiskit-nature0.7.0
System information
Python version3.10.11
Python compilerGCC 12.2.1 20221121 (Red Hat 12.2.1-4)
Python buildmain, Apr 5 2023 00:00:00
OSLinux
CPUs8
Memory (Gb)62.48404312133789
Wed Jun 07 10:45:42 2023 CEST
" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "

This code is a part of Qiskit

© Copyright IBM 2017, 2023.

This code is licensed under the Apache License, Version 2.0. You may
obtain a copy of this license in the LICENSE.txt file in the root directory
of this source tree or at http://www.apache.org/licenses/LICENSE-2.0.

Any modifications or derivative works of this code must retain this
copyright notice, and modified files need to carry a notice indicating
that they have been altered from the originals.

" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "import tutorial_magics\n", "\n", "%qiskit_version_table\n", "%qiskit_copyright" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.10.11" } }, "nbformat": 4, "nbformat_minor": 2 }