|
| 1 | +r"""Quantum Defect Embedding Theory (QDET) |
| 2 | +========================================= |
| 3 | +Many interesting problems in quantum chemistry and materials science feature a strongly correlated |
| 4 | +region embedded within a larger environment. Example of such systems include point defects in |
| 5 | +materials [#Galli]_, active site of catalysts [#SJRLee]_ and surface phenomenon such as adsorption |
| 6 | +[#Gagliardi]_. Such systems can be accurately simulated with **embedding theories**, which effectively |
| 7 | +capture the strong electronic correlations in the active region with high accuracy, while accounting |
| 8 | +for the environment in a more approximate manner. |
| 9 | +
|
| 10 | +In this demo, we show how to implement quantum defect embedding theory (QDET) [#Galli]_. This method |
| 11 | +has been successfully applied to study systems such as defects in calcium oxide [#Galli]_ and to calculate |
| 12 | +excitations of the negatively charged nitrogen-vacancy defect in diamond [#Galli2]_. QDET can be used to calculate ground states, |
| 13 | +excited states and dynamic properties of materials. These make QDET a powerful method for affordable quantum simulation |
| 14 | +of materials. Another important advantage |
| 15 | +of QDET is the compatibility of the method with quantum algorithms as we explain in the following |
| 16 | +sections. |
| 17 | +
|
| 18 | +.. figure:: ../_static/demo_thumbnails/opengraph_demo_thumbnails/OGthumbnail_how_to_build_qdet_hamiltonian.png |
| 19 | + :align: center |
| 20 | + :width: 70% |
| 21 | + :target: javascript:void(0) |
| 22 | +
|
| 23 | +""" |
| 24 | + |
| 25 | +############################################# |
| 26 | +# The main component of a QDET simulation is to construct an effective Hamiltonian that describes |
| 27 | +# the impurity subsystem and its interaction with the environment. In second quantization, the |
| 28 | +# effective Hamiltonian can be represented in terms of electronic creation, :math:`a^{\dagger}`, and |
| 29 | +# annihilation , :math:`a`, operators as |
| 30 | +# |
| 31 | +# .. math:: |
| 32 | +# |
| 33 | +# H^{eff} = \sum_{ij} t_{ij}^{eff}a_i^{\dagger}a_j + \frac{1}{2}\sum_{ijkl} v_{ijkl}^{eff}a_i^{\dagger}a_{j}^{\dagger}a_ka_l, |
| 34 | +# |
| 35 | +# where :math:`t_{ij}^{eff}` and :math:`v_{ijkl}^{eff}` represent the effective one-body and |
| 36 | +# two-body integrals, respectively, and the indices :math:`ijkl` span over the orbitals inside the impurity. |
| 37 | +# This Hamiltonian describes a simplified representation of the quantum system that is more |
| 38 | +# computationally tractable, while properly capturing the essential physics of the problem. |
| 39 | +# |
| 40 | +# Implementation |
| 41 | +# -------------- |
| 42 | +# A QDET simulation typically starts by obtaining a mean field approximation of the whole system |
| 43 | +# using efficient quantum chemistry methods such as density functional theory (DFT). These |
| 44 | +# calculations provide a set of orbitals that are partitioned into **impurity** and **bath** orbitals. |
| 45 | +# The effective Hamiltonian is constructed from the impurity orbitals and is subsequently solved |
| 46 | +# by using either a high accuracy classical method or a quantum algorithm. Let's implement these |
| 47 | +# steps for an example! |
| 48 | +# |
| 49 | +# Mean field calculations |
| 50 | +# ^^^^^^^^^^^^^^^^^^^^^^^ |
| 51 | +# We implement QDET to compute the excitation energies of a negatively charged nitrogen-vacancy |
| 52 | +# defect in diamond [#Galli2]_. We use DFT to obtain a mean field description of the whole system. |
| 53 | +# The DFT calculations are performed with the `QUANTUM ESPRESSO <https://www.quantum-espresso.org/>`_ |
| 54 | +# package. This requires downloading pseudopotentials [#Modji]_ for each atomic species |
| 55 | +# in the system from the QUANTUM ESPRESSO `database <https://www.quantum-espresso.org/pseudopotentials/>`_. |
| 56 | +# To prepare our system, the necessary carbon and nitrogen pseudopotentials can be downloaded by executing |
| 57 | +# the following commands through the terminal or command prompt: |
| 58 | +# |
| 59 | +# .. code-block:: bash |
| 60 | +# |
| 61 | +# wget -N -q http://www.quantum-simulation.org/potentials/sg15_oncv/upf/C_ONCV_PBE-1.2.upf |
| 62 | +# wget -N -q http://www.quantum-simulation.org/potentials/sg15_oncv/upf/N_ONCV_PBE-1.2.upf |
| 63 | +# |
| 64 | +# Next, we need to create the input file for running QUANTUM ESPRESSO. This contains |
| 65 | +# information about the system and the DFT calculations. More details on how to construct the input |
| 66 | +# file can be found in the QUANTUM ESPRESSO `documentation <https://www.quantum-espresso.org/Doc/INPUT_PW.html>`_ |
| 67 | +# page. For the system taken here, the input file can be downloaded with |
| 68 | +# |
| 69 | +# .. code-block:: bash |
| 70 | +# |
| 71 | +# wget -N -q https://west-code.org/doc/training/nv_diamond_63/pw.in |
| 72 | +# |
| 73 | +# DFT calculations can now be initiated using the ``pw.x`` executable in ``WEST``, taking ``pw.in`` as the input file |
| 74 | +# and directing the output to ``pw.out``. This process is parallelized across 2 cores using ``mpirun``. |
| 75 | +# |
| 76 | +# .. code-block:: bash |
| 77 | +# |
| 78 | +# mpirun -n 2 pw.x -i pw.in > pw.out |
| 79 | +# |
| 80 | +# Identify the impurity |
| 81 | +# ^^^^^^^^^^^^^^^^^^^^^ |
| 82 | +# Once we have obtained the mean field description, we can identify our impurity by finding |
| 83 | +# the states that are localized around the defect region in real space. To do that, we compute the |
| 84 | +# localization factor :math:`L_n` for each state ``n``, defined as: |
| 85 | +# |
| 86 | +# .. math:: |
| 87 | +# |
| 88 | +# L_n = \int_{V \in \Omega} d^3 r |\Psi_n(r)|^2, |
| 89 | +# |
| 90 | +# where :math:`V` is the identified volume including the impurity within the supercell volume |
| 91 | +# :math:`\Omega` and :math:`\Psi` is the wavefunction [#Galli2]_. We will use the |
| 92 | +# `WEST <https://pubs.acs.org/doi/10.1021/ct500958p>`_ program to compute the localization factor. |
| 93 | +# This requires the ``westpp.in`` input file, example for which is shown below. Here, we specify the |
| 94 | +# box parameters within which the localization factor is being computed; the vectors for this box are provided in |
| 95 | +# in atomic units as ``[x_start, x_end, y_start, y_end, z_start, z_end]``. |
| 96 | +# |
| 97 | +# .. rst-class:: sphx-glr-script-out |
| 98 | +# |
| 99 | +# .. code-block:: none |
| 100 | +# |
| 101 | +# westpp_control: |
| 102 | +# westpp_calculation: L # triggers the calculation of the localization factor |
| 103 | +# westpp_range: # defines the range of states to compute the localization factor |
| 104 | +# - 1 # start from the first state |
| 105 | +# - 176 # use all 176 states |
| 106 | +# westpp_box: # specifies the parameter of the box in atomic units for integration |
| 107 | +# - 6.19 |
| 108 | +# - 10.19 |
| 109 | +# - 6.28 |
| 110 | +# - 10.28 |
| 111 | +# - 6.28 |
| 112 | +# - 10.28 |
| 113 | +# |
| 114 | +# The calculation can now be performed by running the ``westpp.x`` executable from WEST using mpirun to |
| 115 | +# parallelize it across two cores. |
| 116 | +# |
| 117 | +# .. code-block:: bash |
| 118 | +# |
| 119 | +# mpirun -n 2 westpp.x -i westpp.in > westpp.out |
| 120 | +# |
| 121 | +# This creates the file ``westpp.json`` which contains the information we need here. Since |
| 122 | +# computational resources required to run the calculation are large, for the purpose of this tutorial we just |
| 123 | +# download a pre-computed file with: |
| 124 | +# |
| 125 | +# .. code-block:: bash |
| 126 | +# |
| 127 | +# mkdir -p west.westpp.save |
| 128 | +# wget -N -q https://west-code.org/doc/training/nv_diamond_63/box_westpp.json -O west.westpp.save/westpp.json |
| 129 | +# |
| 130 | +# We can plot the computed localization factor for each of the states: |
| 131 | +# |
| 132 | +# .. code-block:: bash |
| 133 | +# |
| 134 | +# import json |
| 135 | +# import numpy as np |
| 136 | +# import matplotlib.pyplot as plt |
| 137 | +# |
| 138 | +# with open('west.westpp.save/westpp.json','r') as f: |
| 139 | +# data = json.load(f) |
| 140 | +# |
| 141 | +# y = np.array(data['output']['L']['K000001']['local_factor'],dtype='f8') |
| 142 | +# x = np.array([i+1 for i in range(y.shape[0])]) |
| 143 | +# |
| 144 | +# plt.plot(x,y,'o') |
| 145 | +# plt.axhline(y=0.08,linestyle='--',color='red') |
| 146 | +# |
| 147 | +# plt.xlabel('Kohn-Sham orbital index') |
| 148 | +# plt.ylabel('Localization factor') |
| 149 | +# |
| 150 | +# plt.show() |
| 151 | +# |
| 152 | +# |
| 153 | +# .. figure:: ../_static/demonstration_assets/qdet/localization.jpeg |
| 154 | +# :align: center |
| 155 | +# :width: 70% |
| 156 | +# :target: javascript:void(0) |
| 157 | +# |
| 158 | +# From this plot, it is easy to see that the orbitals can be categorized as orbitals with low and |
| 159 | +# high localization factors. For the purpose of defining an impurity, we need highly localized |
| 160 | +# orbitals, so we set a cutoff of :math:`0.06`, illustrated by the red dashed line, and choose the orbitals that have a localization |
| 161 | +# factor larger than :math:`0.06`. We'll use these orbitals for the calculation of the parameters |
| 162 | +# needed to construct the effective Hamiltonian. |
| 163 | +# |
| 164 | +# Electronic Integrals |
| 165 | +# ^^^^^^^^^^^^^^^^^^^^ |
| 166 | +# The next step in QDET is to define the effective one-body and two-body integrals for the impurity. |
| 167 | +# The effective two-body integrals :math:`v^{eff}` are computed first as matrix elements of the |
| 168 | +# partially screened static Coulomb potential :math:`W_0^{R}`. |
| 169 | +# |
| 170 | +# .. math:: |
| 171 | +# |
| 172 | +# v_{ijkl}^{eff} = [W_0^{R}]_{ijkl}, |
| 173 | +# |
| 174 | +# where :math:`W_0^R` results from screening the bare Coulomb potential :math:`v` with the reduced |
| 175 | +# polarizability :math:`P_0^R = P - P_{imp}`, where :math:`P` is the system's polarizability and |
| 176 | +# :math:`P_{imp}` is the impurity's polarizability. Since solving the effective Hamiltonian |
| 177 | +# accounts for the exchange and correlation interactions between the active electrons, we remove |
| 178 | +# these interactions from the Kohn-Sham Hamiltonian :math:`H^{KS}` to avoid double counting them. |
| 179 | +# |
| 180 | +# The one-body term :math:`t^{eff}` is obtained by subtracting from the Kohn-Sham Hamiltonian the |
| 181 | +# double-counting term accounting for electrostatic and exchange-correlation interactions in the |
| 182 | +# active space. |
| 183 | +# |
| 184 | +# .. math:: |
| 185 | +# |
| 186 | +# t_{ij}^{eff} = H_{ij}^{KS} - t_{ij}^{dc}. |
| 187 | +# |
| 188 | +# We use the WEST program to compute these parameters. WEST will first compute the |
| 189 | +# quasiparticle energies, then the partially screened Coulomb potential, and finally |
| 190 | +# the parameters of the effective Hamiltonian. The software offers several execution modes through the |
| 191 | +# :code:`wfreq_calculation` keyword. We choose ``XWGQH`` to ensure the full workflow is executed, which computes |
| 192 | +# both the quasiparticle corrections and the final parameters for the QDET effective Hamiltonian. The input file |
| 193 | +# for such a calculation is shown below: |
| 194 | +# |
| 195 | +# .. rst-class:: sphx-glr-script-out |
| 196 | +# |
| 197 | +# .. code-block:: none |
| 198 | +# |
| 199 | +# wstat_control: |
| 200 | +# wstat_calculation: S # starts the calculation from scratch |
| 201 | +# n_pdep_eigen: 512 # number of eigenpotentials; matches number of electrons |
| 202 | +# trev_pdep: 0.00001 # convergence threshold for eigenvalues |
| 203 | +# |
| 204 | +# wfreq_control: |
| 205 | +# wfreq_calculation: XWGQH # compute quasiparticle corrections and Hamiltonian params |
| 206 | +# macropol_calculation: C # include long-wavelength limit for condensed systems |
| 207 | +# l_enable_off_diagonal: true # calculate off-diagonal elements of G_0-W_0 self-energy |
| 208 | +# n_pdep_eigen_to_use: 512 # number of PDEP eigenvectors to be used |
| 209 | +# qp_bands: [87,122,123,126,127,128] # impurity orbitals |
| 210 | +# n_refreq: 300 # number of frequencies on the real axis |
| 211 | +# ecut_refreq: 2.0 # cutoff for the real frequencies |
| 212 | +# |
| 213 | +# We can now execute the calculation with: |
| 214 | +# |
| 215 | +# .. code-block:: bash |
| 216 | +# |
| 217 | +# mpirun -n 2 wfreq.x -i wfreq.in > wfreq.out |
| 218 | +# |
| 219 | +# This calculation takes some time and requires computational resources, therefore we download a |
| 220 | +# pre-computed output file with |
| 221 | +# |
| 222 | +# .. code-block:: bash |
| 223 | +# |
| 224 | +# mkdir -p west.wfreq.save |
| 225 | +# wget -N -q https://west-code.org/doc/training/nv_diamond_63/wfreq.json -O west.wfreq.save/wfreq.json |
| 226 | +# |
| 227 | +# This output file contains all the information we need to construct the effective Hamiltonian. |
| 228 | +# |
| 229 | +# Effective Hamiltonian |
| 230 | +# ^^^^^^^^^^^^^^^^^^^^^ |
| 231 | +# We now construct the effective Hamiltonian by importing the electron integral results and using |
| 232 | +# WEST: |
| 233 | +# |
| 234 | +# .. code-block:: python |
| 235 | +# |
| 236 | +# from westpy.qdet import QDETResult |
| 237 | +# |
| 238 | +# effective_hamiltonian = QDETResult(filename='west.wfreq.save/wfreq.json') |
| 239 | +# |
| 240 | +# The effective Hamiltonian can be solved using a high level method such as the full configuration |
| 241 | +# interaction (FCI) algorithm from WEST as: |
| 242 | +# |
| 243 | +# .. code-block:: python |
| 244 | +# |
| 245 | +# solution = effective_hamiltonian.solve() |
| 246 | +# |
| 247 | +# Using :code:`solve()` prints the excitation energies, spin multiplicity and relative occupation of |
| 248 | +# the active orbitals. |
| 249 | +# |
| 250 | +# .. rst-class:: sphx-glr-script-out |
| 251 | +# |
| 252 | +# .. code-block:: none |
| 253 | +# |
| 254 | +# ====================================================================== |
| 255 | +# Building effective Hamiltonian... |
| 256 | +# nspin: 1 |
| 257 | +# occupations: [[2. 2. 2. 2. 1. 1.]] |
| 258 | +# ===================================================================== |
| 259 | +# diag[1RDM - 1RDM(GS)] |
| 260 | +# E [eV] char 87 122 123 126 127 128 |
| 261 | +# 0 0.000 3- 0.000 0.000 0.000 0.000 0.000 0.000 |
| 262 | +# 1 0.436 1- -0.001 -0.009 -0.018 -0.067 0.004 0.091 |
| 263 | +# 2 0.436 1- -0.001 -0.009 -0.018 -0.067 0.092 0.002 |
| 264 | +# 3 1.251 1- -0.002 -0.019 -0.023 -0.067 0.054 0.057 |
| 265 | +# 4 1.939 3- -0.003 -0.010 -0.127 -0.860 1.000 0.000 |
| 266 | +# 5 1.940 3- -0.003 -0.010 -0.127 -0.860 0.000 1.000 |
| 267 | +# 6 2.935 1- -0.000 -0.032 -0.043 -0.855 0.929 0.002 |
| 268 | +# 7 2.936 1- -0.000 -0.032 -0.043 -0.855 0.002 0.929 |
| 269 | +# 8 4.661 1- -0.006 -0.054 -0.188 -1.672 0.960 0.960 |
| 270 | +# 9 5.080 3- -0.014 -0.698 -0.213 -0.075 1.000 0.000 |
| 271 | +# ---------------------------------------------------------------------- |
| 272 | +# |
| 273 | +# The solution object is a dictionary that contains information about the FCI eigenstates of the |
| 274 | +# system, which includes various excitation energies, spin multiplicities, eigenvectors etc. |
| 275 | +# More importantly, while FCI handles small embedded effective Hamiltonians with ease, it quickly |
| 276 | +# hits a wall with larger impurities. This is precisely where quantum computing steps in, offering |
| 277 | +# the scalability needed to tackle such complex systems. The first step to solving these effective |
| 278 | +# Hamiltonians via quantum algorithms in PennyLane, is to convert them to qubit Hamiltonians. |
| 279 | +# |
| 280 | +# Quantum Simulation |
| 281 | +# ^^^^^^^^^^^^^^^^^^ |
| 282 | +# We now map the effective Hamiltonian to the qubit basis. Note that the two-electron integrals obtained |
| 283 | +# before are represented in Chemists' notation and need to be converted to a notation |
| 284 | +# that is compatible with PennyLane. Here's how to construct the `qubit Hamiltonian <https://pennylane.ai/qml/demos/tutorial_fermionic_operators>`__: |
| 285 | +# |
| 286 | +# .. code-block:: python |
| 287 | +# |
| 288 | +# from pennylane.qchem import one_particle, two_particle, observable |
| 289 | +# import numpy as np |
| 290 | +# |
| 291 | +# effective_hamiltonian = QDETResult(filename="west.wfreq.save/wfreq.json") |
| 292 | +# |
| 293 | +# one_e, two_e = effective_hamiltonian.h1e, effective_hamiltonian.eri |
| 294 | +# |
| 295 | +# t = one_particle(one_e[0]) |
| 296 | +# v = two_particle(np.swapaxes(two_e[0][0], 1, 3)) |
| 297 | +# qubit_op = observable([t, v], mapping="jordan_wigner") |
| 298 | +# print("Qubit Hamiltonian: ", qubit_op) |
| 299 | +# |
| 300 | +# .. rst-class:: sphx-glr-script-out |
| 301 | +# |
| 302 | +# .. code-block:: none |
| 303 | +# |
| 304 | +# Qubit Hamiltonian: (2.40331309905556+0j) * I(0) + (-0.12208093833046951+0j) * Z(0) + (-0.12208093833046951+0j) * Z(1) + (-0.003330743747901097+0j) * (Y(0) @ Z(1) @ Y(2)) + ... |
| 305 | +# |
| 306 | +# We can use this Hamiltonian in a quantum algorithm such as `quantum phase estimation (QPE) <https://pennylane.ai/qml/demos/tutorial_qpe>`__. |
| 307 | +# You can compare the results and verify that the computed energies from quantum algorithm |
| 308 | +# match those that we obtained before. |
| 309 | +# |
| 310 | +# Conclusion |
| 311 | +# ---------- |
| 312 | +# Quantum defect embedding theory is a novel framework for simulating strongly correlated |
| 313 | +# quantum systems and has been successfully used for studying defects in solids. Applicability of |
| 314 | +# QDET however is not limited to defects: it can be used for other systems where a strongly |
| 315 | +# correlated subsystem is embedded in a weakly correlated environment. Additionally, QDET is able to |
| 316 | +# correct the interaction double counting issue within the active space faced by a variety of |
| 317 | +# other embedding theories. The Green's function based formulation of QDET ensures |
| 318 | +# exact removal of double counting corrections at the `GW <https://en.wikipedia.org/wiki/GW_approximation>`__ |
| 319 | +# level of theory thus removing the |
| 320 | +# approximation present in the initial DFT based formulation. This formulation also helps to capture |
| 321 | +# the response properties and provides access to excited state properties. Another major advantage |
| 322 | +# of QDET is the ease with which it can be used with quantum computers in a hybrid framework [#Baker]_. |
| 323 | +# Given its ability to access excited states and its hybrid quantum compatibility, potential future extensions |
| 324 | +# could involve applying QDET to the wider field of spectroscopy, such as investigating |
| 325 | +# core-level excited states relevant to `X-ray Absorption Spectroscopy (XAS) <https://pennylane.ai/qml/demos/tutorial_xas>`__. |
| 326 | +# In conclusion, QDET is a powerful embedding approach for simulating complex quantum systems. |
| 327 | +# |
| 328 | +# References |
| 329 | +# ---------- |
| 330 | +# |
| 331 | +# .. [#Galli] |
| 332 | +# J. Davidsson, M. Onizhuk, C. Vorwerk, G. Galli, |
| 333 | +# "Discovery of atomic clock-like spin defects in simple oxides from first principles", |
| 334 | +# `arXiv:2302.07523 <https://arxiv.org/abs/2302.07523>`__. |
| 335 | +# |
| 336 | +# .. [#SJRLee] |
| 337 | +# S. J. R. Lee, F. Ding, F. R. Manby, T. F. Miller III, |
| 338 | +# "Analytical Gradients for Projection-Based Wavefunction-in-DFT Embedding" |
| 339 | +# `arXiv:1903.05830 <https://arxiv.org/abs/1903.05830>`__. |
| 340 | +# |
| 341 | +# .. [#Gagliardi] |
| 342 | +# A. Mitra, Matthew Hermes, M. Cho, V. Agarawal, L. Gagliardi, |
| 343 | +# "Periodic Density Matrix Embedding for CO Adsorption on the MgO(001)Surface" |
| 344 | +# `J. Phys. Chem. Lett. 2022, 13, 7483 <https://pubs.acs.org/doi/10.1021/acs.jpclett.2c01915>`__. |
| 345 | +# |
| 346 | +# .. [#Galli2] |
| 347 | +# N. Sheng, C. Vorwerk, M. Govoni, G. Galli, |
| 348 | +# "Green's function formulation of quantum defect embedding theory", |
| 349 | +# `arXiv:2203.05493 <https://arxiv.org/abs/2203.05493>`__. |
| 350 | +# |
| 351 | +# .. [#Modji] |
| 352 | +# M. S. Zini, A. Delgado, *et al.*, |
| 353 | +# "Quantum simulation of battery materials using ionic pseudopotentials" |
| 354 | +# `arXiv:2302.07981 <https://arxiv.org/abs/2302.07981>`__. |
| 355 | +# |
| 356 | +# .. [#Baker] |
| 357 | +# Jack S. Baker, Pablo A. M. Casares, *et al.*, |
| 358 | +# "Simulating optically-active spin defects with a quantum computer" |
| 359 | +# `arXiv:2405.13115 <https://arxiv.org/abs/2405.13115>`__. |
| 360 | +# |
0 commit comments