# Nonadiabatic to Adiabatic Transition of Electron Transfer in Colloidal
Quantum Dot Molecules

Bokang Hou  Department of Chemistry, University of California, Berkeley,
California 94720, United States;  Michael Thoss  Institute of Physics, Albert-
Ludwigs University Freiburg, Hermann-Herder-Straße 3, 79104 Freiburg, Germany;
Uri Banin  Institute of Chemistry and the Center for Nanoscience and
Nanotechnology, The Hebrew University of Jerusalem, 91904 Jerusalem, Israel
Eran Rabani  Department of Chemistry, University of California, Berkeley,
California 94720, United States;  Materials Sciences Division, Lawrence
Berkeley National Laboratory, Berkeley, California 94720, United States;  The
Raymond and Beverly Sackler Center of Computational Molecular and Materials
Science, Tel Aviv University, Tel Aviv 69978, Israel;

###### Abstract

Electron transfer is an important and fundamental process in chemistry,
biology and physics, and has received significant attention in recent years.
Perhaps one of the most intriguing questions concerns with the realization of
the transitions between nonadiabatic and adiabatic regimes of electron
transfer, as the coupling (hybridization) energy, $J$, between the donor and
acceptor is varied. Here, using colloidal quantum dot molecules, a new class
of coupled quantum dot dimers, we computationally demonstrate how the
hybridization energy between the donor and acceptor quantum dots can be tuned
by simply changing the neck dimensions and/or the quantum dot size. This
provides a handle to tune the electron transfer from the nonadiabatic over-
damped Marcus regime to the coherent adiabatic regime in a single system,
without changing the reorganization energy, $\lambda$, or the typical phonon
frequency, $\omega_{c}$. We develop an atomistic model to account for several
donor and acceptor states and how they couple to the lattice vibrations, and
utilize the Ehrenfest mean-field mixed quantum-classical method to describe
the charge transfer dynamics as the nonadiabatic parameter, $\gamma$, is
varied. We find that charge transfer rates increase by several orders of
magnitude as the system is driven to the coherent, adiabatic limit, even at
elevated temperatures, and delineate the inter-dot and torsional acoustic
modes that couple most strongly to the charge transfer reaction coordinate.

The theoretical framework to describe charge transfer reactions between a
donor and an acceptor dates back to the seminal work of Marcus.[1] Using a
semi-classical perturbative approach, Marcus derived a relation for the outer-
shell electron transfer rate at high temperature ($T$) in terms of the driving
force ($\Delta\varepsilon$), the reorganization energy ($\lambda$)
characterizing the coupling to nuclear fluctuations, and the
hybridization/electronic couplings between donor and acceptor states ($J$),
assumed to be small (see Fig. 1(a) for a sketch of the relevant energy
scales):

$k_{\rm M}=\frac{\left|J\right|^{2}}{\hbar}\sqrt{\frac{\pi}{\lambda k_{\rm
B}T}}\exp{\frac{-\left(\Delta\varepsilon+\lambda\right)^{2}}{4\lambda k_{\rm
B}T}}.$ (1)

Marcus theory is suitable for the nonadiabatic electron transfer regime, i.e.
for $\gamma\ll 1$, where $\gamma$ is the adiabatic parameter [2] defined by
the ratio of the donor-acceptor coupling, the characteristic phonon bath
frequency ($\omega_{c}$), and the reorganization energy

$\gamma=\frac{|J|^{2}}{2\hbar\omega_{c}}\sqrt{\frac{\pi}{\lambda k_{\rm
B}T}}.$ (2)

Figure 1: (a) Schematic sketch of the potential energy surfaces along the
reaction coordinate [1] $Q$ and the energy scales appearing in Eq. (1). The
solid and dashed lines represent the adiabatic and diabatic potential energy
surfaces, respectively. (b) Hyper-sphere plots of donor (left NC) and acceptor
states (right NC) obtained from the semi-empirical pseudopotential calculation
and FB-localization of the quasi-electron eigenstates in wurtzite CQD dimers.
The energy difference between $1S_{\rm e}$-like ($\ket{D_{1}}$/$\ket{A_{1}}$)
and $1P_{\rm e}$-like ($\ket{D_{2}}$/$\ket{A_{2}}$) orbitals is in the range
of $60\sim 100$ meV, depending on the size of the NC building blocks.

Marcus nonadiabatic electron transfer theory was extended in several different
directions. Jortner and coworkers described the role of quantum nuclear
fluctuations as well as non-parabolicities in the donor and acceptor free
energies on the electron transfer rate.[3] Redfield theory [4] was used to
extend Marcus theory to account for coherences between the donor and acceptor,
observed in photosynthetic complexes and quantum dot islands exhibiting strong
electronic couplings.[5, 6] Unlike the Marcus regime in which the transfer
dynamics decay exponentially and can be characterized by a rate constant, the
population dynamics in this coherent limit oscillate between the donor and
acceptor, eventually relaxing to equilibrium. Zusman developed a framework to
describe the crossover from the Marcus weak hybridization limit to the
adiabatic limit.[7, 8] In the adiabatic limit, the pre-exponential factor
appearing in Marcus’ rate theory (c.f., Eq. (1)) becomes independent of the
electron tunneling probability between the donor and the acceptor, and the
rate constant is given by transition state rate theory:

$k_{\rm TST}=\frac{\omega_{c}}{4}\sqrt{\frac{\lambda}{\pi k_{\rm
B}T}}\exp{\frac{-\left(\Delta\varepsilon+\lambda\right)^{2}}{4\lambda k_{\rm
B}T}}.$ (3)

Experimental manifestation of transition from nonadiabatic to adiabatic charge
transfer is challenging and requires exquisite control over electronic and
vibrational degrees of freedom.[9, 10] Recently, Zhu et al.[11] studied
electron transfer reactions in mixed valance donor-bridge-acceptor complexes.
By changing the length of the bridge connecting the donor and acceptor and the
functional groups on both the donor and the acceptor, they were able to drive
the system from the adiabatic to the nonadiabatic regime. However, their
approach was limited to a either $1$, $2$, or $3$ bridge units, covering a
narrow range of behaviors.

In this work, we revisit this problem and consider the charge transfer between
two coupled colloidal quantum dot (CQD) nanocrystals (NCs) that are connected
by a neck/bridge (see Fig. 1(b)). Charge transfer in such systems is
particularly interesting due to the flexibility in designing the donor and
acceptor states by, for example, changing the width of the neck
($D_{\text{neck}}$) between two NCs and/or the diameter of each NC cores
($D_{\text{core}}$), as well as controlling the shell thickness
($D_{\text{shell}}$).[12, 13] By continuously varying these parameters, the
hybridization between the donor and the acceptor, $J$, can be tuned across a
wide range of values while at the same time the reorganization energy,
$\lambda$, and the typical vibrations frequency, $\omega_{c}$, change
slightly. Furthermore, the hybridization energies can be tuned to be larger
than the thermal energy at room temperature. This is quite distinct from the
behavior of molecular junction, where a change in the hybridization energy is
often accompanied by a change to the reorganization energy and the vibrational
frequency. This is because control over these parameters in molecular
junctions is achieved by either changing the donor/acceptor molecules, or by
extending the bridge connecting them, which result in changes in the other
parameters as well.

From a theoretical/computational perspective, this system poses several
challenges, particularly with respect to the dimensions and number of valance
electrons. Therefore, our approach to describe charge transfer is based on a
model Hamiltonian which is parameterized by first-principle and semiempirical
calculations. The total Hamiltonian can be divided into a sum of three terms,
$H=H_{S}+H_{B}+H_{SB}$, where $H_{S}$ describes the electronic system (donor
and acceptor states and their hybridization), $H_{B}$ is the Hamiltonian for
the nuclear degrees of freedom (DOF) of the QD dimer (nuclear vibrations), and
the interaction between the electronic system and the nuclear vibrations,
approximated to the lowest order in the nuclear DOF,[14] is described by
$H_{SB}$. The three terms are given by:

$\displaystyle
H_{S}=\sum_{n\in\mathcal{D},\mathcal{A}}\varepsilon_{n}\ket{\phi_{n}}\bra{\phi_{n}}+\sum_{\begin{subarray}{c}n\in\mathcal{D}\\\
m\in\mathcal{A}\end{subarray}}J_{nm}\ket{\phi_{n}}\bra{\phi_{m}}+h.c.$ (4)
$\displaystyle
H_{B}=\sum_{\alpha}\frac{1}{2}{P}_{\alpha}^{2}+U(Q_{1},Q_{2},\dots)$
$\displaystyle
H_{SB}=\sum_{\alpha}\sum_{\begin{subarray}{c}n,m\in\mathcal{D},\mathcal{A}\end{subarray}}\ket{\phi_{n}}\bra{\phi_{m}}V_{nm}^{\alpha}{Q}_{\alpha}.$

In the above equations, $\varepsilon_{n}$ is the energy of state
$\ket{\phi_{n}}$ ($n\in\mathcal{D},\mathcal{A}$ with
$\mathcal{D}=\\{D_{1},D_{2}\dots\\}$, $\mathcal{A}=\\{A_{1},A_{2}\dots\\}$)
and $J_{nm}$ is the hybridization term between the donor state
$\ket{\phi_{n}}$ and acceptor state $\ket{\phi_{m}}$. To obtain the donor and
acceptor states, we used the semi-empirical pseudopotential model [15, 16] to
describe the quasi-electron Hamiltonian ($\hat{h}_{\rm QP}$) and the filter
diagonalization technique [17, 18] to calculate the eigenstates
($\ket{\psi_{n}}$) of the dimer near the bottom of the conduction band. The
Förster-Boys localization scheme [19, 20] was then used to generated the local
donor and acceptor states ($\ket{\phi_{n}}$) from the eigenstates
($\ket{\psi_{n}}$), yielding $\varepsilon_{n}=\bra{\phi_{n}}\hat{h}_{\rm
QP}\ket{\phi_{n}}$ and $J_{nm}=\bra{\phi_{n}}\hat{h}_{\rm QP}\ket{\phi_{m}}$
for $n\in\mathcal{D},m\in\mathcal{A}$, and otherwise by construction, it is
set to $J_{nm}=0$. We use the Stilinger-Weber [21] potential energy surface to
describe the nuclear degrees of freedom, where $P_{\alpha}$ and $Q_{\alpha}$
are the $\alpha$ mass-weighted vibrational normal mode coordinates, determined
by diagonlizing the Hessian matrix at the equilibrium geometry.[14] The
strength of coupling between states $\ket{\phi_{n}}$ and $\ket{\phi_{m}}$ to
mode $\alpha$ is given by $V_{nm}^{\alpha}$, and is determined directly using
the pseudopotential Hamiltonian[14] (see Supporting Information Eq. (6) for
more information).

Figure 2: Hybridization energy $J$ of the ground donor and acceptor states as
a function of (a) $D_{\text{neck}}$ and (b) $D_{\text{shell}}$ for [100] and
[001] orientation attachments. Solid, dashed, and dotted-dashed lines
correspond to $D_{\text{core}}=2.2$nm, $D_{\text{core}}=3.0$nm,
$D_{\text{core}}=3.9$nm. (c) Adiabatic parameter $\gamma$ vs. $J/\lambda$ for
varying neck-core sizes or shell-core sizes at 300 K. We define $\gamma=1$ as
the crossover from the nonadiabatic to adiabatic electron transfer regimes.
Figure 3: Population dynamics of the donor, $p_{\mathcal{D}}(t)$ (blue), and
the acceptor, $p_{\mathcal{A}}(t)$ (red), at $T=300$K for three regimes: (a)
nonadiabatic, $D_{\rm neck}=2.4$ nm, $D_{\rm core}=2.2$ nm (b) intermediate,
$D_{\rm neck}=3.1$ nm, $D_{\rm core}=3.0$ nm (c) adiabatic, $D_{\rm neck}=3.8$
nm, $D_{\rm core}=3.9$ nm. All structures are in symmetric attachment, leading
to $p_{\mathcal{D}}=p_{\mathcal{A}}=\frac{1}{2}$ in the long-time limit. The
dashed lines in (a) are the population generated from master equation with
rates computed using Marcus theory.

Several low lying donor and acceptor states for the attachment orientation
$[100]$ are shown in Fig. 1(b). The $[100]$ attachment results in a symmetric
distribution of the charge density, similar to a homo-nuclear diatomic
molecule, with an ”atomic-like” basis of an effective mass particle-in-a-
sphere The two lowest energy donor/acceptor states are mainly comprised of
$1S_{\rm e}$-like orbitals, while higher lying states show $1P_{\rm e}$-like
character (either along or perpendicular to the dimer axis). For the systems
considered in this work, the energy gap between $1S_{\rm e}$-like and $1P_{\rm
e}$-like orbitals ranges from $60$ to $100$ meV. The symmetry is broken along
the $[001]$ orientation attachment, resulting in a small energy bias and an
asymmetric charge distribution, consistent with the behavior of hetero-nuclear
diatomic molecules (see Fig. S1 for more information).

In Fig. 2 panels (a) and (b) we plot the hybridization energy between the
ground donor and acceptor states, $J=J_{D_{1}A_{1}}$, as a function of the
neck ($D_{\rm neck}$) and shell ($D_{\rm shell}$) widths, respectively, with
different NC core diameters ($D_{\rm core}$). As expected, increasing the neck
width and core diameters or decreasing the shell thickness, results in an
exponential increase of the magnitude of the hybridization energy.
Furthermore, the $[100]$ attachment (blue curves) shows a much steeper
dependence on $D_{\rm neck}$ and $D_{\rm shell}$, as a result of the larger
core-to-core distance for this orientation. In Fig. 2 (c) we plot the
adiabatic parameter, $\gamma$, as a function of the hybridization energy for
the same set of neck and shell dimensions shown in panels (a) and (b). We find
a crossover from the nonadiabatic to the adiabatic electron transfer regimes
as $J\rightarrow\lambda$. Since the characteristic vibrational frequency and
the reorganization energy depend weakly on the dimer geometry, we find that
the crossover mainly is affected by the hybridization energy. Note that for
this set of building block monomers, the electron transfer in the $[001]$
orientation attachment is in the adiabatic limit regardless of the neck sizes
or shell thickness, while it can be tuned from the nonadiabatic to the
adiabatic limits for the $[100]$ orientation attachment. The dimensions of CQD
dimers in Fig. 2 are summarized in the Supporting Information, tables S1 and
S2.

Figure 4: (a) Diagonal spectral densities $S_{D_{1}D_{1}}(\omega)$ and (b)
off-diagonal spectral densities $S_{D_{1}A_{1}}(\omega)$ for the same systems
studied in Fig. 3. The blue, orange and green lines correspond to Fig. 3 panel
(a), (b) and (c). The inset shows the phonon density of state (PDOS). The two
most important modes for contributing to dephasing are the acoustic inter-dot
vibrational mode ($\omega_{1}=0.13$ THz) shown in panel (c) and the torsional
($\omega_{2}=0.58$ THz) mode shown in panel (d). The arrows indicate the
motion of each atom for each mode.

In Fig. 3 we plots the donor and acceptor populations as a function of time in
the nonadiabatic-Marcus regime ($\gamma\ll 1$), the intermediate regime
($\gamma\approx 1$), and the adiabatic regime ($\gamma\gg 1$). We used the
Ehrenfest mean-field mixed quantum-classical method [22, 23, 24] to describe
the dynamics in all three regimes. As discussed below (c.f., Fig. 4), the
electron is mainly coupled to the low frequency acoustic modes, for which the
classical limit is adequate ($\hbar\omega_{c}\ll k_{\rm B}T$). For $\gamma\ll
1$ in Fig. 3 (a) we also compare the mean-field results to a master equation
(due to the presence of multiple donor and acceptor states) with transition
rates obtained from Marcus theory. The donor ($p_{\mathcal{D}}(t)$) and
acceptor ($p_{\mathcal{A}}(t)$) populations are given as projections onto the
donor and acceptor Hilbert spaces, respectively. Individual state populations
corresponding to the results shown in Fig. 3 are shown in the Supporting
Information Fig. S2 (population dynamics for all other structures studied in
this work are shown in Fig. S3 to S6).

For weak hybridization between the donor and acceptor (i.e., for small $D_{\rm
neck}$ and/or $D_{\rm shell}$), the population dynamics are characterized by
an over-damped exponential decay shown in Fig. 3 (a), with a decay rate that
approximately matches the Marcus rate between the ground donor and ground
acceptor states. For intermediate values of the shell and neck thicknesses,
i.e., for $\gamma\approx 1$, the population dynamics show underdamped coherent
oscillations, with a Rabi frequency that matches predominately the ground
state donor-acceptor hybridization (Fig. 3(b)). As the neck thickness is
further increased and/or the shell thickness is further decreased, the Rabi
period shortens, signifying the increase in the hybridization energy, as shown
in Fig. 3 (c).

The agreement between the mean field theory and the perturbative Marcus result
for $\gamma\ll 1$ is consistent with the so called ”average classical
limit”,[25] where the dynamics of the nuclear degrees of freedom are governed
by the arithmetic average donor/acceptor Hamiltonian.[26] The average
classical limit has been motivated by the analysis of the Wigner form of the
quantum mechanical expression for the relevant time-correlation function,[27]
suggesting that the average Hamiltonian (similar to mean field) provides the
most accurate approximation to the fully quantum mechanical results.[26] For
large values of $\gamma$, where the hybridization energy is much larger than
the reorganization energy ($J\gg\lambda$) and the system is in the weak
electron-phonon coupling limit, the average mean force on the nuclei is
similar to the force for each diabatic potential energy surfaces, and the mean
field dynamics accurately reproduce the many-body solution.111Quantum
mechanical test calculations using the multiconfiguration time-dependent
Hartree (MCTDH) method show good agreement with the results of the Ehrenfest
method.

Figure 5: (a) Electron transfer rates, $k$, plotted as a function of
$J/\lambda$ at $T=300$K for all structures studies in this work. The solid
black curve is a fit to a connection formula, $k^{-1}=k_{\rm M}^{-1}+k_{\rm
dp}^{-1}$ with an empirical form for $k_{\rm dp}\propto\sqrt{J}$ shown as
dashed line. We also show the $J\rightarrow 0$ Marcus limit, $k_{\rm M}^{-1}$.
(b) and (c): Electron transfer rates plotted as a function of the neck and
shell thicknesses, respectively. Different core sizes are shown in solid
($D_{\rm core}=2.2$nm), dashed ($D_{\rm core}=3.0$nm), and dot-dashed ($D_{\rm
core}=3.9$nm) lines.

We also find that as the dephasing rate increases with increasing
hybridization energy, $J$, deep in the adiabatic limit ($\gamma\gg 1$). This
seems to be opposite to the behavior expected for the spin-boson model,[29]
and can be traced to the inclusion of off-diagonal electron-phonon couplings
terms in our model Hamiltonian ($V_{nm}^{\alpha}$). As shown in the Supporting
Information Fig. S7, the behavior of the population dynamics and the dephasing
rates are consistent with the spin-boson model when the off-diagonal coupling
terms are turned off ($V_{n\neq m}^{\alpha}=0$). In addition, the presence of
several donor and several acceptor states also affects the dephasing rates,
particularly for $\gamma\gg 1$ where the $1P_{e}$-like donor/acceptor orbitals
play a significant role in the charge transfer dynamics due to a smaller
energy gap between $1S_{e}$ and $1P_{e}$-like states.

To further analyze the dephasing dynamics and delineate the modes that most
strongly affect the electron transfer rates, we define the spectral density
characterizing the electron-phonon interaction between states $\ket{\phi_{n}}$
and $\ket{\phi_{m}}$ ($n,m\in\mathcal{D},{A}$) as

$S_{nm}(\omega)=\pi\sum_{\alpha}\omega_{\alpha}\lambda^{\alpha}_{nm}\delta(\omega-\omega_{\alpha})$
(5)

where
$\lambda^{\alpha}_{nm}=\frac{1}{2}\left(\frac{V_{nm}^{\alpha}}{\omega_{\alpha}}\right)^{2}$
is the reorganization energy for mode $\alpha$. In Fig. 4 panels (a) and (b)
we plot the diagonal and off-diagonal spectral densities between the ground
states of the donor and acceptor ($\ket{\phi_{D_{1}}}$ and
$\ket{\phi_{A_{1}}}$) for the same systems shown in Fig. 3. The spectral
densities are very structured with stronger coupling to the low-frequency
acoustic modes ($\omega/2\pi<1.5$ THz) and weaker coupling to the high-
frequency optical modes ($\omega/2\pi>4$ THz). We find that the overall
magnitude of the diagonal spectral densities (panel (a)) decrease with
increasing $\gamma$ while the off-diagonal spectral densities behave the
opposite. Thus, setting the non-diagonal coupling terms to zero leads to a
decrease in the dephasing rates with increasing hybridization energies, which
is not the case observed in Fig. 3, where the dephasing times are governed by
the off-diagonal couplings. In addition, we identify the modes that contribute
the most to the dephasing dynamics (panels (c) and (d) of Fig. 4). The lower
frequency mode, $\omega_{1}$, involves a breathing motion while the higher
frequency mode, $\omega_{2}$ involves a torsional motion. These vibrational
modes are important to localize the charge and to facilitate the transfer
dynamics.

In Fig. 5 we provide a summary of the electron transfer rate constants across
the nonadiabatic to adiabatic transition, and show how they depend on the neck
and shell dimensions. In the adiabatic regime, we estimate the electron
transfer rate constant by fitting the decay of the envelop of the donor
population (see SI). Panel (a) of Fig. 5 shows a scattered plot of the rate
constants calculated for all the dimers considered in this work as a function
of the ground state hybridization energy, $J$, at room temperature. As can be
seen, the rate constants follows the $J^{2}$-Marcus theory dependence as
$J\rightarrow 0$. For larger values of $J$, the population dynamics are no
longer characterized by an over-damped exponential decay (unlike the Zusman
limit), and the electron transfer is dominated by the dephasing times, which
show a weaker dependence on $J$. The first order dependence of the dephasing
rate on $J$ can be derived for the simple two-level spin-boson model, and is
given by[30] $k_{\rm dp}=2S(2J)\coth{(J/k_{\rm B}T)}$, where $S(\omega)$ is
the spectral density defined above. For the multi-state model considered in
this work, the dephasing rate depends on other factors discussed previously,
such as the off-diagonal spectral densities and contributions from higher-
lying states, and increases with $J$ rather than decrease for the standard
two-state model. The solid black line in Fig. 5(a) is a fit to a connection
formula for the total electron transfer rate, given by $k^{-1}=k_{\rm
M}^{-1}+k_{\rm dp}^{-1}$ with an empirical form for $k_{\rm
dp}\propto\sqrt{J}$.

Fig. 5 panels (b) and (c) show the dependence of $k$ on the neck width,
$D_{\rm neck}$, and shell thickness, $D_{\rm shell}$, for different core
diameters $D_{\rm core}$. We find that across the range of neck and shell
dimensions that can be varied experimentally,[13] the electron transfer rate
constants can be tuned over a broad range of values, from $\approx 10$ ps to
$100$ fs. This suggests that faster, room temperature, electron transfer
devices require thinner shells, wider necks, and larger cores.

In conclusion, we outlined a theoretical framework to calculate charge
transfer between two colloidal quantum dots connected by a bridge. Our
approach used an atomistic model to describe the electronic structure, the
normal modes of the donor/acceptor system, the electron-phonon couplings, and
used a mixed quantum-classical mean-field Ehrenfest method to describe the
time evolution of density matrix of the system. We showed that by increasing
the neck size, increasing the core size, and/or decreasing the shell thickness
of the quantum dot building blocks, it is possible to control the
hybridization energy and tune the system from the Marcus nonadiabatic, slow
electron transfer regime to a coherent, adiabatic limit, even at elevated
temperatures, with electron transfer times that are orders of magnitude faster
and scale mildly with the hybridization energy. In all regimes, the charge
transfer dynamics are mainly governed by the coupling of the excess electron
to inter-dot and torsional acoustic modes.

###### Acknowledgements.

This work was supported by the NSF–BSF International Collaboration in the
Division of Materials Research program, NSF grant number DMR-2026741 and BSF
grant number 2020618. Methods used in this work were provided by the Center
for Computational Study of Excited State Phenomena in Energy Materials
(C2SEPEM), which is funded by the U.S. Department of Energy, Office of
Science, Basic Energy Sciences, Materials Sciences and Engineering Division,
via contract no. DE-AC02- 05CH11231, as part of the Computational Materials
Sciences Program. Computational resources were provided in part by the
National Energy Research Scientific Computing Center (NERSC), a U.S.
Department of Energy Office of Science User Facility operated under contract
no. DE-AC02- 05CH11231. D.J. acknowledges the support of the Computational
Science Graduate Fellowship from the U.S. Department of Energy under grant no.
DE-SC0019323.

## Reference

## References

  * Marcus [1956] R. A. Marcus, J. Chem. Phys. 24, 966 (1956).
  * Newton and Sutin [1984] M. D. Newton and N. Sutin, Annu. Rev. Phys. Chem. 35, 437 (1984).
  * Jortner [1976] J. Jortner, J. Chem. Phys. 64, 4860 (1976).
  * Redfield [1965] A. Redfield, in  _Advances in Magnetic and Optical Resonance_, Vol. 1 (Elsevier, 1965) pp. 1–32.
  * Mohseni  _et al._ [2014] M. Mohseni, A. Aspuru-Guzik, P. Rebentrost, A. Shabani, S. Lloyd, S. F. Huelga, and M. B. Plenio, in  _Quantum Effects in Biology_, edited by G. S. Engel, M. B. Plenio, M. Mohseni, and Y. Omar (Cambridge University Press, Cambridge, 2014) pp. 159–176.
  * Heinrich  _et al._ [2021] A. J. Heinrich, W. D. Oliver, L. M. K. Vandersypen, A. Ardavan, R. Sessoli, D. Loss, A. B. Jayich, J. Fernandez-Rossier, A. Laucht, and A. Morello, Nat. Nanotechnol. 16, 1318 (2021).
  * Zusman [1980] L. Zusman, Chem. Phys. 49, 295 (1980).
  * Zusman [1983] L. Zusman, Chem. Phys. 80, 29 (1983).
  * Demadis, Hartshorn, and Meyer [2001] K. D. Demadis, C. M. Hartshorn, and T. J. Meyer, Chem. Rev. 101, 2655 (2001).
  * Fiebig  _et al._ [2001] T. Fiebig, K. Stock, S. Lochbrunner, and E. Riedle, Chem. Phys. Lett. 345, 8 (2001).
  * Zhu  _et al._ [2021] G. Y. Zhu, Y. Qin, M. Meng, S. Mallick, H. Gao, X. Chen, T. Cheng, Y. N. Tan, X. Xiao, M. J. Han, M. F. Sun, and C. Y. Liu, Nat. Commun. 12, 456 (2021).
  * Cui  _et al._ [2019] J. Cui, Y. E. Panfil, S. Koley, D. Shamalia, N. Waiskopf, S. Remennik, I. Popov, M. Oded, and U. Banin, Nat. Commun. 10, 5401 (2019).
  * Cui  _et al._ [2021] J. Cui, S. Koley, Y. E. Panfil, A. Levi, Y. Ossia, N. Waiskopf, S. Remennik, M. Oded, and U. Banin, J. Am. Chem. Soc. 143, 19816 (2021).
  * Jasrasaria and Rabani [2021] D. Jasrasaria and E. Rabani, Nano Lett. 21, 8741 (2021).
  * Wang and Zunger [1994] L. W. Wang and A. Zunger, J. Phys. Chem. 98, 2158 (1994).
  * Rabani  _et al._ [1999] E. Rabani, B. Hetényi, B. J. Berne, and L. E. Brus, J. Chem. Phys. 110, 5355 (1999).
  * Wall and Neuhauser [1995] M. R. Wall and D. Neuhauser, J. Chem. Phys. 102, 8011 (1995).
  * Toledo and Rabani [2002] S. Toledo and E. Rabani, J. Comput. Phys. 180, 256 (2002).
  * Foster and Boys [1960] J. M. Foster and S. F. Boys, Rev. Mod. Phys. 32, 300 (1960).
  * Kleier  _et al._ [1974] D. A. Kleier, T. A. Halgren, J. H. Hall, and W. N. Lipscomb, J. Chem. Phys. 61, 3905 (1974).
  * Zhou  _et al._ [2013] X. W. Zhou, D. K. Ward, J. E. Martin, F. B. van Swol, J. L. Cruz-Campa, and D. Zubia, Phys. Rev. B 88, 085309 (2013).
  * Crespo-Otero and Barbatti [2018] R. Crespo-Otero and M. Barbatti, Chemical Reviews 118, 7026 (2018).
  * McLachlan [1964] A. McLachlan, Mol. Phys. 8, 39 (1964).
  * Nijjar, Jankowska, and Prezhdo [2019] P. Nijjar, J. Jankowska, and O. V. Prezhdo, J. Chem. Phys. 150, 204124 (2019).
  * Egorov, Rabani, and Berne [1998] S. A. Egorov, E. Rabani, and B. J. Berne, J. Chem. Phys. 108, 1407 (1998).
  * Egorov, Rabani, and Berne [1999] S. A. Egorov, E. Rabani, and B. J. Berne, J. Chem. Phys. 110, 5238 (1999).
  * Shemetulskis and Loring [1992] N. E. Shemetulskis and R. F. Loring, J. Chem. Phys. 97, 1217 (1992).
  * Note [1] Quantum mechanical test calculations using the multiconfiguration time-dependent Hartree (MCTDH) method show good agreement with the results of the Ehrenfest method.
  * Thoss, Wang, and Miller [2001] M. Thoss, H. Wang, and W. H. Miller, J. Chem. Phys. 115, 2991 (2001).
  * Breuer and Petruccione [2002] H.-P. Breuer and F. Petruccione,  _The Theory of Open Quantum Systems_ (Oxford University Press, 2002).