A general suite for benchmarking nonadiabatic dynamics

In this tutorial we provide a general suite for benchmarking nonadiabatic dynamics described in ref ¹ together with data (this link is the collection) of all figures. Among all the dynamics tested, the robustness of nonadiabatic field (NaF)²³ is highlighted.

For a full publication list of the generalized constraint coordinate-momentum phase space (CPS) and NaF please see this page.

For an illustration of the nonadiabatic force on the surface hopping dynamics please see this page (reported in ref²).

1. LVCM (Linear vibronic coupling models)

Figures_and_data

The Hamiltonian of a LVCM (linear vibronic coupling model) reads

where

Note that the relation between the dimensionless weighted normal-mode coordinate/momentum, $\{{\overline{R}}_k,{\overline{P}}_k\}$$\{\bar{R}_k,\bar{P}_k\}$, and the corresponding canonical (mass-weighted) coordinate/momentum, $\{R_k,P_k\}$$\{R_k,P_k\}$, in the diabatic representation reads

The parameters of LVCMs in the suite are listed as follows, except the parameters of the 7-state 39-mode LVCM of Thymine in ref ⁴ which are kindly provided by Martha Yaghoubi Jouybari and Fabrizio Santoro.

Table 1. Parameters of 2-state 3-mode LVCM of Pyrazine (unit: eV) ⁵

Table 2. Parameters of 24-mode LVCM of Pyrazine (unit: eV)⁶

In our simulation of the pyrazine models, the second diabatic (electronic) state is initially occupied, and the nuclear initial condition is set to be the Wigner distribution of the vibrational ground state,

The numerically exact results in our data are produced by MCTDH⁷.

Table 3. Parameters of 3-state 2-mode LVCM of $\mathrm{Cr}(\mathrm{CO}{)}_5$$\mathrm{Cr(CO)_5}$ (unit: eV)⁸

In our test of the $\mathrm{Cr}(\mathrm{CO}{)}_5$$\mathrm{Cr(CO)_5}$ model, the second diabatic (electronic) state is initially occupied, and the nuclear initial condition is set to be the Wigner distribution of the following Wigner distribution,

where $r_1=0,r_2=14.3514,{\alpha }_1=0.4501$$r_1 = 0, r_2=14.3514, \alpha_1 = 0.4501$, and ${\alpha }_2=0.4586$$\alpha_2=0.4586$. The numerically exact results in our data are produced by MCTDH, taken from ref ⁸.

In our test of the Thymine model, the $\pi {\pi }^{\ast }2$$\pi\pi^*2$ state is initially occupied with each nuclear mode in its vibrational ground state, i.e., eq (5). The numerically exact results in our data are produced by ML-MCTDH, taken from ref ⁴.

2. Tully models

Figures_and_data

The single avoided crossing (SAC) model⁹ reads

where $A=0.01,B=1.6,C=0.005,D=1$$A=0.01, B=1.6, C=0.005, D=1$ (in a.u.). The dual avoided crossing (DAC) model⁹ reads

with $A=0.01,B=1.6,C=0.005,D=1.0$$A=0.01, B=1.6, C=0.005, D=1.0$ (in a.u.). The extended coupling region (ECR) model ⁹ reads

where $E_0=-0.0006,B=0.9,C=0.1$$E_0 = -0.0006, B=0.9, C=0.1$ (in a.u.), and $h(R)$$h(R)$ denotes the Heaviside step function. In our data, the nuclear DOF is with mass $m=2000$$m=2000$ (in a.u.). The electronic ground state is selected as the initial state, while the initial condition for the nuclear DOF follows the Wigner distribution

with width parameter $\alpha =1$$\alpha=1$. The coordinate center $R_0$$R_0$ for SAC, DAC, and ECR models is set to -3.8, -10.0, and -13.0, respectively. The initial momentum $P_0$$P_0$ is scanned and is used as the x-axis in our figures.

The exact results in our data are provided by DVR¹⁰.

3. 3-state photodissociation models

Figures_and_data

The potential matrix elements of anharmonic 3-state photodissociation models¹¹ reads

The parameters are listed as follows,

Table 4. Parameters of the 3-state photodissociation models (unit: a.u.)¹¹

The initial electronic configuration is set to the first diabatic state, and the nuclear DOF is sampled from the Wigner distribution of the ground state,

where $m=20000$$m=20000$ a.u. is the mass of the nuclear DOF, $\omega =0.005$$\omega=0.005$ a.u. is the vibrational frequency of the ground state, and the center $R_e$$R_e$ is listed in Table 4.

The exact results in our data are produced by DVR¹⁰.

4. 2-state photodissociation model

Figures_and_data

The potential matrix elements of the anharmonic 2-state photodissociation model¹¹¹ reads

The parameters are listed as follows,

Table 5. Parameters of the 2-state photodissociation model (unit: a.u.)¹¹¹

The initial condition and mass are identical to those of Model 2 of the 3-state photodissociation models.

The exact results in our data are produced by DVR¹⁰.

5. Spin-boson models

Figures_and_data

The Hamiltonian of a spin-boson model reads¹²

We use the discretization scheme for the Ohmic spectral density¹² with the Kondo parameter $\alpha$$\alpha$ and the cut-off frequency ${\omega }_c$$\omega_c$, which results in the following expressions¹³¹⁴ with $N_b$$N_b$ discrete modes:

We employ four specific spin-boson models with $\epsilon =\mathrm{\Delta }=1$$\varepsilon=\Delta=1$ at $\beta =5$$\beta = 5$. The parameters¹⁵ of the bath DOFs are $\alpha \in \{0.1,0.4\}$$\alpha\in\{0.1,0.4\}$ and ${\omega }_c\in \{1,2.5\}$$\omega_c\in\{1,2.5\}$ with $N_b=300$$N_b=300$.

The system occupies the first diabatic electronic state, and the nuclear DOFs are sampled from their equilibrium Wigner density. The numerically exact results in our data¹⁵ are calculated by eHEOM¹⁶.

6. Site-exciton models

Figures_and_data

The Hamiltonian of site-exciton models reads¹⁷

where ${\hat{H}}_b=\sum _{n,j}\frac{1}{2}({\hat{P}}_{nj}^2+{\omega }_{nj}^2{\hat{R}}_{nj}^2)$$\hat{H}_{b}=\sum_{n, j} \frac{1}{2}\left(\hat{P}_{n j}^{2}+\omega_{n j}^{2} \hat{R}_{n j}^{2}\right)$ and ${\hat{H}}_{s-b}=\sum _{n,j}{c}_{nj}{\hat{R}}_{nj}|n⟩⟨n|$$\hat{H}_{s-b}=\sum_{n, j} c_{n j} \hat{R}_{nj}|n\rangle\langle n|$. We employ the Debye spectral density

for each state, where $\lambda$$\lambda$ and ${\omega }_c$$\omega_c$ denote the bath reorganization energy and the characteristic frequency, respectively. The corresponding discretization scheme¹⁸¹⁹²⁰ is

For the 7-state FMO model²¹, the system Hamiltonian reads

The bath reorganization energy is $\lambda =35{\mathrm{cm}}^{-1}$$\lambda=35 \mathrm{~cm}^{-1}$, the characteristic frequency is ${\omega }_c=106.14{\mathrm{cm}}^{-1}$$\omega_{c}=106.14 \mathrm{~cm}^{-1}$, and $N_b=50$$N_b = 50$ for each state.

For 3-state singlet-fission model of pentacene²²²³, the system Hamiltonian reads

The bath reorganization energy is $\lambda =0.1$$\lambda=0.1$ eV, the characteristic frequency is ${\omega }_c=0.18$$\omega_c=0.18$ eV, and $N_b=200$$N_b = 200$ for each state.

In these two site-exciton models, the bath DOFs are sampled from their equilibrium Wigner distribution, and the first diabatic electronic state is initially occupied.

The numerically exact results in our data¹ are produced by TD-DMRG²⁴ in our FMO model, and by HEOM²⁵ in our SF model.

7. Atom-in-cavity model

Figures_and_data

The total Hamiltonian for the atom-in-cavity models²⁶²⁷ can be decomposed into three parts

where

The quantity ${\mu }_{nm}={\mu }_{mn}$$\mu_{nm}=\mu_{mn}$ denotes the transition dipole moment between the $n$$n$-th and $m$$m$-th atomic energy levels. The coupling strength between the atom and the $j$$j$-th optical field mode is given by

where $L=236200$$L=236200$ a.u. represents the cavity length, $r_0=L/2$$r_0=L/2$ corresponds to the cavity center, and ${\epsilon }_0$$\varepsilon_0$ is the vacuum permittivity. The mode frequency is determined by ${\omega }_j=j\pi c/L$$\omega_j=j\pi c/L$, with $c=137.036\text{ a.u.}$$c=137.036 \textrm{\ a.u.}$ to be the light speed in vacuum. We employ 400 standing wave modes and model the hydrogen atom system as a three-level system. The energy levels and transition dipole moments (in atomic units) are: ${\epsilon }_1=-0.6738$$\varepsilon_1=-0.6738$, ${\epsilon }_2=-0.2798$$\varepsilon_2=-0.2798$, ${\epsilon }_3=-0.1547$$\varepsilon_3=-0.1547$, ${\mu }_{12}=-1.034$$\mu_{12}=-1.034$, ${\mu }_{23}=-2.536$$\mu_{23}=-2.536$, and ${\mu }_{13}=0$$\mu_{13}=0$. A reduced two-level system is also considered, where only the lowest two energy levels of the hydrogen atom is retained.

Initially the highest atomic state is occupied, with each optical field mode in their optical vacuum state with Wigner distribution reading

The exact results are produced by truncated configuration interaction and are taken from refs²⁶²⁷.

8. Holstein model

Figures_and_data

The Hamiltonian of the 1-dim Holstein model²⁸²⁹ is given by

where the electronic system is described by a nearest-neighbor tight-binding model:

Here, ${\hat{a}}_n^{†}$$\hat a_n^\dagger$ and ${\hat{a}}_n$$\hat a_n$ denote the fermionic creation and annihilation operators at the $n$$n$-th lattice site, while $E_n$$E_n$ and $V$$V$ represent the on-site energy and the coupling between neighboring sites, respectively. For a periodic chain with a finite number of lattice sites $L$$L$, one typically sets $E_n=0$$E_n=0$ and assumes ${\hat{a}}_{n+L}={\hat{a}}_n$$\hat a_{n+L}=\hat a_n$.

The phonon Hamiltonian describes a bosonic bath environment:

The term ${\hat{H}}_c$$\hat H_c$ encompasses all electron-phonon coupling effects:

Considering the low carrier concentration regime, where the system contains only a single electron. In this case, the fermionic operators act within the single-particle Fock space, and the system Hamiltonian becomes isomorphic to a multistate Hamiltonian³⁰. Consequently, the fermionic operators can be expressed in terms of state representations:

or equivalently,

We set $V=0.083\mathrm{eV},$$V = 0.083 \mathrm{\ eV},$$d=R=7.2\mathrm{Å}$$d=R=7.2\mathrm{\ Å}$ and $L=21$$L=21$³¹³², where the unmentioned parameters are used in the calculation of correlation function. Other parameters are listed as follows.

Table 6. Parameters of the Holstein model in refs³¹³²

We guide readers to Section S5-S6 of Supporting Information of ref¹ for details of initial condition and correlation function. The exact results in our data are taken from ref³¹, produced by TD-DMRG.

9. References