{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Building Molecular Hamiltonian\n",
    "*Copyright (c) 2021 Institute for Quantum Computing, Baidu Inc. All Rights Reserved.*"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Overview\n",
    "\n",
    "In this tutorial, we will demonstrate how to use Paddle Quantum's *qchem* module to build valid Hamiltonian for simulating chemical molecules on a quantum computer. We will go step by step how to build the second quantized Hamiltonian from a molecular structure and how to transform it to a set of Pauli matrices. \n",
    "\n",
    "Hamiltonian is a physical quantity related to the total energy of a physical system. In general, it can be represented as \n",
    "\n",
    "$$\n",
    "\\hat{H}=\\hat{T}+\\hat{V},\\tag{1}\n",
    "$$\n",
    "\n",
    "where $\\hat{T}$ is the kinetic energy and $\\hat{V}$ is the potential energy. Hamiltonian is useful for various quantum algorithms, such as [variational quantum eigensolver](./VQE_EN.ipynb) and [Hamiltonian Simulation with Product Formula](./HamiltonianSimulation_EN.ipynb).\n",
    "\n",
    "When trying to solve a chemistry problem with quantum mechanics, we also need to write down a Hamiltonian that describes the chemical system involved in the problem. Starting from this Hamiltonian, we can, in principle, calculate the ground state and excited states, and use the information to further explore all the physical properties of the quantum system. The dominant Hamiltonian of electronic problems has the form\n",
    "\n",
    "$$\n",
    "\\hat{H}=\\sum_{i=1}^N\\left(-\\frac{1}{2}\\nabla_{x_i}^2\\right)+\\sum_{i=1}^N\\sum_{j< i}\\frac{1}{|x_i-x_j|}-\\sum_{i=1}^N\\sum_{I=1}^M\\frac{Z_I}{|x_i-R_I|},\\tag{2}\n",
    "$$\n",
    "\n",
    "when we use [atomic units](https://en.wikipedia.org/wiki/Hartree_atomic_units). Our electronic problem contains $N$ electrons and $M$ nucleus. We use $x_i$ to denote position of the $i$-th electron, and use $R_I$ to denote position of the $I$-th nuclei. \n",
    "\n",
    "This tutorial will have the following parts. Let's first talk about how to construct a molecule in *qchem*. After that, we will briefly describe how to calculate [Hartree Fock](https://en.wikipedia.org/wiki/Hartree%E2%80%93Fock_method) single particle orbitals by calling external quantum chemistry within Paddle Quantum. Next, we show how we can obtain molecular Hamiltonian using *qchem*'s `Molecule` class."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Defining the molecular structure\n",
    "In this example, we show how to construct water molecule from its chemical formula and coordinates of atoms. \n",
    "\n",
    "![h2o.png](figures/buildingmolecule-fig-h2o.png)\n",
    "\n",
    "Within Paddle Quantum, we specify the atom as a list whose first element is the atomic symbol and the second element is another list that contains its Cartesian coordinate. The molecule is thus a bigger list composed of atoms' list.\n",
    "\n",
    "**Note: As to the environment setting, please refer to [README.md](https://github.com/PaddlePaddle/Quantum/blob/master/README.md).**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Eliminate noisy python warnings\n",
    "import warnings\n",
    "warnings.filterwarnings(\"ignore\")\n",
    "import logging\n",
    "logging.basicConfig(level=logging.INFO)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# in Angstrom\n",
    "h2o_coords = [(\"H\", [-0.02111417,0.8350417,1.47688078]),  # H stands for hydrogen element in water\n",
    "                        (\"O\", [0.0, 0.0, 0.0]),                     # O stands for oxygen element in water\n",
    "                        (\"H\", [-0.00201087,0.45191737,-0.27300254])]"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Calculate Hartree Fock orbitals\n",
    "Hartree Fock method uses the [Slater determinant](https://en.wikipedia.org/wiki/Slater_determinant) to represent the $N$-electron wavefunction. It could provide us with a set of single particle orbitals which are often taken as input to more advanced quantum chemistry methods. \n",
    "\n",
    "Currently, Paddle Quantum uses PySCF [1] as its quantum chemistry engine, we will support more quantum chemistry toolkits, such as psi4, in the future (NOTE: PySCF currently only support Mac and Linux platform). We could use the `PySCFDriver` provided in `qchem` module to manage the quantum chemistry calculation and get the necessary information about the molecule. It takes molecular structure, total molecular charge, and spin multiplicity as its major inputs. We can also control the precision of Hartree Fock calculation by setting different [basis set](https://en.wikipedia.org/wiki/Basis_set_(chemistry)) for `basis` keyword. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:root:\n",
      "#######################################\n",
      "SCF Calculation (Classical)\n",
      "#######################################\n",
      "INFO:root:Basis: sto-3g\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "converged SCF energy = -73.9677038774737\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:root:SCF energy: -73.96770.\n"
     ]
    }
   ],
   "source": [
    "from paddle_quantum.qchem import PySCFDriver\n",
    "\n",
    "driver = PySCFDriver()\n",
    "driver.load_molecule(\n",
    "    atom=h2o_coords,         # Geometry of H2O molecule\n",
    "    basis=\"sto-3g\",          # Basis set for quantum chemistry calculation\n",
    "    multiplicity=1,          # Spin multiplicity for molecule, since the total spin of H2O is S=0，its spin multiplicity is 2S+1=1\n",
    "    charge=0,                 # Total charge of molecule, since H2O is charge neutral, its charge=0\n",
    "    unit=\"Angstrom\"\n",
    ")\n",
    "driver.run_scf()             # Perform Hartree Fock calculation"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Molecular Hamiltonian in second quantization form\n",
    "When we study many electron quantum systems, it's often convenient to write Hamiltonian at the beginning of this tutorial in [second quantization](https://en.wikipedia.org/wiki/Second_quantization) representation \n",
    "\n",
    "$$\n",
    "\\hat{H}=\\sum_{p,q}h_{pq}\\hat{c}^{\\dagger}_p\\hat{c}_q+\\frac{1}{2}\\sum_{p,q,r,s}v_{pqrs}\\hat{c}^{\\dagger}_p\\hat{c}^{\\dagger}_q\\hat{c}_r\\hat{c}_s,\\tag{3}$$\n",
    "\n",
    "where $p$, $q$, $r$ and $s$ are Hartree Fock orbitals computed in the previous section. $\\hat{c}^{\\dagger}_p$ and $\\hat{c}_q$ are creation and annihilation operations, respectively. The two coefficients $h_{pq}$ and $v_{pqrs}$ are called molecular integrals, and can be obtained from `PySCFDriver` in the following way."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-3.2911e+01  5.5623e-01  2.8755e-01 -5.1653e-15 -7.4568e-02 -9.4552e-02 -2.8670e-01]\n",
      " [ 5.5623e-01 -8.0729e+00 -4.0904e-02  1.5532e-15  1.7890e-01  3.5048e-01  1.3460e+00]\n",
      " [ 2.8755e-01 -4.0904e-02 -7.3355e+00  6.8611e-16  4.1911e-01  5.2109e-01 -7.0928e-01]\n",
      " [-5.1650e-15  1.5617e-15  6.7472e-16 -7.5108e+00  5.4096e-16  1.7947e-15  4.5448e-16]\n",
      " [-7.4568e-02  1.7890e-01  4.1911e-01  5.3561e-16 -5.7849e+00  2.0887e+00 -1.2427e-01]\n",
      " [-9.4552e-02  3.5048e-01  5.2109e-01  1.8235e-15  2.0887e+00 -5.0803e+00 -1.3967e-02]\n",
      " [-2.8670e-01  1.3460e+00 -7.0928e-01  4.8475e-16 -1.2427e-01 -1.3967e-02 -5.0218e+00]]\n"
     ]
    }
   ],
   "source": [
    "import numpy as np \n",
    "np.set_printoptions(precision=4, linewidth=150)\n",
    "\n",
    "hpq = driver.get_onebody_tensor(\"int1e_kin\") + driver.get_onebody_tensor(\"int1e_nuc\")\n",
    "vpqrs = driver.get_twobody_tensor()\n",
    "assert np.shape(hpq)==(7, 7)             # H2O has 7 orbitals when using STO-3G basis.\n",
    "assert np.shape(vpqrs)==(7, 7, 7, 7)\n",
    "\n",
    "print(hpq)\n",
    "# print(vpqrs)  # vpqrs is a large tensor，for brevity, we comment this line by default, interested users can uncomment this line to see it."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Most of the time, we don't need to extract those integrals and assemble the Hamiltonian manually, *qchem* has already ecapsulated this procedure in `Molecule` class, we can get the Hamiltonian of water molecule from it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:root:\n",
      "#######################################\n",
      "SCF Calculation (Classical)\n",
      "#######################################\n",
      "INFO:root:Basis: sto-3g\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "converged SCF energy = -73.9677038774737\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:root:SCF energy: -73.96770.\n"
     ]
    }
   ],
   "source": [
    "# Build `Molecule` for water molecule.\n",
    "from paddle_quantum.qchem import Molecule\n",
    "\n",
    "mol = Molecule(\n",
    "    geometry=h2o_coords,\n",
    "    basis=\"sto-3g\",\n",
    "    driver=PySCFDriver()\n",
    ")\n",
    "\n",
    "\n",
    "H_of_water = mol.get_molecular_hamiltonian()"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since quantum computing only allows operation based on Pauli operators $\\hat{\\sigma}_x$, $\\hat{\\sigma}_y$, $\\hat{\\sigma}_z$, we have to transform Hamiltonian in above expression to a form represented by Pauli operators. This can be achieved via [Jordan-Wigner 变换](https://en.wikipedia.org/wiki/Jordan%E2%80%93Wigner_transformation) and we have already integrated it in `get_molecular_hamiltonian` method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "There are 2110 terms in H2O Hamiltonian in total.\n",
      "The first 10 terms are \n",
      " -44.808115297466266 I\n",
      "12.537584025014317 Z0\n",
      "12.537584025014313 Z1\n",
      "1.8115178684479893 Z2\n",
      "1.8115178684479893 Z3\n",
      "1.4546840922727915 Z4\n",
      "1.4546840922727915 Z5\n",
      "1.4299903414012243 Z6\n",
      "1.429990341401224 Z7\n",
      "1.0849294605602529 Z8\n"
     ]
    }
   ],
   "source": [
    "print('There are', H_of_water.n_terms, 'terms in H2O Hamiltonian in total.')\n",
    "print('The first 10 terms are \\n', H_of_water[:10])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Great! Now you know how to build a proper Hamiltonian from a given molecular structure, let's move further and see how to use [variational quantum eigensolver](./VQE_EN.ipynb) (VQE) to determine the ground state of hydrogen molecule.\n"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "## References\n",
    "\n",
    "[1] [PySCF: Python-based Simulations of Chemistry Framework](https://pyscf.org/)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "modellib",
   "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.7.15"
  },
  "vscode": {
   "interpreter": {
    "hash": "8f24120f890011f53feb4ed62c47961d8565ec1de8b7cb23548c15bd6da8f2d2"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}