Note

Cette page a été générée à partir de tutorials/chemistry/02_vibronic_structure.ipynb.

Exécuter en mode interactif dans le IBM Quantum lab.

Structure vibronique

Introduction

L’hamiltonien moléculaire est

\[\mathcal{H} = - \sum_I \frac{\nabla_{R_I}^2}{M_I} - \sum_i \frac{\nabla_{r_i}^2}{m_e} - \sum_I\sum_i \frac{Z_I e^2}{|R_I-r_i|} + \sum_i \sum_{j>i} \frac{e^2}{|r_i-r_j|} + \sum_I\sum_{J>I} \frac{Z_I Z_J e^2}{|R_I-R_J|}\]

Because the nuclei are much heavier than the electrons they do not move on the same time scale and therefore, the behavior of nuclei and electrons can be decoupled. This is the Born-Oppenheimer approximation.

Dans le cadre de l’approximation de Born-Oppenheimer, une fonction d’onde moléculaire est factorisée en tant que produit d’une partie électronique, qui est la solution de l’équation électronique de Schroedinger, et d’une solution vibro-rotationnelle, qui est la solution de l’équation de Schroedinger nucléaire dans la surface d’énergie potentielle (PSE) générée par l’échantillonnage des valeurs propres de l’équation de Schroedinger électronique pour différentes géométries.

The nuclear Schroedinger equation is usually solved in two steps, in analogy with its electronic counterpart. A single-particle basis (the basis functions are called, in this case, modals) is obtained either by the harmonic approximation applied to the PES or from a vibrational self-consistent field (VSCF) calculation. Vibrational anharmonic correlations are added a-posteriori with perturbative or variational approaches. The latter include Vibrational Configuration Interaction (VCI) and Vibrational Coupled Cluster (VCC) for highly-accurate anharmonic energies. The main advantage of VCI and VCC over alternative approaches (such as perturbation theories) is that their accuracy can be systematically improved towards the complete basis set limit for a given PES. However, their applicability is limited to small molecules with up to about 10 atoms due to their unfavorable scaling with system size.

Pour résoudre le problème d’échelle, nous aimerions utiliser des algorithmes quantiques.

L’équation de Schroedinger nucléaire est

\[\mathcal{H}_{\text{vib}} |\Psi_{n}\rangle = E_{n} |\Psi_{n}\rangle\]

On appelle l’Hamiltonien Watson (en négligeant les termes de couplage vibro-rotationnel)

\[\mathcal{H}_\text{vib}(Q_1, \ldots, Q_L) =-\frac{1}{2} \sum_ { l= 1 } ^{L} \frac{ \partial ^ 2 } { \partial Q_l ^ 2 } + V (Q_1, \ldots, Q_L)\]

où \(Q_l\) sont les coordonnées normales pondérées par la masse harmonique.

\(\mathcal{H}_\text{vib}\) doit être mappé à un opérateur qui agit sur les états d’un ensemble donné de \(N_q\) qubits afin de calculer ses fonctions propres sur le matériel quantique. Dans les calculs de structure électronique, la cartographie est réalisée par l’expression de l’hamiltonien électronique non relativiste dans la seconde quantification, \textit{i.e.}, en le projetant sur l’ensemble complet des vecteurs d’occupation antisymétriques (ONV) générés par un ensemble d’orbitales donné (fini). Pour coder l’hamiltonien vibrationnel dans un second opérateur de quantification analogique, nous étendons le potentiel \(V (Q_1, \ldots, Q_L)\) avec l’extension \(n\)-body comme suit:

\[V(Q_1, \ldots, Q_L) = V_0 + \sum_{l=1}^L V^{[l]}(Q_l) + \sum_{l<m}^L V^{[l,m]}(Q_l, Q_m) + \sum_{l<m<n}^L V^{[l,m,n]}(Q_l, Q_m, Q_n) + \ldots\]

where \(V_0\) is the electronic energy of the reference geometry, the one-mode term \(V^{[l]}(Q_l)\) represents the variation of the PES upon change of the \(l\)-th normal coordinate from the equilibrium position. Similarly, the two-body potential \(V^{[l,m]}(Q_l, Q_m)\) represents the change in the exact PES upon a simultaneous displacement along the \(l\)-th and \(m\)-th coordinates. Often, including terms up to three-body in the \(L\)-body expansion is sufficient to obtain an accuracy of about 1~cm\(^{-1}\). We highlight that the many-body expansion of the potential operator defining the Watson Hamiltonian contains arbitrarily high coupling terms. This is a crucial difference compared to the non-relativistic electronic-structure Hamiltonian that contains only pairwise interactions.

Une seconde forme de quantification flexible du Hamiltonien Watson est obtenue dans le cadre de la représentation dite n-mode. Supposons que chaque mode \(l\) est décrit par un ensemble de base \(N_l\)-set \(S_l\) défini comme suit :

\[\mathcal{S}_l = \{ \phi_1^{(l)} (Q_l) , \ldots , \phi_{N_l}^{(l)} (Q_l) \} \, .\]

La fonction d’onde \(n\) en mode peut être étendue dans la base du produit \(\mathcal{S} = \otimes_{i=1}^L \mathcal{S}_i\) comme l’extension suivante de type CI:

\[|\Psi\rangle = \sum_{k_1=1}^{N_1} \cdots \sum_{k_L=1}^{N_L} C_{k_1,\ldots,k_L} \phi_{k_1}^{(1)}(Q_1) \cdots \phi_{k_L}^{(L)}(Q_L) \, ,\]

La fonction de base à plusieurs corps :math:` phi_{k_1}^ { (1) } (Q_ 1) cdots phi_{k_L}^ { (L) } (Q_L) ` sont encodées dans la seconde quantification dite mode :math:` n ` en tant que vecteurs de numéro d’occupation (ONV) comme suit:

\[\phi_{k_1}(Q_ 1) \cdots \phi_{k_L}(Q_L) \equiv | 0_1 \cdots 1_{k_1} \cdots 0_{N_1}, 0_1 \cdots 1_{k_2} \cdots 0_{N_2}, \cdots, 0_1 \cdots 1_{k_L} \cdots 0_{N_L}\rangle \,.\]

L’ONV défini ci-dessus est donc le produit des ONV spécifiques au mode :math:` L”, chacun décrivant un mode individuel. Puisque chaque mode est décrit par une seule et unique fonction de base, l’occupation de chaque ONV spécifique est une occupation. D’un point de vue théorique, chaque mode peut être interprété comme une quasi-particule distincte (définie comme des phonons en physique de l’état solide). La distinction découle du fait que le PSE n’est pas invariant par la permutation de deux modes, aussi dans ce cas, contrairement à l’interaction coulombienne entre deux particules égales. De ce point de vue, une molécule peut être interprétée comme une collection de particules indifférenciables de :math: »L » qui interagissent par l’intermédiaire de l’opérateur PES.

Basé sur cette représentation de deuxième quantification, nous introduisons une paire d’opérateurs de création et d’annihilation par mode \(l\) \textit{and} par fonction de base \(k_l\) définie comme:

\[\begin{split}\start{aligned} a_{k_l}^\dagger |\cdots, 0_1 \cdots 0_{k_l} \cdots 0_{N_l}, \cdots\rangle &= | \cdots, 0_1 \cdots 1_{k_l} \cdots 0_{N_l}, \cdots\rangle \\ a_{k_l}^\dagger | \cdots, 0_1 \cdots 1_{k_l} \cdots 0_{N_l}, \cdots\rangle &= 0 \\ a_{k_l} | \cdots, 0_1 \cdots 1_{k_l} \cdots 0_{N_l}, \cdots\rangle &= | \cdots, 0_1 \cdots 0_{k_l} \cdots 0_{N_l}, \cdots\rangle \\ a_{k_l} | \cdots, 0_1 \cdots 0_{k_l} \cdots 0_{N_l}, \cdots\rangle &= 0 \\ \end{aligned}\end{split}\]

avec

\[\begin{split}\start{aligned} \left[ a_{k_l}^\dagger, a_{h_m}^\dagger \right] &= 0 \\ \left[ a_{k_l}, a_{h_m} \right] &= 0 \\ \left[ a_{k_l}^\dagger, a_{h_m} \right] &= \delta_{l,m} \, , \delta_{k_l,h_m} \end{aligned}\end{split}\]

La deuxième forme de quantification est obtenue en exprimant le potentiel comme

\[\begin{split}\begin{aligned} \mathcal{H}_\text{vib}^{SQ} =& \sum_{l=1}^L \sum_{k_l,h_l}^{N_l} \langle \phi_{k_l} | T(Q_l) + V^{[l]}(Q_l) | \\phi_{h_l} \rangle a_{k_l}^+ a_{h_l} \\ +& \sum_{l<m}^L \sum_{k_l,h_l}^{N_l} \sum_{k_m,h_m}^{N_m} \langle \phi_{k_l} \phi_{k_m} | V^<unk> l, ]}(Q_l, Q_m) | \phi_{h_l} \phi_{h_m} \rangle a_{k_l}^+ a_{k_m}^+ a_{h_l} a_{h_m} + \cdots \end{aligned}\end{split}\]

Nous soulignons ici la différence entre les opérateurs définis ci-desus et la structure électronique. Tout d’abord, comme nous l’avons déjà mentionné, le potentiel contient (en principe) des termes de couplage à trois corps et plus, qui conduisent à des cordes avec six (ou plus) opérateurs de seconde quantification. De plus, l’Hamiltonien conserve le nombre de particules pour chaque mode, comme on peut le voir du fait que le nombre d’opérateurs de création et d’annihilation pour un mode donné est le même à chaque terme. Néanmoins, différents modes sont couplés par des termes à deux corps (et plus) contenant des opérateurs de SQ appartenant à différents modes :math:` l ` et :math:` m `.

Référence: Ollitrault, Pauline J., et al., arXiv:2003.12578 (2020).

Calcul du potentiel électronique

Résoudre l’ESE pour différentes configurations nucléaires afin d’obtenir la fonction PES :math:` V (Q_1, lpoints, Q_L) `. Jusqu’à présent, Qiskit donne la possibilité d’approcher le PSE avec un champ de force quartique.

\[V(Q_1, \ldots, Q_L) = \frac{1}{2} \sum_{ij} k_{ij} Q_i Q_j + \frac{1}{6} \sum_{ijk} k_{ijk} Q_i Q_j Q_k + \frac{1}{16} \sum_{ijkl} k_{ijkl} Q_i Q_j Q_k Q_l\]

L’avantage de cette forme pour le PSE est que les champs de force anharmonique (:math:` k_{ij}`, :math:` k_{ijk}`, :math:` k_{ijkl}`) peuvent être calculés par des approches de différences finies. Pour les méthodes pour lesquelles l’énergie nucléaire Hessian peut être calculée de manière analytique avec des méthodes théorique de réponse (telles que HF et DFT), le champ de force quartique peut être calculé par différenciation semi-numérique comme suit:

\[k_{ijk} = \frac{H_{ij}(+\delta Q_k) - H_{ij}(-\delta Q_k)}{2\delta Q_k}\]

et

\[k_{ijkl} = \frac{H_{ij}(+\delta Q_k+\delta Q_l) - H_{ij}(+\delta Q_k-\delta Q_l) -H_{ij}(-\delta Q_k+\delta Q_l) + H_{ij}(-\delta Q_k+\delta Q_l)} {4\delta Q_k \delta Q_l}\]

Une telle procédure numérique est mise en œuvre, par exemple, dans la suite Gaussienne de programmes.

En pratique, cela peut être fait avec Qiskit en utilisant le GaussianForceDriver.

from qiskit.chemistry.drivers import GaussianForcesDriver
# if you ran Gaussian elsewhere and already have the output file
driver = GaussianForcesDriver(logfile='aux_files/CO2_freq_B3LYP_ccpVDZ.log')

# if you want to run the Gaussian job from Qiskit
# driver = GaussianForcesDriver(
#                 ['#p B3LYP/6-31g Freq=(Anharm) Int=Ultrafine SCF=VeryTight',
#                  '',
#                  'CO2 geometry optimization B3LYP/6-31g',
#                  '',
#                  '0 1',
#                  'C  -0.848629  2.067624  0.160992',
#                  'O   0.098816  2.655801 -0.159738',
#                  'O  -1.796073  1.479446  0.481721',
#                  '',
#                  ''

Schématiser avec un qubit Hamiltonien

Maintenant que nous avons une approximation pour le potentiel, nous devons écrire le Hamiltonien en deuxième quantification. Pour cela, nous devons sélectionner une base modale pour calculer les intégraux d’un corps unique \(\langle\phi_{k_i}| V(Q_i) | \phi_{h_i} \rangle\), intégraux à deux corps \(\langle\phi_{k_i} \phi_{k_j}| V(Q_i,Q_j) | \\phi_{h_i} \phi_{h_j} \rangle\)…

Dans le cas le plus simple, les fonctions :math:` phi sont les fonctions propres d’oscillateur harmonique pour chaque mode. Le principal avantage de ce choix réside dans le fait que les intégrales d’un PES exprimé en tant qu’expansion de Taylor sont faciles à calculer avec cette base. Une routine pour le calcul de ces intégrales est implémentée dans Qiskit.

L’opérateur bosonique, :math:` mathcal{H}_text{vib}^{SQ}`, est alors créé et doit être associé à un opérateur de qubit. Le schéma direct introduit dans la première section de ce tutoriel peut être utilisé comme suit:

from qiskit.chemistry.transformations import (BosonicTransformation,
                                              BosonicTransformationType, # type of modal basis
                                              BosonicQubitMappingType) # type of boson to qubit mapping

bosonic_transformation = BosonicTransformation(qubit_mapping=BosonicQubitMappingType.DIRECT,
                                               transformation_type=BosonicTransformationType.HARMONIC,
                                               basis_size=2,
                                               truncation=2)

Dans la cellule précédente, nous avons défini une transformation bosonique afin d’exprimer le Hamiltonien dans la base modale harmonique, avec 2 modals par mode avec le potentiel tronqué à l’ordre 2 et le boson ‘direct’ à qubit cartographié. Le calcul est alors exécuté comme suit :

qubit_op, _ = bosonic_transformation.transform(driver)
print(qubit_op.print_details())

print('\nThe total number of modes in this system is {}'.format(bosonic_transformation.num_modes))

IIIIIIII        (5077.2365606250005+0j)
IIIIIIIZ        (-601.9000340625+0j)
IIIIIIZI        (-1810.6538965625004+0j)
IIIIIZII        (-342.04261812500005+0j)
IIIIXXII        (-84.41897986727395+0j)
IIIIYYII        (-84.41897986727395+0j)
IIIIZIII        (-1027.4430731250004+0j)
IIIZIIII        (-155.20426031250003+0j)
IIZIIIII        (-467.50585156250014+0j)
IZIIIIII        (-157.22939093750003+0j)
ZIIIIIII        (-472.28961593750006+0j)
IIIIIZIZ        (1.2599134375000005+0j)
IIIIIZZI        (3.779740312500002+0j)
IIIIXXIZ        (22.272716326623776+0j)
IIIIYYIZ        (22.272716326623776+0j)
IIIIXXZI        (66.81814897987134+0j)
IIIIYYZI        (66.81814897987134+0j)
IIIIZIIZ        (3.779740312500002+0j)
IIIIZIZI        (11.339220937500007+0j)
IIIZIIIZ        (-1.5962606250000009+0j)
IIIZIIZI        (-4.788781875000002+0j)
IIIZIZII        (-0.6414307812500002+0j)
IIIZXXII        (-5.411241592930711+0j)
IIIZYYII        (-5.411241592930711+0j)
IIIZZIII        (-1.9242923437500008+0j)
IIXXIIXX        (-0.5021409375000001+0j)
IIYYIIXX        (-0.5021409375000001+0j)
IIXXIIYY        (-0.5021409375000001+0j)
IIYYIIYY        (-0.5021409375000001+0j)
IIZIIIIZ        (-4.788781875000002+0j)
IIZIIIZI        (-14.366345625000008+0j)
IIZIIZII        (-1.924292343750001+0j)
IIZIXXII        (-16.233724778792137+0j)
IIZIYYII        (-16.233724778792137+0j)
IIZIZIII        (-5.772877031250004+0j)
IZIIIIIZ        (-1.1489604687500004+0j)
IZIIIIZI        (-3.4468814062500015+0j)
IZIIIZII        (-0.42099484375000024+0j)
IZIIXXII        (-1.6031887335286772+0j)
IZIIYYII        (-1.6031887335286772+0j)
IZIIZIII        (-1.2629845312500008+0j)
IZIZIIII        (-0.13775546875000005+0j)
IZZIIIII        (-0.4132664062500002+0j)
XXIIIIXX        (0.8980418750000001+0j)
YYIIIIXX        (0.8980418750000001+0j)
XXIIIIYY        (0.8980418750000001+0j)
YYIIIIYY        (0.8980418750000001+0j)
XXXXIIII        (1.986637812500001+0j)
YYXXIIII        (1.986637812500001+0j)
XXYYIIII        (1.986637812500001+0j)
YYYYIIII        (1.986637812500001+0j)
ZIIIIIIZ        (-3.4468814062500015+0j)
ZIIIIIZI        (-10.340644218750004+0j)
ZIIIIZII        (-1.2629845312500005+0j)
ZIIIXXII        (-4.809566200586032+0j)
ZIIIYYII        (-4.809566200586032+0j)
ZIIIZIII        (-3.7889535937500023+0j)
ZIIZIIII        (-0.4132664062500002+0j)
ZIZIIIII        (-1.2397992187500007+0j)


The total number of modes in this system is 4

Pour avoir un nombre différent de modaux par mode:

bosonic_transformation_2 = BosonicTransformation(qubit_mapping=BosonicQubitMappingType.DIRECT,
                                               transformation_type=BosonicTransformationType.HARMONIC,
                                               basis_size=[0,0,4,2],
                                               truncation=2)
qubit_op_2, _ = bosonic_transformation_2.transform(driver)
print(qubit_op_2.print_details())

IIIIII  (3309.8257831250007+0j)
IIIIIZ  (-164.15502593750006+0j)
IIIXIX  (1.8979435436158592+0j)
IIIYIY  (1.8979435436158592+0j)
IIIIZI  (-494.3581484375001+0j)
IIXIXI  (6.378706032974414+0j)
IIYIYI  (6.378706032974414+0j)
IIIZII  (-828.3474121875001+0j)
IIZIII  (-1166.1228171875002+0j)
IZIIII  (-161.85614656250007+0j)
ZIIIII  (-486.16988281250013+0j)
IZIIIZ  (-0.13775546875000005+0j)
IZIXIX  (0.19481565219731314+0j)
IZIYIY  (0.19481565219731314+0j)
IZIIZI  (-0.4132664062500002+0j)
IZXIXI  (0.3374306077154137+0j)
IZYIYI  (0.3374306077154137+0j)
IZIZII  (-0.6887773437500002+0j)
IZZIII  (-0.9642882812500004+0j)
XXIIXX  (1.986637812500001+0j)
YYIIXX  (1.986637812500001+0j)
XXIIYY  (1.986637812500001+0j)
YYIIYY  (1.986637812500001+0j)
XXXIIX  (0.8288402231450827+0j)
YYXIIX  (0.8288402231450827+0j)
XXYIIY  (0.8288402231450827+0j)
YYYIIY  (0.8288402231450827+0j)
XXIXXI  (4.245123515804727+0j)
YYIXXI  (4.245123515804727+0j)
XXIYYI  (4.245123515804727+0j)
YYIYYI  (4.245123515804727+0j)
XXXXII  (6.957428881323849+0j)
YYXXII  (6.957428881323849+0j)
XXYYII  (6.957428881323849+0j)
YYYYII  (6.957428881323849+0j)
ZIIIIZ  (-0.4132664062500002+0j)
ZIIXIX  (0.5844469565919395+0j)
ZIIYIY  (0.5844469565919395+0j)
ZIIIZI  (-1.2397992187500007+0j)
ZIXIXI  (1.0122918231462412+0j)
ZIYIYI  (1.0122918231462412+0j)
ZIIZII  (-2.066332031250001+0j)
ZIZIII  (-2.8928648437500013+0j)

Maintenant que l’Hamiltonien est prêt, il peut être utilisé dans un algorithme quantique pour trouver des informations sur la structure électronique vibronique de la molécule correspondante. Consultez nos tutoriels sur le calcul de l’état du sol (Ground State Calculation) et le calcul des états excités (Excited State Calculation) pour en savoir davantage sur la façon de le faire dans Qiskit !

import qiskit.tools.jupyter
%qiskit_version_table
%qiskit_copyright

Version Information

Qiskit Software	Version
Qiskit	0.23.0
Terra	0.16.0
Aer	0.7.0
Ignis	0.5.0
Aqua	0.8.0
IBM Q Provider	0.11.0
System information
Python	3.7.4 (default, Aug 13 2019, 15:17:50) [Clang 4.0.1 (tags/RELEASE_401/final)]
OS	Darwin
CPUs	6
Memory (Gb)	16.0
Wed Oct 21 12:22:51 2020 CEST

This code is a part of Qiskit

© Copyright IBM 2017, 2020.

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.