^{1}, Sebastian Fernandez-Alberti

^{2}, Adrian E. Roitberg

^{3}and Sergei Tretiak

^{1}

### Abstract

Within the fewest switches surface hopping (FSSH) formulation, a swarm of independent trajectories is propagated and the equations of motion for the quantum coefficients are evolved coherently along each independent nuclear trajectory. That is, the phase factors, or quantum amplitudes, are retained. At a region of strong coupling, a trajectory can branch into multiple wavepackets. Directly following a hop, the two wavepackets remain in a region of nonadiabatic coupling and continue exchanging population. After these wavepackets have sufficiently separated in phase space, they should begin to evolve independently from one another, the process known as decoherence. Decoherence is not accounted for in the standard surface hopping algorithm and leads to internal inconsistency. FSSH is designed to ensure that at any time, the fraction of classical trajectories evolving on each quantum state is equal to the average quantum probability for that state. However, in many systems this internal consistency requirement is violated. Treating decoherence is an inherent problem that can be addressed by implementing some form of decoherence correction to the standard FSSH algorithm. In this study, we have implemented two forms of the instantaneous decoherence procedure where coefficients are reinitialized following hops. We also test the energy-based decoherence correction (EDC) scheme proposed by Granucci et al. and a related version where the form of the decoherence time is taken from Truhlar's Coherent Switching with Decay of Mixing method. The sensitivity of the EDC results to changes in parameters is also evaluated. The application of these computationally inexpensive ad hoc methods is demonstrated in the simulation of nonradiative relaxation in two conjugated oligomer systems, specifically poly-phenylene vinylene and poly-phenylene ethynylene. We find that methods that have been used successfully for treating small systems do not necessarily translate to large polyatomic systems and their success depends on the particular system under study.

T.N. and S.T. acknowledge support of Directed Research and Development Fund at Los Alamos National Laboratory (LANL). A.E.R. and S.F.-A. acknowledge supported from CONICET, UNQ, ANPC_{ I }T (PICT-2010-2375), NSF Grant Nos. CHE-0239120 and CHE-0808910. Los Alamos National Laboratory is operated by Los Alamos National Security, LLC, for the National Nuclear Security Administration of the U.S. Department of Energy under Contract No. DE-AC52-06NA25396. We acknowledge support of Center for Integrated Nanotechnology (CINT) and Center for Nonlinear Studies (CNLS).

I. INTRODUCTION

II. THEORETICAL METHODOLOGY

A. Standard Tully surface hopping

B. Instantaneous decoherence

C. Energy-based decoherence correction

D. Internal consistency

III. RESULTS AND DISCUSSION

A. Molecular dynamics simulations

B. Effective energy gaps

C. Relaxation rates and internal consistency

D. Energy-based decoherence correction

E. Electronic wavepackets

IV. CONCLUSIONS

### Key Topics

- Band gap
- 29.0
- Excited states
- 28.0
- Non adiabatic couplings
- 17.0
- Absorption spectra
- 11.0
- Surface states
- 9.0

## Figures

Cartoon depicting the evolution of the quantum amplitudes c i during nonradiative relaxation. The wavefunction is initialized as a pure state. The red state symbolizes the “current state.” (a) During the standard Tully procedure, the wavepacket continues to broaden but remains at relatively high energy even as the system transitions to lower energy states. (b) Following instantaneous decoherence the wavepacket is allowed to broaden, but the coefficients are reinitialized after each hop so that the new current state has a probability of unity. In this way, the wavepacket center follows the relaxation to lower energy. (c) In energy-based decoherence, the wavepacket broadens slightly after each timestep and the coefficients are rescaled to narrow the wavepacket before the next timestep. The damping factor is related to the energy separation from the current state. Population removed from other states is deposited into the current state so that the center of the wavepacket follows the relaxation to lower energy.

Cartoon depicting the evolution of the quantum amplitudes c i during nonradiative relaxation. The wavefunction is initialized as a pure state. The red state symbolizes the “current state.” (a) During the standard Tully procedure, the wavepacket continues to broaden but remains at relatively high energy even as the system transitions to lower energy states. (b) Following instantaneous decoherence the wavepacket is allowed to broaden, but the coefficients are reinitialized after each hop so that the new current state has a probability of unity. In this way, the wavepacket center follows the relaxation to lower energy. (c) In energy-based decoherence, the wavepacket broadens slightly after each timestep and the coefficients are rescaled to narrow the wavepacket before the next timestep. The damping factor is related to the energy separation from the current state. Population removed from other states is deposited into the current state so that the center of the wavepacket follows the relaxation to lower energy.

The model molecular systems studied with the NA-ESMD approach: A 3-ring oligomer of PPV (PPV3), and a system composed of meta-substituted linear PPE segments of 2-, 3-, and 4-rings (2-3-4 PPE).

The model molecular systems studied with the NA-ESMD approach: A 3-ring oligomer of PPV (PPV3), and a system composed of meta-substituted linear PPE segments of 2-, 3-, and 4-rings (2-3-4 PPE).

(Top) Equilibrium absorption spectrum for PPV3 showing absorbance from the ground state (S0 → S n ), the lowest excited state (S1 → S n ) and equilibrium density of states. (Bottom) Density of excited states for the lowest 15 excited states at T = 300 K computed using all initial ground state configurations. S m is the state with the largest oscillator strength from S1, and the excitation is performed using a simulated laser pulse centered at λ = 245 nm with FWHM = 100 fs.

(Top) Equilibrium absorption spectrum for PPV3 showing absorbance from the ground state (S0 → S n ), the lowest excited state (S1 → S n ) and equilibrium density of states. (Bottom) Density of excited states for the lowest 15 excited states at T = 300 K computed using all initial ground state configurations. S m is the state with the largest oscillator strength from S1, and the excitation is performed using a simulated laser pulse centered at λ = 245 nm with FWHM = 100 fs.

(Top) Equilibrium absorption spectrum for 2-3-4 PPE showing absorbance from the ground state (S0 → S n ) and equilibrium density of states. (Bottom) Density of excited states for the lowest 6 excited states at T = 300 K computed using all initial ground state configurations. The excitation is performed using a simulated laser pulse centered at λ = 346 nm with FWHM = 100 fs.

(Top) Equilibrium absorption spectrum for 2-3-4 PPE showing absorbance from the ground state (S0 → S n ) and equilibrium density of states. (Bottom) Density of excited states for the lowest 6 excited states at T = 300 K computed using all initial ground state configurations. The excitation is performed using a simulated laser pulse centered at λ = 346 nm with FWHM = 100 fs.

Histograms of the energy gap between S1 and S2 are shown for PPV3 (top) and 2-3-4 PPE (bottom), while the system is evolving on S1 ( ) and S2 ( ) as well as the energy gap for the final effective hop to S1 ( ) during Standard Tully dynamics and for G-EDC. Energy gaps are smaller in 2-3-4 PPE and the final effective hop corresponds to the most probably energy gap for S2. In PPV3, energy gaps are large, and the final effective hop can only be made from lower energies, which are infrequent in the S2 spectrum.

Histograms of the energy gap between S1 and S2 are shown for PPV3 (top) and 2-3-4 PPE (bottom), while the system is evolving on S1 ( ) and S2 ( ) as well as the energy gap for the final effective hop to S1 ( ) during Standard Tully dynamics and for G-EDC. Energy gaps are smaller in 2-3-4 PPE and the final effective hop corresponds to the most probably energy gap for S2. In PPV3, energy gaps are large, and the final effective hop can only be made from lower energies, which are infrequent in the S2 spectrum.

(Top) Combined population rise of the two lowest energy states (S1 + S2) and (bottom) population rise of the lowest energy S1 state. Results are shown for ID-S, ID-A, and G-EDC using the recommended parameters (C = 1; E0 = 0.1 hartree) compared to the Standard Tully algorithm. Populations for the classical system are depicted as solid lines while the corresponding dashed lines represent the average quantum probabilities.

(Top) Combined population rise of the two lowest energy states (S1 + S2) and (bottom) population rise of the lowest energy S1 state. Results are shown for ID-S, ID-A, and G-EDC using the recommended parameters (C = 1; E0 = 0.1 hartree) compared to the Standard Tully algorithm. Populations for the classical system are depicted as solid lines while the corresponding dashed lines represent the average quantum probabilities.

Comparison of (top) G-EDC and (bottom) T-EDC populations using the recommended parameters (C = 1; E0 = 0.1 hartree). For PPV3, the populations are shown for S1, S2, and the initial S m state. Populations for S1, S2, and S3 are shown for 2-3-4 PPE. Classical populations and average quantum probabilities are depicted as solid and dashed lines, respectively.

Comparison of (top) G-EDC and (bottom) T-EDC populations using the recommended parameters (C = 1; E0 = 0.1 hartree). For PPV3, the populations are shown for S1, S2, and the initial S m state. Populations for S1, S2, and S3 are shown for 2-3-4 PPE. Classical populations and average quantum probabilities are depicted as solid and dashed lines, respectively.

Combined population rise for the lowest two excited states (S1 + S2) using G-EDC algorithm with varying parameters (C, E0) for (top) PPV3 and (bottom) 2-3-4 PPE. Classical populations are shown as solid lines and the corresponding dashed lines represent the average quantum probabilities.

Combined population rise for the lowest two excited states (S1 + S2) using G-EDC algorithm with varying parameters (C, E0) for (top) PPV3 and (bottom) 2-3-4 PPE. Classical populations are shown as solid lines and the corresponding dashed lines represent the average quantum probabilities.

Quantum wavepacket evolution for relaxation in PPV3 following excitation to S m compared for the different decoherence methods. The height of each point corresponds to the population of the state with a given energy, where the energy is plotted as the energy difference with respect to the current running state, α. Wavepackets are plotted at 2 fs intervals for the first 100 fs of the dynamics.

Quantum wavepacket evolution for relaxation in PPV3 following excitation to S m compared for the different decoherence methods. The height of each point corresponds to the population of the state with a given energy, where the energy is plotted as the energy difference with respect to the current running state, α. Wavepackets are plotted at 2 fs intervals for the first 100 fs of the dynamics.

Quantum wavepacket evolution for relaxation in 2-3-4 PPE following excitation to S n compared for the different decoherence methods. The height of each point corresponds to the population of the state with a given energy, where the energy is plotted as the energy difference with respect to the current running state, α. Wavepackets are plotted at 1 fs intervals for the first 50 fs of the dynamics.

Quantum wavepacket evolution for relaxation in 2-3-4 PPE following excitation to S n compared for the different decoherence methods. The height of each point corresponds to the population of the state with a given energy, where the energy is plotted as the energy difference with respect to the current running state, α. Wavepackets are plotted at 1 fs intervals for the first 50 fs of the dynamics.

Time evolution of the ensemble averaged ratio of coefficients for the current state, α, and the state directly below in energy, α − 1, is compared for the different decoherence methods in (top) PPV3 and (bottom) 2-3-4 PPE.

Time evolution of the ensemble averaged ratio of coefficients for the current state, α, and the state directly below in energy, α − 1, is compared for the different decoherence methods in (top) PPV3 and (bottom) 2-3-4 PPE.

## Tables

Energy gaps between states S1 and S2 for the AM1 optimized structures and averaged during dynamics while running on S1, S2, and the final effective hop to S1.

Energy gaps between states S1 and S2 for the AM1 optimized structures and averaged during dynamics while running on S1, S2, and the final effective hop to S1.

Article metrics loading...

Full text loading...

Commenting has been disabled for this content