University of Southern California Dissertations and Theses

Adiabatic and non-adiabatic molecular dynamics in nanoscale systems: theory and applications

(USC Thesis Other)

Adiabatic and non-adiabatic molecular dynamics in nanoscale systems: theory and applications

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content ADIABATIC AND NON-ADIABATIC MOLECULAR DYNAMICS IN NANOSCALE SYSTEMS: THEORY AND APPLICATIONS by Parmeet Nijjar A Dissertation Presented to the FACULTY OF THE USC GRADUATE SCHOOL UNIVERSITY OF SOUTHERN CALIFORNIA In Partial Fulfillment of the Requirements for the Degree DOCTOR OF PHILOSOPHY (CHEMICAL PHYSICS) May 2019 Copyright 2019 Parmeet Nijjar Acknowledgements It is such a happy coincidence that Professor Oleg Prezhdo moved to USC in the Fall of 2014, when I was deciding which research group to join. My Ph.D. journey was such a pleasant experience under Oleg’s supervision. Oleg is not only a great research advisor who taught me how to approach complex scientific problems in systematic and creative ways, he is also an effective mentor, who taught me the value of perseverance and challenging oneself in pursuit of excellence. Oleg, thank you for all your support. I am very grateful to you. I also want to thank all the members of Oleg’s group, Dr. Alexey Akimov, Dr. Ben Nebgen, Dr. Joanna Jankowska, Dr. Saugata Pal, Dr. Yisiang Wang, Andrew Sifain, Dhara Trivedi, and Linqiu Li for being a great support system throughout my years in graduate school. With your help, quantum chemistry doesn’t seem all that intimidating anymore. Thank you! I owe a great debt to my undergraduate Professors at IIT Kharagpur, especially Dr. Sanjoy Bandopadhyay, who introduced me to computational chemistry. That first course I took in Computational Chemistry taught by him sparked my interest in the field and led to this journey of pursuing the Ph.D. I am also indebted to Dr. Vannajan Sanghiran Lee, who gave me the invaluable experience of working in a research group in an international setting and paved the way for my career at USC. Thanks to Dr. Srabani Taraphder and ii Dr. Swagata Dasgupta for their contagious passion for chemistry that motivated and inspired me to follow the path to a Ph.D. in the subject. Thanks to Professors Jahan Dawlaty, Curt Wittig, and Anna Krylov who have taught and advised me at USC. Professor Wittig’s encouragement in the form of lemon squares, brownies and chocolate strawberries in class and his constant support and advise through- out made graduate school so much easier. Dr. Wittig, you are a great mentor and friend, and I will continue to seek your mentorship even after graduation. I have been fortunate to have such a good support system in graduate school and I recognize and appreciate it. A special thanks to the members of my qualifying exam and thesis committees: Professors Oleg Prezhdo, Curt Wittig, Surya Prakash, Jahan Dawlaty, Aichiro Nakano, and Moh El-Naggar. I appreciate you making the time given your busy schedules. Thank you! I had the opportunity to work with Dr. Sergei Tretiak and Dr. Tammie Nelson at Los Alamos National Laboratory. A significant portion of the work presented in this thesis, chapter 4 was completed under their guidance. Thank you for giving me the opportunity to experience working in a prestigious national organization and devoting so much of your time to help me. I have learned a lot from both of you. I have met some of my dearest friends in graduate school, among others, Buddhinie, Deepak, Archith, Bibek, Debanjan, Amruta, Natalie, Arman, Bailey, Andrew, and Bren- dan Gifford. I cannot thank you enough for all the help, guidance, and fun we had to- gether that helped me tackle everything these last four years. Thank you for being my family at USC and filling my years here with laughter, dinners, hiking and many more joyful memories. A big thanks goes to my family. My mother, Kulwinder Kaur, taught me to dream big but more importantly, to go get those dreams. Thanks to my father, Jaswant Singh, iii for providing us with everything from food and shelter and always encouraging me to achieve my goals. I could not imagine being who I am and where I am today if it was not for my parents’ encouragement and faith in me. My sister, Jaspreet Jagpal, never failed to give me the good advice and words of encouragement whenever I needed a push to surmount the obstacles along the way. My brother, Dilraj always inspires me by his achievements and with his encouraging words. Thank you guys for always being there. Last but not the least, thank you to Sidharth Madan, for all his love and support. iv Table of Contents Acknowledgements ii List of Tables vii List of Figures viii Abstract xv Chapter 1: Introduction 1 1.1 Adiabatic Molecular Dynamics . . . . . . . . . . . . . . . . . . . . . . 2 1.1.1 Born Oppenheimer Approximation . . . . . . . . . . . . . . . . 3 1.2 Non-adiabatic Molecular Dynamics . . . . . . . . . . . . . . . . . . . 5 1.2.1 Ehrenfest Dynamics . . . . . . . . . . . . . . . . . . . . . . . . 6 1.2.2 Fewest Switches Surface Hopping . . . . . . . . . . . . . . . . 10 1.2.3 Decoherence . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 Chapter 1 References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 Chapter 2: Atomistic Analysis of Room Temperature Quantum Coherence in Two-Dimensional CdSe Nanostructures 21 2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 2.2 Theoretical Methods and Simulation details . . . . . . . . . . . . . . . 24 2.3 Results and discussion . . . . . . . . . . . . . . . . . . . . . . . . . . 27 2.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40 Chapter 2 References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42 Chapter 3: Conversion of He(2 3 S) to He 2 (a 3 + u ) in Helium droplets 47 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47 3.2 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . 49 3.2.1 3-atom linear configurations . . . . . . . . . . . . . . . . . . . 54 3.2.2 4-atom linear configurations . . . . . . . . . . . . . . . . . . . 57 3.2.3 3-atom nonlinear configurations . . . . . . . . . . . . . . . . . 58 3.2.4 10-atom 2D configurations . . . . . . . . . . . . . . . . . . . . 59 3.3 Computational Details . . . . . . . . . . . . . . . . . . . . . . . . . . . 60 v 3.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 3.5 Appendix A: Trajectories for 3-atom linear geometries, 3-atom non- linear geometries and 4 atoms constrained to a straight line. . . . . . . . 63 Chapter 3 References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73 Chapter 4: The Trivial Crossing Problem of FSSH 76 4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76 4.2 Theoretical Methodology . . . . . . . . . . . . . . . . . . . . . . . . . 79 4.2.1 Non-adiabatic Excited State Molecular Dynamics . . . . . . . . 79 4.2.2 Transition Density Analysis . . . . . . . . . . . . . . . . . . . 80 4.2.3 Simulation Details . . . . . . . . . . . . . . . . . . . . . . . . 81 4.3 Results and discussion . . . . . . . . . . . . . . . . . . . . . . . . . . 82 4.3.1 Absorption Spectra . . . . . . . . . . . . . . . . . . . . . . . . 82 4.3.2 Photoinduced dynamics of non-interacting polyfluorene segments 84 4.3.3 Photoinduced dynamics of interacting spirobifluorene . . . . . . 88 4.3.4 Relaxation time-scales . . . . . . . . . . . . . . . . . . . . . . 92 4.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93 4.5 Appendix A: Natural Transition Orbitals for S 1 - S 10 excited states of spirobifluorene at 300 K . . . . . . . . . . . . . . . . . . . . . . . . . . 94 Chapter 4 References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97 Chapter 5: Ehrenfest Dynamics with Decoherence and Detailed Balance 100 5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100 5.2 Theory: Formulation of the Ehrenfest method with decoherence and de- tailed balance corrections . . . . . . . . . . . . . . . . . . . . . . . . . 104 5.2.1 The Ehrenfest method . . . . . . . . . . . . . . . . . . . . . . . 104 5.2.2 Coherence penalty functional . . . . . . . . . . . . . . . . . . . 108 5.2.3 Detailed balance and thermodynamic equilibrium . . . . . . . . 110 5.2.4 Pure-dephasing/decoherence time . . . . . . . . . . . . . . . . 112 5.2.5 Simulation details . . . . . . . . . . . . . . . . . . . . . . . . . 113 5.3 Results: Charge injection, cooling and recombination in a metallic-tip/porphyrin system . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114 5.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128 Chapter 5 References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 130 Chapter 6: Future Work 134 Chapter 6 References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 138 vi List of Tables 2.1 Canonically Averaged Energy, E of Heavy (HX) and Light (LX) Hole Excitations in the CdSe Nanoplatelet a . . . . . . . . . . . . . . . . . . 28 2.2 Canonically Averaged Heavy (HX) and Light (LX) Hole Excitations in the CdSe/CdZnS Core/Shell Nanoplatelet . . . . . . . . . . . . . . . . 30 2.3 Canonically Averaged Heavy (HX) and Light (LX) Hole Excitations in the CdSe/CdZnS Core/Shell Nanoplatelets with the S/Se Exchange at the Interface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 4.1 Rise-time of the average population of state S 1 (S a 1 for FSSH, S d 1 for the Min-Cost method) for separated system and spirobifluorene. . . . . . . 92 5.1 Pure-dephasing time (fs) (lower triangle) and time-averaged absolute value of nonadiabatic coupling (meV) (upper triangle) . . . . . . . . . 116 vii List of Figures 1.1 The left panel shows stochastic hops from one surface to the other when passing through a region of strong coupling in multiple trajectories in surface hopping method. The right panel shows that the nuclear tra- jectory starts evolving on a mean surface in Ehrenfest dynamics after passing through the region of strong coupling. . . . . . . . . . . . . . . 10 1.2 Schrodinger’s cat. Decoherence are all those phenomena that result in the cat being either dead or alive. The cat is not in a superposition of the two states anymore. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 2.1 Optimized geometries of (a) CdSe nanoplatelet, (b) CdSe/ CdZnS core/shell nanoplatelet, (c) CdSe/CdZnS core/shell nano-platelet with S/Se exchange, and (d) CdSe/CdZnS core/shell nanoplatelet with Zn/Cd exchange. The bracket in part a identifies a unit layer. The blue, magenta, yellow, and red balls represent Cd, Se, Zn, and S atoms, respectively. The S/Se and Zn/Cd exchanges are shown in parts c and d by arrows. . . . . . . . . . 26 2.2 Projected density of states (DOS) of (a) CdSe and (b) CdSe/CdZnS core/shell nanoplatelets. For CdSe, the blue line indicates the total DOS, whereas the orange and magenta lines represent the projected DOS of Cd and Se atoms, respectively. For the core/shell platelet, the black line rep- resents the total DOS, whereas the blue and red lines indicate the core and shell DOS, respectively. . . . . . . . . . . . . . . . . . . . . . . . 28 2.3 (a) Room temperature pure-dephasing functions for pairs of heavy and light hole excitations in the CdSe nanoplatelet. The red and blue lines represent the HX1-LX1 and HX2-LX3 excitation pairs, whose energy differences most closely match the experiment, Table 5.1. The inset in panel a displays the optical spectrum. The individual transitions shown by the red line are Gaussian-broadened to obtain the black line. (b) Autocorrelation functions of the energy gaps between the same pairs of excitations. (c) Fourier transforms of (i) HX1, (ii) HX2, (iii) LX1, and (iv) LX3 excitation energies. . . . . . . . . . . . . . . . . . . . . . . . 29 viii 2.4 (a) Room temperature pure-dephasing functions for pairs of heavy and light hole excitations in the CdSe/CdZnS core/shell nanoplatelet. The red and blue lines represent the HX1-LX1 and HX1-LX2 excitation pairs, whose energy differences most closely match the experiment, Ta- ble 2.2. The inset in panel a displays the optical spectrum. The individ- ual transitions shown by the red line are Gaussian-broadened to obtain the black line. (b) Autocorrelation functions of the energy gaps between the same pairs of excitations. (c) Fourier transforms of (i) HX1, (ii) LX1, and (iv) LX2 excitation energies. . . . . . . . . . . . . . . . . . 31 2.5 (a) Room temperature pure-dephasing functions obtained by averaging over all pairs of heavy and light hole excitations in CdSe (upper panel), CdSe/CdZnS (middle panel), and CdSe/CdZnS core/shell nanoplatelets with S/Se exchange at the interface (lower panel). (b) Room temperature optical spectrum of the CdSe/CdZnS core/shell nanoplatelet with Zn/Cd exchange at the interface. The individual transitions shown by the red line are Gaussian- broadened to obtain the black line. . . . . . . . . . . 33 2.6 Plots of charge densities of different electronic states that are involved in the transitions defined in Tables 5.1, 2.2 and 2.3 for (a) CdSe, (b) CdSe/CdZnS and (c) CdSe/CdZnS with S/Se exchange, respectively, at 300 K. VBM and CBM refer to valence band maximum and conduc- tion band minimum. For CdSe, the energy of the VB 0 , VB 00 , VB 000 , CB 0 and CB 00 states correspond to VBM-0.08 eV , VBM-0.15 eV , VBM-0.19 eV , CBM+0.12 eV and CBM+0.26 eV , respectively. The HX1, HX2, LX1 and LX3 excitations arise from the following electronic transitions, VB 00 ! CBM, VB 000 ! CBM, VB 00 ! CBM 0 and VB 0 ! CBM 00 . For CdSe/CdZnS, VB 0 , VB 00 and CB 0 represent states with the ener- gies VBM-0.04 eV , VBM-0.19 eV and CBM+0.19 eV , respectively. The HX1, LX1 and LX2 excitations correspond to the VB 0 ! CBM, VB 00 ! CBM, VB 0 ! CBM 0 transitions, respectively. For CdSe/CdZnS with S/Se exchange, VB 0 , VB 00 , VB 000 and CB 0 represent states with the ener- gies VBM-0.05 eV , VBM-0.2 eV and CBM+0.10 eV , respectively. The HX1, HX2, HX3, LX1 and LX2 excitations correspond to the following transitions VB 0 ! CBM, VB 00 ! CBM, VB 0 ! CBM 0 , VB 000 ! CBM and VB 000 ! CBM 0 , respectively. . . . . . . . . . . . . . . . . . . . . . 37 ix 2.7 (a) Room temperature pure-dephasing functions for pairs of heavy and light hole excitations in the CdSe/CdZnS core/shell nanoplatelet with S/Se exchange at the interface. Red, blue, magenta, orange, light green, and dashed dark green lines represent HX1-LX1, HX1-LX2, HX2-LX1, HX2-LX2, HX3-LX1, and HX3-LX2 excitation pairs, respectively, see Table 2.3. The inset in panel a displays the optical spectrum. The individual transitions shown by the red line are Gaussian-broadened to obtain the black line. (b) Autocorrelation functions of the energy gaps between the same pairs of excitations, maintaining the color code used for the dephasing functions. (c) Fourier transforms of (i) HX1, (ii) HX2, (iii) HX3, (iv) LX1, and (v) LX2 excitation energies. . . . . . . . . . . 38 2.8 Schematic showing the coupling of the HX and LX excitons to the same phonon modes leading to long quantum coherence timescales. . . . . . 40 3.1 The a 3 + u and c 3 + g potential energy curves correlate to He* + He. The vertical scale for the a 3 + u barrier region is expanded in the insert. . . . 50 3.2 (a) 3 atoms in a linear configuration. (b) 4 atoms in a linear configu- ration. (c) 3-atom bent configurations included 0 = 30 , 60 , and 90 , and a T-shaped starting geometry. Initial interatomic distances satisfy R 0 2 R 0 1 and R 0 2 R 0 3 . . . . . . . . . . . . . . . . . . . . . . . . . . 51 3.3 1D, 3-atom lowest adiabat slices: NTOs are shown for 2.5 and 5.0 ˚ A (note dashed lines in top panel). The NTOs show holes and electrons on the left and right as indicated in the top panel; atom labels (1,2,3) are as per Figure 3.2a. As with MOs, red and blue orbital colors only have meaning within an NTO (phases of1). Orbital cutoff is chosen for viewing convenience. . . . . . . . . . . . . . . . . . . . . 53 3.4 (a) is a view of energy versus R 1 and R 2 for 3 atoms. (b) shows a top view. The rectangle contains the (R 0 1 ; R 0 2 ) combinations used in the trajectories. Gradients enable one to see how trajectories are launched from different (R 0 1 ; R 0 2 ). The arrows originate at the points where the gradients are evaluated. . . . . . . . . . . . . . . . . . . . . . . . . . . 54 3.5 The smooth curve is the calculated radial distribution function g(r) for helium at 1.21 K. Points show experimental data. Retrieved from refer- ence 35. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55 3.6 Trajectories for (a) linear 3-atom cases (Figure 3.2a), and (b) linear 4- atom cases (Figure 3.2b). Note that R 2 (t) oscillation transfers to R 3 (t) in the first two panels. (c) Snapshots of a trajectory for 3-atoms with R 0 1 = 3.1 ˚ A, R 0 2 = 2.5 ˚ A, and 0 = 90 (Figure 3.2c). Smaller initial distances (R 0 2 ) are characteristic of the rising part of g(r). Larger initial distances (R 0 1 ) are characteristic of the peak region of g(r). . 56 3.7 Snapshots of a trajectory for 10-atoms showing exciton hopping. . . . . 59 x 3.8 Trajectories for linear 3-atom cases for R 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.6 - 2.8 ˚ A. The R 2 pair does not always remain intact, but it always passes through the region of the potential minimum as noted by the downward thrusts of the red traces. . . . . . . . . . . . . . . . . . . . . . . . . . 63 3.9 Trajectories for 2D 3-atom geometries with = 30 . The initial distances lie in the range R 0 1 = 3.0 - 3.4 ˚ A and R 0 2 = 2.6 - 2.8 ˚ A. Results are similar to those for the corresponding 1D 3-atom cases. . . . . . . . . . . . . . 64 3.10 Trajectories for 2D 3-atom geometries with = 60 . The initial dis- tances lie in the range R 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.6 - 2.8 ˚ A. Results show He 2 formation takes place in more initial configurations than the corresponding 1D 3-atom cases (Figure 3.8) as well as the 2D 30 ge- ometries (Figure 3.9). . . . . . . . . . . . . . . . . . . . . . . . . . . . 65 3.11 Trajectories for 2D 3-atom geometries with = 90 . The initial distances lie in the rangeR 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.4 - 2.6 ˚ A. The diatom pair, even when it initially localizes the exciton, is discouraged from bonding because of repulsion and the angular momentum it receives by being pushed from outside its center-of-mass. . . . . . . . . . . . . . . . . . 66 3.11 Trajectories for 2D 3-atom geometries with = 90 . The initial distances lie in the rangeR 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.7 - 2.9 ˚ A. The diatom pair, even when it initially localizes the exciton, is discouraged from bonding because of repulsion and the angular momentum it receives by being pushed from outside its center-of-mass. . . . . . . . . . . . . . . . . . 67 3.12 Trajectories for 2D 3-atom geometries with T-shaped initial configura- tions. The initial distances lie in the range R 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.4 - 2.6 ˚ A. For R 0 2 = 2.4 - 2.6 ˚ A, initial localization of the exciton at the diatom is retained and the third atom separates from He 2 . . . . . . . . . 68 3.12 Trajectories for 2D 3-atom geometries with T-shaped initial configura- tions. The initial distances lie in the range R 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.7 - 2.9 ˚ A. For R 0 2 = 2.7 - 2.9 ˚ A, He 2 did not survive in any of the trajectories. 69 3.13 Trajectories for symmetric (R 0 1 = R 0 3 ) arrangements of 4-atoms arranged in a linear fashion. The initial distances lie in the range R 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.6 - 3.0 ˚ A. The middle two atoms pass through the region of the potential minimum. For R 0 2 = 2.6 ˚ A, 2.7 ˚ A, all R 0 1 values yielded He 2 ; for R 0 2 = 2.8 ˚ A, R 0 1 = 3.0 - 3.3 ˚ A yielded He 2 ; for R 0 2 = 2.9 ˚ A, R 0 1 = 3.0 - 3.2 ˚ A yielded He 2 ; and for R 0 2 = 3.0 ˚ A, R 0 1 = 3.0 ˚ A and 3.1 ˚ A yielded He 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70 xi 3.13 Trajectories for symmetric (R 0 1 = R 0 3 ) arrangements of 4-atoms arranged in a linear fashion. The initial distances lie in the range R 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.6 - 3.0 ˚ A. The middle two atoms pass through the region of the potential minimum. For R 0 2 = 2.6 ˚ A, 2.7 ˚ A, all R 0 1 values yielded He 2 ; for R 0 2 = 2.8 ˚ A, R 0 1 = 3.0 - 3.3 ˚ A yielded He 2 ; for R 0 2 = 2.9 ˚ A, R 0 1 = 3.0 - 3.2 ˚ A yielded He 2 ; and for R 0 2 = 3.0 ˚ A, R 0 1 = 3.0 ˚ A and 3.1 ˚ A yielded He 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71 3.14 Trajectories for non-symmetric arrangements of 4-atoms arranged in a linear fashion. R 0 1 is fixed at 3.4 ˚ A, the other initial distances lie in the range R 0 2 = 2.6 - 2.9 ˚ A and R 0 3 = 2.8 - 3.3 ˚ A. Exciton transfer took place in 10 of the 32 initial configurations. . . . . . . . . . . . . . . . . . . . 72 4.1 (a) Chemical structure of spirobifluorene (top) and the corresponding structure of polyfluorene segments separated by 100 nm (bottom). (b) Natural transition orbitals for the S 0 ! S 1 and S 0 ! S 2 transitions at the optimized ground-state geometry for spirobifluorene. (c) Simulated ab- sorption spectra with contributions of different excited states for spirob- ifluorene (solid lines) and separated polyfluorenes (dashed lines) at 300 K. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83 4.2 Simulation results for polyfluorene segments separated by 100 nm. Columns 1 and 2 correspond to (i) standard FSSH, and (ii) Min-Cost methods, re- spectively. (a) Average population on each electronic surface as a func- tion of time obtained from the fraction of trajectories in each state. (b) Time-dependence of the average fraction of TD localized on segment 1 [( gf ) 2 A > 0:5 at t = 1000 fs] and segment 2. (c) Time dependence of the average participation number. (d) Histogram of the number of TD lo- calization changes between the polyfluorene segments during dynamics. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85 4.3 (a) Absolute value of the non-adiabatic coupling terms (NACT) between states (S 1 , S 2 ), (S 2 , S 3 ) and (S 1 , S 3 ) averaged over all trajectories, (b) Relative probability of changing TD localization between segments as a function of time, for polyfluorene segments separated by 100 nm. Columns 1 and 2 correspond to (i) standard FSSH, and (ii) Min-Cost methods, respectively. . . . . . . . . . . . . . . . . . . . . . . . . . . 87 4.4 Simulation results for spirobifluorene. Columns 1 and 2 correspond to (i) standard FSSH, and (ii) Min-Cost method, respectively. (a) Aver- age population on each electronic surface as a function of time obtained from the fraction of trajectories in each state. (b) Time-dependence of the average fraction of TD localized on segment 1 [( gf ) 2 A > 0:5 at t = 1000 fs] and segment 2. (c) Time dependence of the average participa- tion number. (d) Histogram of the number of TD localization changes between the polyfluorene segments during dynamics. . . . . . . . . . . 89 xii 4.5 (a) Absolute value of the non-adiabatic coupling terms (NACT) between states (S 1 , S 2 ), (S 2 , S 3 ) and (S 1 , S 3 ) averaged over all trajectories, (b) Relative probability of changing TD localization between segments as a function of time, for spirobifluorene. Columns 1 and 2 correspond to (i) standard FSSH, and (ii) Min-Cost methods, respectively. . . . . . . . . 91 5.1 (a) Top view on the tip-porphyrin system. (b) The studied nonadiabatic processes and charge densities of the corresponding states. . . . . . . . 115 5.2 (a) Phonon-driven evolution of the electronic state energies, S0-S4. (b) Energy gaps between excited states, S1-S4. (c) Absolute value of nona- diabatic coupling. (d) Un-normalized autocorrelation functions. (e) In- fluence spectra. (f) Pure-dephasing functions for transitions between ad- jacent electronic states. The energy gap between states S1 and S2 drops below k B T=0.025eV at room temperature (part b). Transition between this pair of states exhibits the largest nonadiabatic coupling (part c) and couples to high frequency vibrations (part e), arising from the porphyrin molecule that supports these states. The energy gap fluctuation is largest for the S0-S1 state pair (part d) which exhibits fastest dephasing (part f). 117 5.3 Dynamics of electronic state populations displayed up to 1 ns (left pan- els) and focused on the first 100 ps (right panels) calculated within the Ehrenfest-DDB method (upper panels) and the DISH method (lower panels) for the partial electron injection process (S4!S3). Both meth- ods produce similar timescales, with Ehrenfest-DBB marginally slower than DISH. Whereas the DISH method shows almost complete relax- ation, the Ehrenfest-DDB scheme stops at about 6% short of complete relaxation. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118 5.4 Dynamics of electronic state populations displayed up to 1 ns (left pan- els) and focused on the first 100 ps (right panels) calculated within the Ehrenfest-DDB method (upper panels) and the DISH scheme (lower panels) for the complete electron injection process (S4!S3!S2). Al- though the initial relaxation rate is marginally slower with Ehrenfest- DDB as compared to DISH, electronic populations equilibrate at a sim- ilar state with both methods. . . . . . . . . . . . . . . . . . . . . . . . 121 xiii 5.5 Dynamics of electronic state populations displayed up to 1 ns (left pan- els) and focused on the first 100 ps (right panels) calculated within the Ehrenfest-DDB scheme (upper panels) and the DISH scheme (lower panels) for the combined injection and cooling processes (S4!S1). The initial dynamics are similar and slightly slower with Ehrenfest-DDB compared to DISH. At the long-time limit the populations of states S1 and S2 reach different values. The different behavior of the two methods can be attributed to the facts that the energy gap between these states drops below kBT (Figure 5.2b) and that the nonadiabatic coupling is large (Figure 5.2c). . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122 5.6 Dynamics of electronic state populations displayed up to 1 ns (left pan- els) and focused on the first 100 ps (right panels) calculated within the Ehrenfest-DDB scheme (upper panels) and the DISH scheme (lower panels) just for the pair of states S1 and S2 that exhibit different long- time limits during the S4!S1 dynamics shown in Figure 5.5. The dif- ferences arise because the S1-S2 energy gap drops below k B T (Figure 5.2b) and the nonadiabatic coupling is large (Figure 5.2c). . . . . . . . 124 5.7 Dynamics of electronic state populations displayed up to 1 ns (left pan- els) and focused on the first 100 ps (right panels) calculated within the Ehrenfest-DDB scheme (upper panels) and the DISH scheme (lower panels) for the complete set of processes, including charge injection, cooling and recombination (S4!S0). The initial dynamics are similar and slightly slower in Ehrenfest-DDB compared to DISH. At long times, the Ehrenfest-DDB dynamics slows down considerably, with the larges differences observed for the populations of states S0, S1 and S2. See Figures 5.5 and 5.6 for additional insights. . . . . . . . . . . . . . . . 125 5.8 Dynamics of electronic state populations displayed up to 1 ns (left pan- els) and focused on the first 100 ps (right panels) calculated within the Ehrenfest-DDB scheme (upper panels) and the DISH scheme (lower panels) for the recombination process (S1!S0). Both DISH and Ehrenfest- DDB require more than 1 ns to reach equilibrium, owing to the large en- ergy gap (Figure 5.2a) and small nonadiabatic coupling (Figure 5.2c and Table 5.1) between the states. Ehrenfest-DDB approaches equilibrium faster with two states (S1!S0) than with five states (S4!S0), Figure 5.7. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126 xiv Abstract Adiabatic as well as non-adiabatic molecular dynamics methods are essential tools to simulate many important physical processes. Molecular dynamics simulations pro- vide mechanistic information that often cannot be detected by experiment alone and also serve as a predictive tool to simulate new research. However, sophisticated mod- eling of ultrafast processes, especially those involving excited states, is intellectually challenging and computationally demanding. The work presented here provides several diverse examples of challenging gas phase and condensed phase ultrafast photoinduced processes. For each system, we offer a detailed account of the computational modeling scheme employed for an accurate description of the chemical and physical properties of interest. At the same time, we highlight the limitations of the particular dynamics methods used. In Chapter 1, we start with an overview of the adiabatic and the non-adiabatic dy- namics. This is followed by a description of the primary non-adiabatic dynamics meth- ods, the Ehrenfest dynamics and the Fewest Switches Surface Hopping method. We also describe the major limitations of the two approaches, especially the neglect of decoher- ence effects among others. xv This is followed by an example of the application of adiabatic molecular dynam- ics in Chapter 2 to study the atomistic origin of room temperature coherence between the light-hole and heavy-hole excitations in two dimensional CdSe and CdSe/CdZnS core/shell nanoplatelets, synthesized and studied recently by optical spectroscopy. The molecular dynamics calculations allowed us to mimic closely the experimental systems and measurements, and to elucidate the role of defects that are inevitably present in ex- periment. Chapter 3 provides another application of quantum-classical molecular dy- namics to study the conversion of He(2 3 S) to He 2 (a 3 + u ) in large helium nanodroplets. He 2 is a likely participant in the production of He + 4 on surfaces of large helium droplets which shows an anomalously intense peak in mass spectra of large helium droplets. Ab initio molecular dynamics with highly accurate and robust equation-of-motion coupled- cluster (EOM-CC) theory, showed rapid and efficient conversion of He to He 2 in sys- tems composed of small numbers of Helium atoms. The Born-Oppenheimer dynamics is followed by the non-adiabatic dynamics with electronic transitions in chapters 4 and 5. Chapter 4 analyzes the fewest switches surface hopping (FSSH) algorithm and brings attention to one of it’s limitations, the trivial crossing problem. It provides a comparison of the standard fewest switches sur- face hopping (FSSH) algorithm with the Min-Cost adiabatic state identity tracking al- gorithm, which treats the trivial crossing problem of FSSH. The performance of the two methods is investigated for the internal conversion and energy transfer dynamics of a weakly coupled spiro-linked polyfluorene (spirobifluorene) system. The results re- veal that FSSH fails to identify trivial crossings, therefore, the system remains on the upper adiabatic surface at such crossings and reaches the lowest excited state only via non-adiabatic transitions resulting in a longer relaxation time-scale as compared to the Min-Cost method. xvi Chapter 5 addresses the other primary mixed quantum classical dynamics method, the Ehrenfest method. Ehrenfest method does not produce the Boltzmann equilibrium populations of the quantum states at a finite temperature. Additionally, both the Ehren- fest MD and the FSSH techniques are fully quantum coherent, they do not account for decoherence. In this chapter, we present a semiclassical approach for non-adiabatic molecular dynamics based on the Ehrenfest method with corrections for decoherence effects and detailed balance. The decoherence effects are described via a coherence penalty functional that drives dynamics away from regions in Hilbert space character- ized by large values of coherences. Detailed balance is described by modification of the off-diagonal matrix elements with a quantum correction factor used in semiclassical approximations to quantum time-correlation functions. Both decoherence and detailed balance corrections introduce non-linear terms to the Schr odinger equation. At the same time, the simplicity of the Ehrenfest method, which involves fully deterministic dynam- ics and requires a single trajectory for each initial condition, is preserved. In contrast, surface hopping techniques contain stochastic elements and require averaging over mul- tiple realization of stochastic process for each initial condition. The method is applied to a hybrid nanoscale system consisting of a fluorophore molecule and an STM tip studied experimentally. By investigating current-induced charge injection, cooling and recom- bination, we demonstrate that the Ehrenfest-decoherence-detailed-balance (Ehrenfest- DDB) method produces timescales that are similar to those obtained with the deco- herence induced surface hopping (DISH), which is a popular nonadiabatic molecular dynamics technique used with many condensed matter and nanoscale systems. The Ehrenfest-DDB technique provides efficient means to study quantum dynamics in large, realistic systems. xvii Finally, in chapter 6, we outline future directions and discuss further investigations of the developed computational scheme to establish the range of problems it can be applied to with useful accuracy. xviii Chapter 1: Introduction Molecular Dynamics is an indispensable tool for studying many important phenomena in all fields of science including physics, mathematics, chemistry, biology, and materials sciences. The increasing availability of faster computers (high performance computing) and simultaneous development of more efficient algorithms has made numerical simula- tions prevalent in virtually every field of research. The molecular dynamics simulations are used to study processes ranging from protein folding 1, 2 with explicit solvent inter- actions where the quantum effects are neglected to phenomena like electron transfer, 3, 4 charge transport 5, 6 and exciton dynamical processes 7–16 in relatively smaller systems where quantum mechanical behavior is central to the description of the process. With the advent of ab initio molecular dynamics, very accurate potential energy surfaces can be computed on the fly to study dynamical processes in systems of the order of 10 2 atoms. The development of techniques like Surface Hopping 17 and Ehrenfest dynam- ics 18 has enabled ab initio simulations of complex ultrafast processes with electronic transitions in condensed phase systems. Classical Path approximation allows simula- tion of non-adiabatic processes in realistic condensed phase systems of interest albeit the potentials are computed first and the nuclear sub-system does not receive feedback from the electronic configuration at every time step. In the following chapters, we will explore the application of molecular dynamics simulations to study the adiabatic dynamics in 2D CdSe nanoplatelets, the adiabatic 1 dynamics on a triplet potential energy surface of small Helium systems, as well as pho- toinduced dynamics, in particular, nonradiative relaxation in a spirobifluorene molecule. We will begin with a derivation of the theoretical framework for both adiabatic and non- adiabatic dynamics as well as outline the major approximations used. 1.1 Adiabatic Molecular Dynamics The typical starting point in the study of quantum chemical systems is the following non-relativistic Schr odinger equation H i (r;R) =E i i (r;R) (1.1) wherer = r 1 ;r 2 ;:::r n andR = R 1 ;R 2 ;:::R N represent the electron and nuclear coor- dinates respectively. The Hamiltonian of a molecular system consisting of n electrons and N nuclei is a sum of the kinetic energy of the electrons,T e , and the nuclei,T n , and the coulomb terms including the electron-nuclear attraction, V en , the electron-electron repulsion,V ee , and the nuclear-nuclear repulsion,V nn . H = n X i=1 1 2 r 2 i N X a=1 1 2M a r 2 a X a;i Z a ~ R a ~ r i + X i<j 1 j~ r i ~ r j j + X a<b Z a Z b ~ R a ~ R b (1.2) whereM a andZ a are the mass and the charge of thea-th nucleus respectively. The solution of equation 1.2 yields wavefunctions that depend on the coordinates of all electrons and nuclei. However, this equation is analytically solvable only for the H atom and can be solved very accurately for the lightest molecule, H 2 using the variational approach. For other chemical systems, the Born Oppenheimer approximation must be used. 2 1.1.1 Born Oppenheimer Approximation The Born Oppenheimer approximation 19, 20 separates the electronic and nuclear motions allowing to compute the total wavefunction of the system in two less complicated con- secutive steps. The proton mass is 3 orders of magnitude larger than the electron mass. So, it is reasonable to assume that the electronic motion will be much faster than the nu- clear motion and the electronic subsystem adjusts almost instantaneously to the changes in the nuclear subsystem. Following this assumption, we fix the nuclear geometry and ignore the nuclear kinetic energy term and solve the remaining electronic SE: H el i (r;R) =U i (R) i (r;R) (1.3) H el =T e +V en +V ee +V nn (1.4) The resulting solution (wavefunctions and energies) dependsparametrically on nuclear coordinates. At each nuclear geometry, the lowest energy solution of the electronic SE, equation 1.3, corresponds to the ground state energy and the higher energy solutions correspond to the electronically excited states. In the next step, the nuclear kinetic energy,T n , is reintroduced and the full problem, equation 1.1, is solved in terms of the solution of the electronic SE, equation 1.3. To do that, the total wavefunction, (r;R), is recast as a product of the electronic and nuclear wavefunctions. i(total) = X (r;R) i (R) (1.5) 3 where , the electronic wavefunctions with parametric dependence on R are used as a basis and the nuclear wavefunctions, i (R) are expansion coefficients to be determined. Substituting ansatz 1.5 into equation 1.1: (H el +T n ) X (r;R) i (R) =E i X (r;R) i (R) (1.6) Multiplying the above equation by on the left and integrating over the electronic coordinates, we obtain the following equation, (T n +U (R)E i ) i = X a 1 2M a X h jr 2 a i r i + 2 X h jr a i r r a i ! (1.7) The orthonormality of the electronic wavefunctions is also used in the derivation of the equation 1.7. The diagonal parts of the terms on the right, called the adiabatic correction is a correction to the electronic energy due to nuclear motion. The non-diagonal parts of the terms are called the non-adiabatic coupling terms which couple the nuclear dynamics on different electronic states. In the adiabatic approximation, we assume that the electronic states are well separated and the terms on the right of equation 1.7 are negligible. This results in a simpler eigenproblem for the nuclei moving on a single Born Oppenheimer electronic potential surface surface (PES). T n i + (U (R)E i ) i = 0 (1.8) This way, we can find the nuclear eigenstates for each electronic PES and the total wavefunction is given by i = i i (R) (1.9) 4 Adiabatic molecular dynamics with the PES calculated “on the fly” has been used to study ultrafast dynamics in organic molecules. 21–24 The trajectories can be propagated along the ground state or the excited state potential energy surfaces. The nuclei are propagated classically using the Newton’s equation of motion with the forces from the electronic potential at each time-step. The situation becomes more complex when the electron and the nuclear dynamics become non-adiabatic. When the coupling terms on the right of equation 1.7 are large, the adiabatic approximation is no longer valid. These terms can be large when the electronic wavefunction depends strongly on the nuclear geometry. For example, when the nuclear velocity is high, the non-adiabatic coupling terms can become very large. In this scenario, the electrons cannot instantaneously adjust to the nuclear motion and non-adiabatic hops between potential energy surfaces become probable. 1.2 Non-adiabatic Molecular Dynamics Non-adiabatic molecular dynamics has become increasingly invaluable for simulating and understanding complex ultrafast processes such as electron transfer, charge trans- port, exciton relaxation and recombination, etc. Computational modeling complements experiment by providing mechanistic information which cannot be determined via mea- surements alone and also serves as a predictive tool to simulate new research. For the simulation of non-adiabatic dynamics, the two most widely used approaches are the Ehrenfest method 18 and the surface hopping method. 17 In both the Ehrenfest and sur- face hopping approaches, the nuclei are treated classically while the electrons are treated quantum mechanically. The transitions (hops) among coupled excited states incorporate feedback between electronic and nuclear subsystems, i.e., the so-called back reaction of 5 the quantum onto the classical system is included. A brief discussion of the Ehrenfest and surface hopping schemes is presented below. 1.2.1 Ehrenfest Dynamics In the Ehrenfest dynamics, the BO approximation is considered; however, the electronic subsystem does not evolve adiabatically on one of the potential energy surfaces. The electrons can move between potential energy surfaces and evolve on linear combination of adiabatic states. The nuclei evolve subject to the mean-field force from the average potential energy surface. Starting with the time-dependent Schr odinger equation i~ @ @t (r;R;t) =H (r;R;t) (1.10) and using the Born-Oppenheimer separation ansatz for the wavefunction, (r;R;t) = X i i (r;R(t)) i (t;R(t)) (1.11) and projecting it onto the corresponding electronic states, one obtains the dynamics of nuclear wavepackets correlated with different electronic states: i~ @ i (t;R(t)) @t = X j ( T nucl + ~ E j (t;R) ij ~ 2 ~ d (1) ij M ~ r~ 2 ~ d (2) ij 2M ) i (t;R(t)) (1.12) where ~ d (n) ij are the nth order non-adiabatic couplings between electronic statesi andj, defined as ~ d (1) ij = h i (r;R(t))j ~ r j (r;R(t))i (1.13a) ~ d (2) ij = h i (r;R(t))jr 2 j (r;R(t))i (1.13b) 6 ~ E i (t;R)) are the eigenvalues of the stationary electronic SE and correspond to the adia- batic potential energy surfaces as explained in the section 1.1 above. The diagonal parts of the term~ 2 ~ d (2) ij 2M in Eq. 1.12 are the adiabatic correction to the electronic energy and can be included in the adiabatic energies as follows: E i (R) = ~ E i (R)~ 2 ~ d (2) ij 2M (1.14) The off-diagonal terms representing the second order NA coupling between electronic states are often neglected under the assumptions that adiabatic wavefunctions are slowly varying with respect to nuclear coordinates. 25–27 However, these terms may become im- portant in some cases, and one must be careful when using this ansatz. 28 After applying these assumptions, Eq. 1.12 simplifies to: i~ @ i (t;R(t)) @t = X j ( (T nucl +E j (t;R)) ij ~ 2 ~ d (1) ij M ~ r ) j (t;R(t)) (1.15) The next approximation involves treating the nuclear degrees of freedom classically, using the following standard quantum-classical correspondence rules. T nucl ~ 2 2M r 2 ! P 2 2M (1.16) ~ 2 ~ d (1) ij M r!i~ ~ d (1) ij ~ P M (1.17) The solution of the TD-SE is expressed on a finite set of adiabatic basis functions defined as the eigenfunctions of the stationary SE, Eq. 5.7,f i g,i = 0;:::;N b 1, whereN b is the basis size: (r;R;t) = N b 1 X i=0 c i (t) i (r;R(t)) (1.18) 7 The notation i (r;R(t)) indicates functional dependence of the electronic wave func- tions, i , on the electronic coordinates,r, and its parametric dependence on the nuclear coordinates,R. It should be noted that with this definition, there exists no explicit func- tional dependence of (r;R;t) on the nuclear positionsR. Using the ansatz, Eq. 1.18 together with Eqs. 1.16 and 1.17, one converts the TD- SE, Eq. 1.15, into the equations of motion for the amplitudes,c i (t): i~ @c i @t = N b 1 X j=0 E i ij i~ ~ d (1) ij ~ P M ! c j (1.19) Solution of Eq. 1.19 describes evolution of the electronic degrees of freedom. The first term represents evolution of the coefficient of an eigenfunction of the Hamiltonian, and the second term couples the electronic states with each other via the nuclear kinetic energy and the non-adiabatic coupling. In the Ehrenfest method, the nuclei move on the average electronic potential energy surface weighted by the corresponding populations, subject to the mean-field force: 29 ~ F mf =h (r;R;t)jr R H el (r;t;R)j (r;R;t)i = * X i c i (t) i (r;R(t))jr R H el (r;t;R)j X j c j (t) j (r;R(t)) + = * X i c i c i i jr R H el j i + * X i;j j6=i c i c j i jr R H el j j + = X i jc i j 2 ~ F i + X i;j c i c j (E i E j )d (1) ij (1.20) where ~ F i = i j @H el @R j i (1.21) 8 is the Hellmann-Feynman force corresponding to the adiabatic electronic state i . The mean-field force is computed on the basis of the contribution of each adiabatic state to the time-dependent electronic wave function. The first term couples the populations of the electronic statesjc i j 2 to the nuclear trajectory, while the second term with the coefficients c i c j includes interferences between these states. The system of coupled differential equations 1.19 and 1.20 constitute the Ehrenfest method. The trajectories generated in the Ehrenfest dynamics are fully deterministic and continuous. Although the Ehrenfest technique is simple and easy to implement, it should be men- tioned that the mean-field approximation is only adequate when either the quantum and classical subsystems are weakly coupled to each other, or when the classical subsys- tem reacts similarly to all quantum states included in the active space. The Ehrenfest approach also cannot collapse the system’s wavefunction on to a pure state. After pass- ing through a region of strong electronic coupling, the nuclear trajectory evolves on a combination of nearby adiabatic wavefunctions even when the coupling between the states becomes negligible. However, interpreting the results of the Ehrenfest dynamics as the average values of the properties and their distributions resolves this deficiency of the Ehrenfest method. The Ehrenfest scheme is fully quantum coherent and does not account for decoherence. Additionally, the standard Ehrenfest method does not satisfy detailed balance between transitions downward and upward in energy, and therefore, it cannot generate thermodynamic equilibrium between quantum and classical subsys- tems. In other words, it does not produce the Boltzmann equilibrium populations of the quantum states at a finite temperature. 9 Figure 1.1: The left panel shows stochastic hops from one surface to the other when passing through a region of strong coupling in multiple trajectories in surface hopping method. The right panel shows that the nuclear trajectory starts evolving on a mean surface in Ehrenfest dynamics after passing through the region of strong coupling. 1.2.2 Fewest Switches Surface Hopping The fewest switches surface hopping algorithm was developed by J. C. Tully in 1990 to incorporate the non-adiabatic effects by including excited adiabatic surfaces in the cal- culations, and allowing the classical subsystem to “hop” between the various coupled states. 17 At each time step, the classical subsystem evolves on a single adiabatic poten- tial energy surface following the Newton equation unlike the mean field in the Ehren- fest technique, and the quantum parts are propagated according to the time-dependent Schr odinger equation. An ensemble of trajectories is propagated to obtain meaningful physical properties. As in the Ehrenfest method, the evolution of the electronic degrees of freedom is described by the equation of motion for the amplitudes, c i (t) given by equation 1.19. The FSSH algorithm allows the system to hop between electronic states at any point in time. The hopping probability from the current electronic state,, to another state,, 10 during the time interval t is determined by the strength of the non-adiabatic coupling between the two states. It is defined as g = t b (t) a (t) (1.22) where a (t) = c (t)c (t) is the electronic population of the adiabatic state, and b (t) =2Re a _ Rd is related to the probability flux, _ a (t) = P 6= b . The switching probability,g , is computed at each integration time-step. If it turns out to be a negative number, it is reset to zero. Once computed, the switching probability is compared to a random number 2 [0; 1] to determine whether the system should remain in the current electronic state or hop to a new state. The state to which the system hops is defined by: X =1 g < +1 X =1 g (1.23) However, if P n =1 g < 1 wheren is the total number of adiabatic states included in the simulation, the system remains in the current state, and there is no hop. No- tice that the FSSH probabilities depend on the flux of state populations, and not the state populations themselves. This minimizes the number of hops, and hence the prefix “fewest switches”. Once a trajectory hops to a new state, the nuclear velocities are rescaled in the di- rection of the non-adiabatic coupling, equation 1.13a to conserve the total energy of the electon-nuclear system. If the energy available in the nuclear coordinates is insufficient to accommodate a hop to a higher energy state, the hop is rejected. Such hops are known as classically ”forbidden” hops. 11 Notice,g =g , meaning transitions up and down in energy are equally prob- able. This violates the detailed balance principle. Rejection of hops to higher energy states and velocity rescaling restores the validity of this fundamental principle. FSSH is a popular and well-tested method for the simulation of nonadiabatic pro- cesses. It has been applied to study a broad range of systems and processes including photodissociation of small molecules, 30 vibrational relaxation in clusters and condensed phase 31 electron transfer, 32 and photoinduced dynamics in organic molecules. 33 The popularity of FSSH has also brought forward its limitations, as demonstrated by the numerous types of corrections and modifications to the standard FSSH algorithm. A detailed discussion of such methods can be found in this review. 34 A major limitation of FSSH is called the trivial crossing problem. In large scale systems, such as extended polyatomic molecules, molecular crystals or assemblies of quantum dots, there may exist electronic states localized in different parts of a system and coupled very weakly. Such non-interacting states may have similar energy and cross during dynamics, leading to trivial crossings in which the NA couplings behave as delta functions of time. Finite time-step numerical propagation is likely to miss such crossings and give erroneous re- sults like unphysical long range energy transfer. Further, FSSH is known to give more accurate results in the adiabatic representation rather than diabatic representation, be- cause the velocity rescaling process differs strongly in the two representations. The rep- resentation dependence of the results may be viewed as an undesirable property. 34 Both the Ehrenfest and the FSSH methods are intrinsically limited by the classical mechani- cal description of the nuclear motion. Many nuclear quantum effects, such as tunneling and zero-point vibrational energy that are important for low temperature dynamics, can- not be easily described by classical trajectories. 35, 36 Last but not least, loss of quantum 12 coherence in the electronic subsystem induced by coupling to nuclear motions has to be addressed. 37–44 1.2.3 Decoherence Decoherence is the loss of quantum coherence when a quantum system is not perfectly isolated, but is in contact with it’s environment. It is the absence of macroscopic quan- tum effects like interference (which is a coherent process) due to interactions between a quantum system and the larger macroscopic environment. Figure 1.2: Schrodinger’s cat. Decoherence are all those phenomena that result in the cat being either dead or alive. The cat is not in a superposition of the two states anymore. For example, consider Schr odinger’s cat in a box whose life depends on the disin- tegration of a single atom. In the quantum mechanical view, the cat is alive and dead at the same time. It is in a superposition of these two quantum states. However, if we open the box to check the state of the cat, the cat is either dead or alive. The process of opening the box and other phenomena like the air molecules in the box, the movement of the cat, etc. could trigger the disintegration of the atom, that leaves the cat either dead 13 or alive. Any real quantum system is never truly isolated from these interactions with the environment. The semiclassical equations of motion derived above are not always consistent in the sense that the electronic state and the nuclear forces correspond correctly to each other. The solution of the semiclassical TD-SE (equation 1.19), c i (t), preserves coherences between electronic wavefunctions, c i c j , indefinitely. According to the fully quantum- mechanical description of the electron-nuclear system, in regions of strong derivative coupling, these coherences should decay to zero and the wavefunction should collapse to an eigenstate. In FSSH, a swarm of trajectories is propagated and, the stochastic hops ensure that the number of trajectories on each surface match the electronic wavefunction. When a trajectory enters a region of strong derivative coupling, the nuclear wavefunctions can bifurcate into more than one wavepackets on different surfaces in exact quantum dy- namics, however, in FSSH, this is not treated correctly. The electronic wavefunction is always propagated exactly along a particular nuclear trajectory, unaware that the total wave packet is now breaking apart. Due to this, the components of the electronic wave- function on different electronic surfaces always remain in phase, and thus, are coherent with each other. Neglect of decoherence is a well-known problem of FSSH. 17, 45–47 The Ehrenfest algorithm is based on the same semiclassical equations derived above from the TD-SE. Thus, it is trivial that the Ehrenfest method also does not account for decoherence effects. In a 2-level system described by the Ehrenfest method, a complete loss of coherence between the pair of electronic states is possible only if the population of one of the states is zero. This thesis consists of three parts. The first part (chapters 2 & 3) presents the ap- plication of adiabatic molecular dynamics in two distinct scenarios. In chapter 2, we 14 employed the molecular dynamics and tight-binding density functional theory to study the atomistic origin of room temperature coherence between the light-hole and heavy- hole excitations in two dimensional CdSe and CdSe/CdZnS core/shell nanoplatelets, synthetized and studied recently by optical spectroscopy. This work was done in col- laboration with Dr. Sougata Pal, Prof. Thomas Frauenheim and Prof. Oleg Prezhdo. Dr. Sougata Pal ran the molecular dynamics simulations; analysis (figures and tables) of the simulation results was performed by myself as well as Dr. Pal. The paper was written jointly by myself and Dr. Pal and edited by Dr. Oleg Prezhdo and Dr. Thomas Frauenheim. In chapter 3, we employed the highly accurate quantum chemistry method, EOM-EE-CCSD to study the mechanism of formation of He + 4 in mass spectra of large helium droplets. This work was done in collaboration with Prof. Curt Wittig, Prof. Anna Krylov, Prof. Andrey Vilesov and Prof. Oleg Prezhdo. The second part (chapter 4) presents an application of non-adiabatic molecular dy- namics to analyze the electron-vibrational energy relaxation following photoexcitation of an extended spirobifluorene molecule using FSSH-based simulations. We highlight the trivial crossing problem of FSSH and provide a comparison of two different mod- ified FSSH schemes with ad hoc corrections to the trivial crossing problem of FSSH. This work was supervised by Dr. Tammie Nelson and Dr. Sergei Tretiak at Los Alamos National Laboratory. The third part (chapter 5) presents the method development and testing of a semi- classical approach for non-adiabatic molecular dynamics based on the Ehrenfest method with corrections for decoherence effects and detailed balance. This work was done in collaboration with Dr. Joanna Jankowska and Prof. Oleg Prezhdo. The three parts of this thesis are connected by the main challenge of these projects: an accurate descrip- tion of condensed phase non-adiabatic dynamics. This study resulted in the development 15 of E-DDB, which is comparably accurate and computationally more efficient than the surface-hopping based counterpart, DISH. 16 Chapter 1 References [1] David A.C Beck and Valerie Daggett, Methods 34, 112 (2004), Investigating Protein Folding, Misfolding and Nonnative States: Experimental and Theoretical Methods. [2] Yinglong Miao, Ferran Feixas, Changsun Eun, and J. Andrew McCammon, Jour- nal of Computational Chemistry 36, 1536 (2015). [3] Anirban Hazra, Alexander V . Soudackov, and Sharon Hammes-Schiffer, The Jour- nal of Physical Chemistry Letters 2, 36 (2011), PMID: 26295211. [4] Puja Goyal and Sharon Hammes-Schiffer, The Journal of Physical Chemistry Let- ters 6, 3515 (2015), PMID: 26275870. [5] Xing Gao, Hua Geng, Qian Peng, Jiajun Ren, Yuanping Yi, Dong Wang, and Zhi- gang Shuai, The Journal of Physical Chemistry C 118, 6631 (2014). [6] Linjun Wang, Oleg V . Prezhdo, and David Beljonne, Phys. Chem. Chem. Phys. 17, 12395 (2015). [7] Jing Huang, Likai Du, Jun Wang, and Zhenggang Lan, The Journal of Physical Chemistry C 119, 7578 (2015). [8] Xing Gao, Qian Peng, Yingli Niu, Dong Wang, and Zhigang Shuai, Phys. Chem. Chem. Phys. 14, 14207 (2012). [9] S. Fernandez-Alberti, Valeria D. Kleiman, S. Tretiak, and Adrian E. Roitberg, The Journal of Physical Chemistry Letters 1, 2699 (2010). [10] Evan G. Buchanan, Joseph R. Gord, and Timothy S. Zwier, The Journal of Physical Chemistry Letters 4, 1644 (2013), PMID: 26282972. [11] Linjun Wang, Yoann Olivier, Oleg V . Prezhdo, and David Beljonne, The Journal of Physical Chemistry Letters 5, 3345 (2014), PMID: 26278443. [12] Alexey V . Akimov and Oleg V . Prezhdo, Journal of the American Chemical Soci- ety 136, 1599 (2014), PMID: 24397723. [13] Mirjana Eckert-Maksi´ c and Ivana Antol, The Journal of Physical Chemistry A 113, 12582 (2009), PMID: 19603774. 17 [14] Haiming Zhu, Ye Yang, Kim Hyeon-Deuk, Marco Califano, Nianhui Song, Youwei Wang, Wenqing Zhang, Oleg V . Prezhdo, and Tianquan Lian, Nano Letters 14, 1263 (2014), PMID: 24359156. [15] Run Long, Niall J. English, and Oleg V . Prezhdo, The Journal of Physical Chem- istry Letters 5, 2941 (2014), PMID: 26278240. [16] Jin Liu and Oleg V . Prezhdo, The Journal of Physical Chemistry Letters 6, 4463 (2015), PMID: 26505613. [17] John C. Tully, The Journal of Chemical Physics 93, 1061 (1990). [18] P. Ehrenfest, Zeitschrift fur Physik 45, 455 (1927). [19] M. Baer, Beyond Born-Oppenheimer: Electronic Nonadiabatic Coupling Terms andConicalIntersections. Wiley, 2006. [20] C.J. Cramer, EssentialsofComputationalChemistry: TheoriesandModels. Wiley, 2005. [21] Anders M. N. Niklasson, Peter Steneteg, Anders Odell, Nicolas Bock, Matt Chal- lacombe, C. J. Tymczak, Erik Holmstr¨ om, Guishan Zheng, and Valery Weber, The Journal of Chemical Physics 130, 214109 (2009). [22] Stoyan Karabunarliev, Martin Baumgarten, Eric R. Bittner, and Klaus M¨ ullen, The Journal of Chemical Physics 113, 11372 (2000). [23] Stoyan Karabunarliev, Martin Baumgarten, and Klaus M¨ ullen, The Journal of Physical Chemistry A 104, 8236 (2000). [24] F Molnar, M Ben-Nun, T.J Martinez, and K Schulten, Journal of Molecular Struc- ture: THEOCHEM 506, 169 (2000). [25] Eyal Neria and Abraham Nitzan, The Journal of Chemical Physics 99, 1109 (1993). [26] Pavel Jungwirth and R. Benny Gerber, The Journal of Chemical Physics 104, 5803 (1996). [27] Todd J. Martinez, M. Ben-Nun, and R. D. Levine, The Journal of Physical Chem- istry 100, 7884 (1996). [28] Alexey V . Akimov and Oleg V . Prezhdo, Journal of Chemical Theory and Compu- tation 10, 789 (2014), PMID: 26580053. [29] Folkmar A. Bornemann, Peter Nettesheim, and Christof Sch¨ utte, The Journal of Chemical Physics 105, 1074 (1996). 18 [30] S. Fernandez Alberti, N. Halberstadt, J. A. Beswick, and J. Echave, The Journal of Chemical Physics 109, 2844 (1998). [31] Gustavo Pierdominici-Sottile, Sebasti´ an Fern´ andez Alberti, and Juliana Palma, Chapter 8 - applications of mixed-quantum/classical trajectories to the study of nu- clear quantum effects in chemical reactions and vibrational relaxation processes, in Combining Quantum Mechanics and Molecular Mechanics. Some Recent Pro- gresses in QM/MM Methods, edited by John R. Sabin and Erkki Br¨ andas, vol- ume 59 of Advances in Quantum Chemistry, pages 247 – 282. Academic Press, 2010. [32] Soo Young Kim and Sharon Hammes-Schiffer, The Journal of Chemical Physics 119, 4389 (2003). [33] Ivana Antol, Mirjana Eckert-Maksi´ c, Mario Barbatti, and Hans Lischka, The Jour- nal of Chemical Physics 127, 234303 (2007). [34] Linjun Wang, Alexey Akimov, and Oleg V . Prezhdo, The Journal of Physical Chemistry Letters 7, 2100 (2016), PMID: 27171314. [35] Neil Shenvi, Sharani Roy, and John C. Tully, The Journal of Chemical Physics 130, 174107 (2009). [36] Tammie Nelson, Sebastian Fernandez-Alberti, Adrian E. Roitberg, and Sergei Tre- tiak, Chemical Physics Letters 590, 208 (2013). [37] Eric R. Bittner and Peter J. Rossky, The Journal of Chemical Physics 103, 8130 (1995). [38] Benjamin J. Schwartz, Eric R. Bittner, Oleg V . Prezhdo, and Peter J. Rossky, The Journal of Chemical Physics 104, 5942 (1996). [39] Oleg V . Prezhdo, The Journal of Chemical Physics 111, 8366 (1999). [40] Amber Jain, Ethan Alguire, and Joseph E. Subotnik, Journal of Chemical Theory and Computation 12, 5256 (2016), PMID: 27715036. [41] Shu Chun Cheng, Chaoyuan Zhu, Kuo Kan Liang, Sheng Hsien Lin, and Don- ald G. Truhlar, The Journal of Chemical Physics 129, 024112 (2008). [42] Ross E. Larsen, Michael J. Bedard-Hearn, and Benjamin J. Schwartz, The Journal of Physical Chemistry B 110, 20055 (2006), PMID: 17020394. [43] Heather M. Jaeger, Sean Fischer, and Oleg V . Prezhdo, The Journal of Chemical Physics 137, 22A545 (2012). 19 [44] Alexey V . Akimov, Run Long, and Oleg V . Prezhdo, The Journal of Chemical Physics 140, 194107 (2014). [45] Joseph E. Subotnik and Neil Shenvi, The Journal of Chemical Physics 134, 024105 (2011). [46] Joseph E. Subotnik and Neil Shenvi, The Journal of Chemical Physics 134, 244114 (2011). [47] Joseph E. Subotnik, The Journal of Physical Chemistry A 115, 12083 (2011), PMID: 21995423. 20 Chapter 2: Atomistic Analysis of Room Temperature Quantum Coherence in Two-Dimensional CdSe Nanostructures 2.1 Introduction An intrinsic feature of any quantum system, coherence has attracted considerable atten- tion due to its possible role in energy and charge transfer in nano-scale and biological systems. 1–15 The development of ultrafast multidimensional coherent spectroscopies has led to the discovery of time-dependent oscillatory signals in biological and man-made light harvesting systems at room-temperature. 6, 16, 17 These oscillations are interpreted as quantum beating due to a coherent superposition of electronic states with well-defined relative phase. 1 The oscillation frequency is dictated by the energy gap between the participating electronic states. The longer this oscillation remains intact, the longer is quantum coherence. 18, 19 The knowledge of the decoherence time is of fundamental importance to light harvesting, opto-electronic, spintronic, lasing, quantum computing and other applications of semiconductor nanocrystals. It sets the timescale for which the operations based on coherent light-matter interactions can occur in a given material. Typically, long coherence times are desired. 2, 6, 8, 11, 20 21 Many factors contribute to the loss of quantum mechanical phase. In bulk materials with continuous electronic spectra, the dephasing is mainly associated with quantum transitions due to carrier-carrier and carrier-phonon scattering, carrier recombination, Auger-type phenomena, etc. 21, 22 In comparison, transitions are slow in quantum sys- tems with discrete electronic spectra, and dephasing happens elastically, without change in population of the participating states. Called pure-dephasing, 23, 24 this process sets a firm time limit for which a coherent system can progress through the Rabi cycle. 25, 26 Decoherence is a consequence of the fact that a quantum system can never be completely isolated. Its interactions with the surrounding environment will eventually destroy phase correlations, leading to decoherence. Both ensemble effects (inhomogeneous broadening) and coupling to the bath (ho- mogeneous broadening) result in fast dephasing, thus contaminating and challenging the detection of electronic coherence. 27 Several reported electronic coherence studies remain controversial highlighting the challenges in experimental detection. To ratio- nalize the long-lived coherences that last for hundreds of femtoseconds, several models invoking strong correlation in the mechanism of decoherence have been proposed. It has been shown that resonant lattice vibrations (phonons) play an essential role. 11, 28, 29 Coupling with the same phonon modes results in a correlated evolution of the excitons, enabling long-lived coherences. However, electron-nuclear correlation is not the only factor that can lead to the development of long-lived coherences. A small magnitude of the fluctuation in the energy gap between the electronic states is generally needed to sustain coherence, 19 and a small energy gap can also lead to a small fluctuation. Inde- pendent bosons lead to the simplest model of excitonic dephasing. 30 The cumulant and quadratic coupling models take into account both real and virtual transitions, and are of- ten used to solve the exciton-phonon problem. 31 If available, an atomistic representation 22 of the electron-phonon system provides the most realistic description of the dephasing processes. 32–38 Colloidal semiconductor nanocrystals exhibit unique optical and electronic prop- erties that can be tuned by variation of size, shape and composition. 39–41 Among the different morphologies, two dimensional nanoplatelets (NPLs) draw particular attention because of their large absorption cross-section and narrow spectral line-widths. 42, 43 The NPL structural pattern leads to strong anisotropic quantum confinement. 42, 44 Owing to their considerably larger intrinsic absorbance, as compared to quantum dots, NPLs are more desirable for light-harvesting applications. 44–46 The absence of ensemble disorder and narrow spectral lineshapes make colloidal NPLs an ideal system to study electronic coherence and gain insight into its role in charge and energy transfer. 47 A significant body of work has been reported on these high quality nano-systems, focusing on coher- ent quantum dynamics of excitons. 1, 48, 49 The recent observation of quantum beats in colloidal CdSe and CdSe/CdZnS core/shell NPLs at room temperature was attributed to a coherent superposition of the heavy hole (HX) and light hole (LX) excitons. 48, 49 The beats were characterized using 2D electronic spectroscopy by Scholes and co-workers. 49 By synthesizing NPLs with precisely controlled thickness, the authors minimized the inhomogeneous broadening and observed dephasing times in the range of 13–24 fs for both CdSe and CdSe/CdZnS core/shell NPLs. The exact knowledge of the system size creates an excellent opportunity for an atomistic simulation, which can reproduce di- rectly the experimental setup. In addition, atomistic modeling takes into consideration interfacial interactions, core/shell interpenetration, defects and other realistic details that are inevitably present in experiment. In this contribution, we elucidate the atomistic origin of the long lived coherence between heavy and light hole excitations in CdSe and CdSe/CdZnS core/shell NPLs, 23 including both pristine systems, and NPLs with cationic and anionic defects caused by cross diffusion at the core/shell interface. 50 Based on density functional theory (DFT), our molecular dynamics simulations directly mimic the experiment of Scholes and co- workers, 49 showing excellent agreement with the measured coherence times. The simu- lations demonstrate that long coherence arises due to small fluctuations in the energy gap between the heavy and light hole excitations, as well as due to a relatively slow decay of the correlation between the heavy and light hole energies. The gap fluctuation is an order of magnitude smaller than thermal energy. Since both heavy and light holes cou- ple to the same phonon modes, their energy fluctuations can synchronize, minimizing coherence loss. The electronic structure calculations demonstrate that the excitations are localized near NPL surfaces. Hence, coupling to surface acoustic phonons constitutes the main source of decoherence, and atom exchange at a core/shell interface provides additional decoherence channels. Surprisingly, exchange of Se and S atoms that form the core and shell valence bands has little effect on the coherence time. On the contrary, exchange of Cd and Zn that contribute to the conduction bands has a major effect on the electronic coherence properties in the relevant energy range. The atomistic under- standing of the origins of the long-lived coherence and the key sources of coherence loss, established for the well-defined core/shell NPLs, should apply to other nanoscale systems. 2.2 Theoretical Methods and Simulation details Constrained by the large number of atoms in the present problem, and aiming to per- form time-domain atomistic simulation with a reasonable computational cost, we carried out all quantum-mechanical calculations, including geometry optimization, electronic 24 structure and adiabatic molecular dynamics (MD), with the self-consistent charge den- sity functional tight binding (SCC-DFTB) method 51–57 as implemented in the DFTB+ code. 51, 58 The parameter set (Slater-Koster files) used in this calculation were exten- sively tested for a variety of Cd-chalcogenide systems. 59 The details of the methods can be found elsewhere. 51–57 The simulations were performed in a cell that was peri- odically replicated along the X and Y directions. 40 ˚ A of vacuum space was imposed between the replicas in the Z direction, resulting in 2D NPLs. The optimized cell lengths were 8.54 ˚ A and 13.92 ˚ A in the X and Y direction, respectively. The structures were fully optimized at 0 K and then heated to 300 K with repeated velocity rescaling. One picosecond micro-canonical MD trajectories were generated for each NPLs using the Verlet algorithm with a 1 fs time step and Hellman-Feynman forces. The CdSe NPLs were five (5) unit layers thick. Each layer contained (2x2) unit cells. The total number of Cd and Se atoms in the simulation cell was 160. The opti- mized structure is shown in Figure 5.1(a), in which the bracket indicates one unit layer thickness. The core part of the CdSe/CdZnS core/shell NPLs was exactly the same. Two (2) unit layers of the CdZnS shell were added on each side of the core, resulting in a total of 9 unit layers, and 288 Cd, Se Zn and S atoms in the simulation cell. Figure 5.1(b) rep- resents the optimized CdSe/CdZnS core/shell NPLs whereas, Figures 5.1(c) and 5.1(d) represent the optimized structures corresponding to the core/shell NPLs with anionic defects (Se/S atoms exchanged) and cationic defects (Cd/Zn atoms exchanged) at the interface, respectively. The arrows indicate each Se/S or Cd/Zn atom pair exchanged at the interface between core and shell. In all cases, we chose the stoichiometric 1010 nonpolar surface derived from the wurtzite bulk, and the NPLs were constructed in ac- cordance with the experimental work by Scholes and co-workers. 49 25 Figure 2.1: Optimized geometries of (a) CdSe nanoplatelet, (b) CdSe/ CdZnS core/shell nanoplatelet, (c) CdSe/CdZnS core/shell nano-platelet with S/Se ex- change, and (d) CdSe/CdZnS core/shell nanoplatelet with Zn/Cd exchange. The bracket in part a identifies a unit layer. The blue, magenta, yellow, and red balls represent Cd, Se, Zn, and S atoms, respectively. The S/Se and Zn/Cd exchanges are shown in parts c and d by arrows. The coherence times between pairs of optical excitations in the NPLs are character- ized by computing and analyzing the pure-dephasing times calculated using the optical response theory. 30 The phonon-induced fluctuations in the energy gap between the states are characterized by the autocorrelation function (ACF). C(t) = E(t)E(0) (2.1) 26 The brackets indicate canonical averaging. The ACF is normalized C norm (t) = E(t)E(0) E(0) 2 (2.2) by its initial value C(0) = E(0) 2 . The square root of this value gives the aver- age fluctuation of the energy gap. The pure-dephasing function is computed using the second-order cumulant expansion to the optical response function. 30 D cumu (t) =exp(g(t)) (2.3) whereg(t) is given by g(t) = t Z 0 d 1 1 Z 0 d 2 C( 2 ) (2.4) The electronic energy gaps used for the computation of the coherence times were ob- tained using the Kohn-Sham orbital energies from the ground state DFTB calculations. The optical spectra were also computed in the Kohn-Sham representation to correspond to the energy levels used in the coherence time calculations. 2.3 Results and discussion Depending on the core and shell thickness, a CdSe/CdZnS hetero-structure shows type- I, type-II or quasi type-II band alignment. 60, 61 The partial density of states (PDOS) of the optimized structures of the CdSe, CdSe/CdZnS core/shell NPLs are given in Figures 2.2(a) and 2.2(b), respectively. The PDOS plot for the CdSe NPLs shows that Se atoms contribute to the valence band maximum (VBM), whereas contribution to the conduction band minimum (CBM) 27 Figure 2.2: Projected density of states (DOS) of (a) CdSe and (b) CdSe/CdZnS core/shell nanoplatelets. For CdSe, the blue line indicates the total DOS, whereas the orange and magenta lines represent the projected DOS of Cd and Se atoms, respectively. For the core/shell platelet, the black line represents the total DOS, whereas the blue and red lines indicate the core and shell DOS, respectively. mainly comes from Cd atoms. The VBM of the core/shell structure arises from the core, while CBM stems from the shell resulting in an overall type-II band alignment. Table 2.1: Canonically Averaged Energy, E of Heavy (HX) and Light (LX) Hole Excitations in the CdSe Nanoplatelet a no. HX E (eV) LX E (eV) E exp (eV) E (eV) (E) 2 (LX-HX) (LX-HX) (meV 2 ) (LX-HX) 1 HX1 2.37 LX1 2.49 0.14 0.12 2.22 2 LX2 2.48 0.11 3 LX3 2.56 0.19 4 HX2 2.41 LX1 2.49 0.07 5 LX2 2.48 0.08 6 LX3 2.56 0.15 2.58 a Experimental, E exp , and calculated, E energy gaps between heavy and light hole excitations, and square gap fluctuation(E) 2 E E 2 . The room-temperature optical spectrum of CdSe NPLs is shown in the insert of Figure 2.3(a), in which the dotted red line indicates individual transitions, while the 28 Figure 2.3: (a) Room temperature pure-dephasing functions for pairs of heavy and light hole excitations in the CdSe nanoplatelet. The red and blue lines represent the HX1-LX1 and HX2-LX3 excitation pairs, whose energy differences most closely match the experiment, Table 5.1. The inset in panel a displays the optical spec- trum. The individual transitions shown by the red line are Gaussian-broadened to obtain the black line. (b) Autocorrelation functions of the energy gaps between the same pairs of excitations. (c) Fourier transforms of (i) HX1, (ii) HX2, (iii) LX1, and (iv) LX3 excitation energies. black curve corresponds to Gaussian broadened peaks. Following the notation of the experimental paper, 49 two different heavy hole excitons (HX), labeled as HX1, HX2 29 with the average energy of 2.37 and 2.41 eV and three different light hole excitons (LX), labeled as LX1, LX2, LX3 with the average energy of 2.46, 2.48 and 2.56 eV , are identified. Table 5.1 summarizes different HX and LX transition energies for the CdSe NPLs. To determine the pure-dephasing time between HX and LX, a total of six pairs of excitons, namely HX1-LX1, HX1-LX2, HX1-LX3, HX2-LX1, HX2-LX2, HX2-LX3 have to be considered. The average energy differences between HX and LX excitons range from 0.34 to 0.15 eV , Table 5.1. We focus on the two pairs, HX1-LX1 and HX2- LX3, that are closest to the experimental data. The energy gaps between these pairs are within 15% from the experimental value of 0.14 eV . 49 The pure-dephasing functions of these two exciton pairs are shown in Figure 2.3(a). The room temperature optical spectrum for the CdSe/CdZnS core/shell hetero- structure is shown in the insert of Figure 2.4(a). Here, there is only one significant transition in the HX peak, labeled as HX1 with an average energy of 2 eV , and three different LX excitons, labeled as LX1, LX2 and LX3 with the average energies of 2.16, 2.20, and 2.23 eV , respectively. Table 2.2: Canonically Averaged Heavy (HX) and Light (LX) Hole Excitations in the CdSe/CdZnS Core/Shell Nanoplatelet no. HX E (eV) LX E (eV) E exp (eV) E (eV) (E) 2 (LX-HX) (LX-HX) (meV 2 ) (LX-HX) 1 HX1 2.00 LX1 2.16 0.18 0.16 1.89 2 LX2 2.20 0.20 5.44 3 LX3 2.23 0.23 Table 2.2 summarizes different HX and LX optical transition energies for the CdSe/CdZnS core/shell NPLs and shows the canonically averaged energy difference 30 Figure 2.4: (a) Room temperature pure-dephasing functions for pairs of heavy and light hole excitations in the CdSe/CdZnS core/shell nanoplatelet. The red and blue lines represent the HX1-LX1 and HX1-LX2 excitation pairs, whose energy differences most closely match the experiment, Table 2.2. The inset in panel a displays the optical spectrum. The individual transitions shown by the red line are Gaussian-broadened to obtain the black line. (b) Autocorrelation functions of the energy gaps between the same pairs of excitations. (c) Fourier transforms of (i) HX1, (ii) LX1, and (iv) LX2 excitation energies. between various excitons pairs. A total of three pairs of HX-LX excitons, namely HX1- LX1, HX1-LX2 and HX1-LX3 are possible, for which the average energy differences 31 are 0.16, 0.20 and 0.23 eV , respectively. The experimental HX-LX average energy dif- ference is found to be 0.18 eV . Hence, we choose the two pairs of excitons, HX1-LX1 and HX1-LX2, as their energy difference is within 15% of the experimental value. The dephasing functions of these pairs of excitons are shown in 2.4(a). The data of Tables 5.1 and 2.2 indicate that the HX and LX excitation energies are lower for the CdSe/CdZnS core/shell NPLs compared to the CdSe NPLs, while the HX- LX splitting is increased. This is in excellent agreement with the experimental results of Scholes and co-workers. 49 This red-shift of the spectra for the core/shell system can be explained by the type-II alignment of the core and shell energy levels, and the primary contribution of the core to the VBM, and the shell to the CBM, Figure 2.2. The room temperature pure-dephasing functions for the selected pairs of excitons in the CdSe and CdSe/CdZnS core/shell NPLs are plotted in Figures 2.3(a) and 2.4(a), respectively. The dephasing times are computed by fitting to the Gaussian, f(t) = Aexp(0:5( t 1 ) 2 ). The dephasing times 18.51 and 17.27 fs computed for the CdSe NPL are in excellent agreement with the experimentally measured time of 163 fs. 49 The dephasing times computed for the CdSe/CdZnS core/shell NPL are 25.5 fs and 11.57 fs, which is also consistent with the experimental value of 243 fs. 49 In addition to calculating the dephasing times for individual pairs of transitions in the HX and LX peaks, we also calculated the dephasing times by averaging the energies all transitions in the HX and LX peaks. The values are 34.31 fs for the CdSe NPL and 30.86 fs for the CdSe/CdZnS core/shell NPL, Figure 2.5(a) upper and middle panel, respectively. According to the second-order cumulant approximation analysis, equations 5.26 and 2.4, the coherence time is determined by the unnormalized ACF, equation 5.25, for the energy gap between the pairs of states forming a coherence superposition. In turn, the unnormalized ACF can be analyzed in terms of its initial value, C(0) = E(0) 2 , 32 Figure 2.5: (a) Room temperature pure-dephasing functions obtained by aver- aging over all pairs of heavy and light hole excitations in CdSe (upper panel), CdSe/CdZnS (middle panel), and CdSe/CdZnS core/shell nanoplatelets with S/Se exchange at the interface (lower panel). (b) Room temperature optical spectrum of the CdSe/CdZnS core/shell nanoplatelet with Zn/Cd exchange at the interface. The individual transitions shown by the red line are Gaussian- broadened to obtain the black line. corresponding to the energy gap fluctuation, and the timescale of decay of the normal- ized ACF, equation 2.2. Coherence is long-lived if the energy gap fluctuation is small and/or if the decay of the normalized ACF is slow. 19 The energy gap fluctuations for the core and core/shell NPLs are given in the last columns of Tables 5.1 and 2.2, while the normalized ACFs are shown in Figures 2.3(b) and 2.4(b). Both factors contribute to the coherences observed in the systems under investigation. The energy gap fluctuations are small, on the order of 1 meV . This is much smaller than thermal energy at room temper- ature, 25 meV . The ACFs decay on the same timescale as the pure-dephasing functions, and exhibit a series of recurrences. The faster decay of the pure-dephasing function for 33 the HX1-LX2 exciton pair, compared to the HX1-LX1 pair in the CdSe/CdZnS NPLs, Figure 2.4, is due to the larger gap fluctuation, Table 2.3, rather than faster ACF decay. In general, the ACFs decay faster in the core/shell hetero-structure than in the core NPL, because a broader range of phonon modes is available in the system composed of more types of atoms. Table 2.3: Canonically Averaged Heavy (HX) and Light (LX) Hole Excitations in the CdSe/CdZnS Core/Shell Nanoplatelets with the S/Se Exchange at the Interface no. HX E (eV) LX E (eV) E (eV) (E) 2 (LX-HX) (meV 2 ) (LX-HX) 1 HX1 1.98 LX1 2.14 0.16 2.26 2 LX2 2.26 0.28 5.68 3 HX2 2.03 LX1 2.14 0.11 1.50 4 LX2 2.26 0.23 5.13 5 HX3 2.10 LX1 2.14 0.04 5.31 6 LX2 2.26 0.16 2.26 In order to identify the phonon modes that couple to the electronic subsystems and induce coherence loss, we computed Fourier transforms (FT) of the excitation energies of the states forming coherent superpositions, Figures 2.3(c) and 2.4(c). These influ- ence (power) spectra show more frequencies for the CdSe/CdZnS NPL than for CdSe. This should be expected, because the core/shell system has a broader range of vibrations, and since the VBM and the CBM reside on core and shell, respectively, such that the ex- citations involve both subsystems. Particularly important for the long-lived coherence, the influence spectra for the HX and LX excitations reveal very similar frequencies. It is true for both CdSe and CdSe/CdZnS NPLs. By coupling to the same phonon modes, the HX and LX excitations exhibit similar fluctuations that maintain correlations for a fairly long time. 34 The phonon peak appearing below 100 cm 1 can be attributed to the discrete longi- tudinal acoustic (LA) mode as well the surface acoustic (SA) mode, whereas the mode around 200 cm 1 is due to the longitudinal optical (LO) phonon of CdSe. 28, 62, 63 The surface optical (SO) phonon modes together with the LO modes that have a signature in the range of 350-400 cm 1 also couple with these excitons. 28, 62, 63 These phonon modes are associated with the surface atoms of the NPLs, indicating that the electronic states have contributions from the surface. The frequencies appearing around 400 cm 1 correspond the first overtone of the LO modes. 28, 62, 63 Figure 2.6 presents charge density distribution of the key orbitals of the CdSe and CdSe/CdZnS NPLs contributing to the HX and LX excitations. The charge densities are mainly localized on the NPL surfaces and extended towards the core in a few cases. The orbital localization rationalizes why the major contribution to the loss of quantum coherence in the electronic subsystems arises from the surface phonon modes with a few bulk modes contributing as well. Given the successful interpretation of the decoherence mechanism of the CdSe and CdSe/CdZnS NPLs, we extended our calculations to structures containing interfacial defects. The primary goal was to understand whether defects, that are inevitably present in experimental systems, have a strong influence on the electronic coherence. Cross diffusion of atoms across interface is common in core/shell materials. We consider two complementary situations, cationic exchange and anionic exchange. The former involves exchange of Cd and Zn atoms, whereas the latter is the exchange of S and Se atoms. The optimized structures for the core/shell systems with the anionic and cationic exchange are shown in Figures 5.1(c) and 5.1(d), respectively. Calculation of the room temperature spectrum of the cationic defect structure, Figure 2.5(b) reveals that such de- fects have a strong influence on the optical properties. The obtained spectrum no longer 35 36 Figure 2.6: Plots of charge densities of different electronic states that are involved in the transitions defined in Tables 5.1, 2.2 and 2.3 for (a) CdSe, (b) CdSe/CdZnS and (c) CdSe/CdZnS with S/Se exchange, respectively, at 300 K. VBM and CBM refer to valence band maximum and conduction band minimum. For CdSe, the energy of the VB 0 , VB 00 , VB 000 , CB 0 and CB 00 states correspond to VBM-0.08 eV , VBM-0.15 eV , VBM-0.19 eV , CBM+0.12 eV and CBM+0.26 eV , respectively. The HX1, HX2, LX1 and LX3 excitations arise from the following electronic transitions, VB 00 ! CBM, VB 000 ! CBM, VB 00 ! CBM 0 and VB 0 ! CBM 00 . For CdSe/CdZnS, VB 0 , VB 00 and CB 0 represent states with the energies VBM-0.04 eV , VBM-0.19 eV and CBM+0.19 eV , respectively. The HX1, LX1 and LX2 excitations correspond to the VB 0 ! CBM, VB 00 ! CBM, VB 0 ! CBM 0 transitions, respectively. For CdSe/CdZnS with S/Se exchange, VB 0 , VB 00 , VB 000 and CB 0 represent states with the energies VBM-0.05 eV , VBM-0.2 eV and CBM+0.10 eV , respectively. The HX1, HX2, HX3, LX1 and LX2 excitations correspond to the following transitions VB 0 ! CBM, VB 00 ! CBM, VB 0 ! CBM 0 , VB 000 ! CBM and VB 000 ! CBM 0 , respectively. shows well resolved peaks in the energy range of interest. The cationic atoms (Cd, Zn) are larger than the anionic atoms (Se, S). As a result, the geometric structure is more sensitive to defects in the larger cationic sub-lattice, and a stronger perturbation in the geometric structure induced by the cationic exchange leads to a stronger perturbation in the optical properties. Thus, the cationic defects are not desirable for sustained excitonic coherence in the CdSe/CdZnS NPLs for the considered excitation energy range. The scenario is completely different for the anionic defect structure. Both excitonic bands are retained and remain well resolved, Figure 2.7(a) insert. Table 2.3 lists dif- ferent HX and LX transitions that are involved in this particular case, along with the corresponding canonically averaged energies, energy differences and their fluctuations. There are three different HX excitons HX1, HX2, HX3 with the average energies of 1.98, 2.03 and 2.10 eV , respectively, and two different LX excitons, LX1 and LX2 and the energies of 2.14 and 2.26 eV . The pure-dephasing functions shown in Figure 2.7(a) demonstrate coherence times within the 11-21 fs range, which only slightly faster than 37 Figure 2.7: (a) Room temperature pure-dephasing functions for pairs of heavy and light hole excitations in the CdSe/CdZnS core/shell nanoplatelet with S/Se ex- change at the interface. Red, blue, magenta, orange, light green, and dashed dark green lines represent HX1-LX1, HX1-LX2, HX2-LX1, HX2-LX2, HX3-LX1, and HX3-LX2 excitation pairs, respectively, see Table 2.3. The inset in panel a dis- plays the optical spectrum. The individual transitions shown by the red line are Gaussian-broadened to obtain the black line. (b) Autocorrelation functions of the energy gaps between the same pairs of excitations, maintaining the color code used for the dephasing functions. (c) Fourier transforms of (i) HX1, (ii) HX2, (iii) HX3, (iv) LX1, and (v) LX2 excitation energies. 38 the 14-25 fs coherence times for the pristine core/shell system, Figure 2.4(a). Hence, the anionic defect has only a minor influence on the electronic coherence within the electronic subsystem. The small influence of the anion exchange at the CdSe/CdZnS interface on the co- herence time is supported by comparison of the energy gap fluctuations, ACF decay, phonon influence spectra and orbital localization to those of the pristine system. The gap fluctuations are still within 1-5 meV , Table 2.3, much less than the thermal energy of 25 meV . The ACFs computed for all excitation pairs in the defected hetero-structure, Figure 2.7(b), decay similarly to those for the pristine core/shell system, Figure 2.4(b). FTs of the fluctuations of the HX1, HX2, HX3, LX2, LX3 electronic transitions, Figure 2.7(c), exhibit frequencies that match those of Figure 2.4(c). The influence spectrum for the LX1 and LX2 excitations is slightly shifted to higher frequencies, panels (iv) and (v) of Figure 2.7(c). The low-frequency phonons below 100 cm 1 correspond to surface acoustic modes, while the phonons peaks around 150-200 cm 1 are due to sur- face optical modes. Introduction of defects disrupts the longer-range acoustic motions, decreasing the low frequency spectral response. At the same time, defect modes and optical modes have similar frequencies, and therefore, the spectral signal in the optical mode range is enhanced. The stronger coupling of the higher frequency phonons to the HX and LX transitions in the CdSe/CdZnS NPL with the S/Se atom exchange explain why decoherence proceeds 10-15% faster compared to the pristine system. Since loss of coherence between the HX and LX transitions is induced primarily by surface phonons, modifications of the bulk core region of the nanoplatelets should not significantly reduce the coherence time. Our calculations suggest that it is particularly important to control defects at nanoplatelet surfaces and core/shell interfaces. 39 2.4 Conclusions To summarize, we have established the atomistic origin of room temperature coher- ence between the light-hole and heavy-hole excitations in two dimensional CdSe and CdSe/CdZnS core/shell nanoplatelets, synthetized and studied recently by optical spec- troscopy. The molecular dynamics and tight-binding density functional theory calcu- lations have allowed us to mimic closely the experimental systems and measurements, and to elucidate the role of defects that are inevitably present in experiment. The CdSe core thickness of the simulated systems matches the experiment exactly, while the shell is smaller in the simulation. Figure 2.8: Schematic showing the coupling of the HX and LX excitons to the same phonon modes leading to long quantum coherence timescales. The computed 18 fs coherence time in the CdSe platelet is in excellent agreement the experimental 163 fs value. The 14-25 fs coherence times computed for the core/shell system agrees well with the measured 243 fs time. The relatively long coherence is 40 maintained due to small energy gap fluctuations, which are an order of magnitude less than thermal energy, and because of a slow decay of correlation between the heavy and light hole energies. Two different scenarios emerge when defects are introduced at the core/shell interface. Anionic defects prove harmless, and a coherence lifetime of the same order of magnitude is obtained. However, the system with cationic defects does not show well resolved peaks in the energy range of interest, indicating that it would be hard to obtain sustained coherence in this case. Coherence is maintained due to correlated evolution of the heavy-hole and light-hole excitons that couple to the same phonon modes, including surface acoustic modes below 100 cm 1 , surface optical at around 180 cm 1 , and to a smaller extent, higher frequency modes in the 350-400 cm 1 range. The charge densities of the key electronic orbitals reveal their localization on nano-system surfaces, rationalizing coupling to the surface phonons. The fundamental principles underlying coherence loss in the CdSe and CdSe/CdZnS nanoplatelets can be used to design other nanoscale systems with long-lived quantum coherence. 41 Chapter 2 References [1] Daniel B. Turner, Yasser Hassan, and Gregory D. Scholes, Nano Letters 12, 880 (2012), PMID: 22201519. [2] Gregory S Engel, Tessa R Calhoun, Elizabeth Read, Tae-Kyu Ahn, Tomas Man- cal, Yuan-Chung Cheng, Robert Blankenship, and Graham R Fleming, 446, 782 (2007). [3] Pascal Gehring, Hatef Sadeghi, Sara Sangtarash, Chit Siong Lau, Junjie Liu, Arzhang Ardavan, Jamie H. Warner, Colin J. Lambert, G. Andrew. D. Briggs, and Jan A. Mol, Nano Letters 16, 4210 (2016), PMID: 27295198. [4] Yuan-Chung Cheng and Graham R. Fleming, The Journal of Physical Chemistry A 112, 4254 (2008), PMID: 18376878. [5] Shuo Dong, Dhara Trivedi, Sabyasachi Chakrabortty, Takayoshi Kobayashi, Yinthai Chan, Oleg V . Prezhdo, and Zhi-Heng Loh, Nano Letters 15, 6875 (2015), PMID: 26359970. [6] Elisabetta Collini, Cathy Y Wong, Krystyna Wilk, Paul M G Curmi, Paul Brumer, and Gregory D Scholes, 463, 644 (2010). [7] Kai Hao, Lixiang Xu, Philipp Nagler, Akshay Singh, Kha Tran, Chandriker Kavir Dass, Christian Sch¨ uller, Tobias Korn, Xiaoqin Li, and Galan Moody, Nano Letters 16, 5109 (2016), PMID: 27428509. [8] Gitt Panitchayangkoon, Dugan Hayes, Kelly A. Fransted, Justin R. Caram, Elad Harel, Jianzhong Wen, Robert E. Blankenship, and Gregory S. Engel, Proceedings of the National Academy of Sciences 107, 12766 (2010). [9] Muhandis Shiddiq, Dorsa Komijani, Yan Duan, Alejandro Gaita-Ari˜ no, Eugenio Coronado, and Stephen Hill, 531, 348 (2016). [10] Andreas T Haedler, Klaus Kreger, Abey Issac, Bernd Wittmann, Milan Kivala, Na- talie Hammer, Jurgen Kohler, Hans-Werner Schmidt, and Richard Hildner, Nature 523, 196 (2015). [11] Daniel B. Turner, Krystyna E. Wilk, Paul M. G. Curmi, and Gregory D. Scholes, The Journal of Physical Chemistry Letters 2, 1904 (2011). [12] Alberto Torres, Robson S. Oliboni, and Luis G. C. Rego, The Journal of Physical Chemistry Letters 6, 4927 (2015), PMID: 26606950. 42 [13] Run Long and Oleg V . Prezhdo, Nano Letters 16, 1996 (2016), PMID: 26882202. [14] Manjin Zhong, Morgan P Hedges, Rose L Ahlefeldt, John G Bartholomew, Sarah E Beavan, Sven M Wittig, Jevon J Longdell, and Matthew J Sellars, Nature 517, 177 (2015). [15] Hong-Guang Duan, R. J. Dwayne Miller, and Michael Thorwart, The Journal of Physical Chemistry Letters 7, 3491 (2016), PMID: 27547995. [16] Sarah Maria Falke, Carlo Andrea Rozzi, Daniele Brida, Margherita Maiuri, Michele Amato, Ephraim Sommer, Antonietta De Sio, Angel Rubio, Giulio Cerullo, Elisa Molinari, and Christoph Lienau, Science 344, 1001 (2014). [17] Carlo Andrea Rozzi, Sarah Maria Falke, Nicola Spallanzani, Angel Rubio, Elisa Molinari, Daniele Brida, Margherita Maiuri, Giulio Cerullo, Heiko Schramm, Jens Christoffers, and Christoph Lienau, Nat Commun 4, 1602 (2013). [18] Ki-Hee Song, Munui Gu, Min-Seok Kim, Hyeok-Jun Kwon, Hanju Rhee, Hogyu Han, and Minhaeng Cho, The Journal of Physical Chemistry Letters 6, 4314 (2015), PMID: 26722967. [19] Alexey V . Akimov and Oleg V . Prezhdo, The Journal of Physical Chemistry Letters 4, 3857 (2013). [20] Niklas Christensson, Harald F. Kauffmann, T˜ onu Pullerits, and Tom´ aˇ s Manˇ cal, The Journal of Physical Chemistry B 116, 7449 (2012), PMID: 22642682. [21] Jin Liu, Amanda J. Neukirch, and Oleg V . Prezhdo, The Journal of Chemical Physics 139, 164303 (2013). [22] E. A. Muljarov, T. Takagahara, and R. Zimmermann, Phys. Rev. Lett. 95, 177405 (2005). [23] J. L. Skinner and D. Hsu, The Journal of Physical Chemistry 90, 4931 (1986). [24] D. Hsu and J.L. Skinner, Journal of Luminescence 37, 331 (1987). [25] Ranojoy Bose, Tao Cai, Kaushik Roy-Choudhury, Glenn S. Solomon, and Edo Waks, 8 (2014). [26] J. Bochmann, M. M¨ ucke, G. Langfahl-Klabes, C. Erbel, B. Weber, H. P. Specht, D. L. Moehring, and G. Rempe, Phys. Rev. Lett. 101, 223601 (2008). [27] A. Tokmakoff and M. D. Fayer, The Journal of Chemical Physics 103, 2810 (1995). 43 [28] Justin R. Caram, Haibin Zheng, Peter D. Dahlberg, Brian S. Rolczynski, Gra- ham B. Griffin, Andrew F. Fidler, Dmitriy S. Dolzhnikov, Dmitri V . Talapin, and Gregory S. Engel, The Journal of Physical Chemistry Letters 5, 196 (2014). [29] A. W. Chin, J. Prior, R. Rosenbach, F. Caycedo-Soler, S. F. Huelga, and M. B. Plenio, Nature Physics 9, 113 (2013). [30] S. (Shaul) Mukamel, Principles of nonlinear optical spectroscopy. Oxford : Ox- ford University Press, 1999. [31] P. Borri, W. Langbein, U. Woggon, V . Stavarache, D. Reuter, and A. D. Wieck, Phys. Rev. B 71, 115328 (2005). [32] Zhenyu Guo, Bradley F. Habenicht, Wan-Zhen Liang, and Oleg V . Prezhdo, Phys. Rev. B 81, 125415 (2010). [33] Bradley F. Habenicht, Hideyuki Kamisaka, Koichi Yamashita, and Oleg V . Prezhdo, Nano Letters 7, 3260 (2007), PMID: 17949045. [34] Bradley F. Habenicht, Oleg N. Kalugin, and Oleg V . Prezhdo, Nano Letters 8, 2510 (2008), PMID: 18646832. [35] Hideyuki Kamisaka, Svetlana V . Kilina, Koichi Yamashita, and Oleg V . Prezhdo, Nano Letters 6, 2295 (2006), PMID: 17034100. [36] Hideyuki Kamisaka, Svetlana V . Kilina, Koichi Yamashita, and Oleg V . Prezhdo, The Journal of Physical Chemistry C 112, 7800 (2008). [37] Jin Liu, Svetlana V . Kilina, Sergei Tretiak, and Oleg V . Prezhdo, ACS Nano 9, 9106 (2015), PMID: 26284384. [38] Angeline B. Madrid, Kim Hyeon-Deuk, Bradley F. Habenicht, and Oleg V . Prezhdo, ACS Nano 3, 2487 (2009), PMID: 19722505. [39] C. B. Murray, D. J. Norris, and M. G. Bawendi, Journal of the American Chemical Society 115, 8706 (1993). [40] Liberato Manna, Delia J Milliron, Andreas Meisel, Erik C Scher, and A Paul Alivisatos, Nat Mater 2, 382 (2003). [41] X. Michalet, F. F. Pinaud, L. A. Bentolila, J. M. Tsay, S. Doose, J. J. Li, G. Sun- daresan, A. M. Wu, S. S. Gambhir, and S. Weiss, Science 307, 538 (2005). [42] Sandrine Ithurria, Mickael Tessier, Benoit Mahler, R P S M Lobo, B Dubertret, and Alexander Efros, 10, 936 (2011). [43] Gregory D Scholes, 10, 906 (2011). 44 [44] Alexander W. Achtstein, Artsiom Antanovich, Anatol Prudnikau, Riccardo Scott, Ulrike Woggon, and Mikhail Artemyev, The Journal of Physical Chemistry C 119, 20156 (2015). [45] M Aagesen and C.B. Sorensen, page 109 (2008). [46] Bingqiang Cao and Weiping Cai, The Journal of Physical Chemistry C 112, 680 (2008). [47] Zhen Wan and Gaoke Zhang, J. Mater. Chem. A 3, 16737 (2015). [48] Micka¨ el D. Tessier, Cl´ ementine Javaux, Ivan Maksimovic, Vincent Loriette, and Benoit Dubertret, ACS Nano 6, 6751 (2012), PMID: 22783952. [49] Elsa Cassette, Ryan Pensack, Benoit Mahler, and Gregory D. Scholes, (2015). [50] Quan Ma, Miguel Isarraraz, Chen S. Wang, Edwin Preciado, Velveth Klee, Sarah Bobek, Koichi Yamaguchi, Emily Li, Patrick Michael Odenthal, Ariana Nguyen, David Barroso, Dezheng Sun, Gretel von Son Palacio, Michael Gomez, Andrew Nguyen, Duy Le, Greg Pawin, John Mann, Tony. F. Heinz, Talat Shahnaz Rahman, and Ludwig Bartels, ACS Nano 8, 4672 (2014), PMID: 24684434. [51] B. Aradi, B. Hourahine, and Th. Frauenheim, The Journal of Physical Chemistry A 111, 5678 (2007), PMID: 17567110. [52] M. Elstner, The Journal of Physical Chemistry A 111, 5614 (2007), PMID: 17564420. [53] M. Elstner, D. Porezag, G. Jungnickel, J. Elsner, M. Haugk, Th. Frauenheim, S. Suhai, and G. Seifert, Phys. Rev. B 58, 7260 (1998). [54] Thomas Kr¨ uger, Marcus Elstner, Peter Schiffels, and Thomas Frauenheim, The Journal of Chemical Physics 122, 114110 (2005). [55] T. A. Niehaus, S. Suhai, F. Della Sala, P. Lugli, M. Elstner, G. Seifert, and Th. Frauenheim, Phys. Rev. B 63, 085108 (2001). [56] D. Porezag, Th. Frauenheim, Th. K¨ ohler, G. Seifert, and R. Kaschner, Phys. Rev. B 51, 12947 (1995). [57] G. Seifert, The Journal of Physical Chemistry A 111, 5609 (2007), PMID: 17439198. [58] Thomas Frauenheim, Gotthard Seifert, Marcus Elstner, Thomas Niehaus, Christof K¨ ohler, Marc Amkreutz, Michael Sternberg, Zolt´ an Hajnal, Aldo Di Carlo, and S´ andor Suhai, Journal of Physics: Condensed Matter 14, 3015 (2002). 45 [59] Sunandan Sarkar, Sougata Pal, Pranab Sarkar, A. L. Rosa, and Th. Frauenheim, Journal of Chemical Theory and Computation 7, 2262 (2011), PMID: 26606495. [60] Haiming Zhu and Tianquan Lian, Energy Environ. Sci. 5, 9406 (2012). [61] Hagai Eshet, Michael Gr¨ unwald, and Eran Rabani, Nano Letters 13, 5880 (2013), PMID: 24215466. [62] Peng Han and Gabriel Bester, Phys. Rev. B 85, 041306 (2012). [63] V olodymyr Dzhagan, Mykhailo Valakh, Nikolai Mel’nik, Olexandra Rayevska, Irina Lokteva, Joanna Kolny-Olesiak, and Dietrich R. T. Zahn, International Jour- nal of Spectroscopy 2012, 6 (2012). 46 Chapter 3: Conversion of He(2 3 S) to He 2 a 3 + u in Helium droplets 3.1 Introduction Superfluid helium nanodroplets serve as a unique platform upon which a broad range of interesting and at times intriguing phenomena have been discovered and examined. 1–6 Indeed, many of the experimental results that have accrued throughout the past sev- eral decades have proven seminal. They have often been rationalized using qualitative models, with the understanding that detailed mechanistic pictures are likely more com- plicated, or perhaps even deviate significantly from the original models. The experimental studies that motivated this work used electron impact to create the lowest energy triplet, He (2 3 S) within helium droplets. 7 Hereafter, He and He 2 (a 3 + u ) within and on bulk helium and helium droplets shall be referred to as He and He 2 respectively. They do not autoionize or dissociate, have long radiative lifetimes, hence, they are often referred to as metastable. Of course, they can interact rapidly with atoms and molecules. The story began with a mass spectrometric study of large superfluid helium droplets by Buchenau et al. 8, 9 Plots of He + n signal intensity versus n displayed an anomaly when 47 electron kinetic energy (eKE) exceeded roughly 40 eV . Rather than decaying monotoni- cally with n, as in small droplets, 8, 10–15 the He + 4 signal was an order of magnitude more intense than expected. The authors concluded that pairs of He atoms created within large droplets are the progenitors of enhanced He + 4 production, as the energy of isolated He is 19.82 eV . They proposed that an incident electron creates He , and then goes on to create a second He within the same droplet, and that He might convert to He 2 on droplet surfaces. Sch obel et al. established the threshold by showing that the slope of He + 4 yield versus eKE increases sharply at 41 eV . 13 The mechanism whereby He converts to He 2 in liquid helium has remained elusive, in large part because of an association barrier of about 500 cm 1 for the gas phase counterpart. We shall explain the He + He! He 2 mechanism by applying a computa- tional approach that combines high level ab initio theory with molecular dynamics (MD) simulations. It enables us to analyze in detail small model systems consisting of a few He atoms, thereby obtaining insights that are germane to larger systems, notably helium droplets. Earlier, Keto et al. 16–18 and Hill et al. 19 had shown that He 2 and He are produced when superfluid helium is bombarded with 160 keV electrons. The He 2 was removed by He 2 + He 2 ! products (though products were not identified). Keto et al. showed that He was removed by bulk helium with a temperature independent 15s exponen- tial decay. 16–18 They also assigned two low resolution spectral features to high a 3 + u vibrational levels (Figure 1), and suggested that such high- levels might be due to He combining with He. Each 160 keV electron resulted in 500 He atoms, and about the same number of He 2 molecules, via cascades, which made interpretation challenging. Brooks et al. confirmed the presence of high vibrational levels ( = 10 12) following bombardment of dense cryogenic gas with 6.5 MeV protons. 20 Northby and coworkers 48 found that, in addition to N = 1, highly rotationally excited He 2 (11 N 29) survives on droplet surfaces for at least 4 ms. 21–24 They concluded that the precursor is He(2 3 P), whose gas phase energy is 20.96 eV . They showed that He was not responsible for the He 2 they observed, and noted that He(2 3 P) correlates without a barrier to theb 3 g state (which radiates to a 3 + u ) as opposed to the barrier on a 3 + u (Figure 1). Their spectra did not include a region where transitions originate from He 2 high vibrational levels. Experiments by von Issendorff et al. 15 and an electronic structure study by Knowles and Murrell 25, 26 concluded that metastable He + 4 ( 4 A 2 ) is responsible for the anomalous He + 4 intensity. Its energy relative to 4 He atoms is 38.93 eV , which is consistent with the energy of two He atoms (39.64 eV) and the 41 eV threshold. Fine et al. found that the time required for metastables on surfaces of large droplets to undergo associative ionization to He + 4 was roughly 10s under their experimental conditions. 27 They argued that excited species (He , He 2 ) travel to droplet surfaces with minimal encounter with one another along the way. One thing is certain: He created at t = 0 exists in an environment rife with complex- ity, meaning questions abound. For example, does He create a bubble, or, alternatively, does it convert to He 2 which then creates a bubble? 28–30 Conversion of He to He 2 is central. In fact, it is relevant to all studies that involve electron impact creation of metastables in superfluid helium. 3.2 Results and Discussion We began the investigation with the potential energy curves of the He + He system, shown in Figure 3.1. The He 2 potential has a bond strength of 14,833 cm 131, 32 and an entrance barrier of 516 cm 1 relative to He + He. The barrier can be understood qualitatively as repulsion between Rydberg electron density and a He atom, countered 49 by attraction between the partially shielded He + core and the He atom. This barrier peaks at a distance of 2.72 ˚ A in agreement with the previous studies. 28, 33 In this region, the potential is quite flat decreasing by only 0.5 cm 1 on both sides of the peak at 2.70 and 2.74 ˚ A. This 516 cm 1 barrier seems insurmountable at T = 0.4 K in helium droplets. 34 That is the essence of the puzzle and we shall explain below how He 2 can be formed efficiently. Figure 3.1: The a 3 + u and c 3 + g potential energy curves correlate to He* + He. The vertical scale for the a 3 + u barrier region is expanded in the insert. To interrogate this issue, we employed a bottom up approach. We started with a simple model in which 3 He atoms are constrained to a straight line, Figure 3.2a. The energy landscape of the lowest triplet adiabat of this system is shown in figure 3.4. We then performed ab initio molecular dynamics simulations on the lowest triplet adiabatic 50 surface using potentials computed “on the fly”. The trajectories were launched from various combinations of initial interatomic distances,R 0 1 andR 0 2 , selected based on the superfluid helium radial distribution function (vide infra, Figure 3.5). 35–37 The initial kinetic energy was set to zero consistent with superfluid conditions. The simulations were run for 500 time-steps (484 fs). Following this simple 1D system, we studied the dynamics in 3 atom non-linear geometries. Various configurations were considered, specifically, 0 = 0 , 30 , 60 , 90 and a T-shaped starting geometry (refer to Figure 3.2c). The model was then extended to 4 atoms constrained to a line and finally a few 10 atom systems confined to a plane were studied. Figure 3.2: (a) 3 atoms in a linear configuration. (b) 4 atoms in a linear configura- tion. (c) 3-atom bent configurations included 0 = 30 , 60 , and 90 , and a T-shaped starting geometry. Initial interatomic distances satisfy R 0 2 R 0 1 and R 0 2 R 0 3 . Natural Transition Orbitals (NTOs) provide an easier understanding of electronic transitions often representing them with just a single dominant pair of electron and hole orbitals. For a deeper understanding of the lowest triplet excitation of the linear 3-atom Helium system, we looked at the NTOs at select configurations. Figure 3.3 shows slices of the lowest adiabat and the dominant pair of NTOs at two configurations (R 1 = 2.5 ˚ A andR 1 = 5 ˚ A) in each potential curve. 51 The potential curves show that there is initial repulsion between atoms 1 and 2 for R 0 1 values in the range 3 - 4 ˚ A. This repulsion pushes atom 1 away, and the force on atom 2 results in translational motion of theR 2 diatom and increased vibrational kinetic energy. The latter also promotes motion toward the He 2 potential minimum as observed in our simulations. This repulsive force from atom 1 can also provide the He 2 moiety with enough energy to result in its dissociation upon its return to the barrier peak region. The initial push toward the He 2 minimum and the return to the outer turning point region is a robust and important result of our linear 3-atom simulations. The NTOs atR 1 = 5 ˚ A indicate that the Rydberg electron prefers atom 1 over theR 2 diatom, whose distances in the figure are 2.7 - 3.0 ˚ A. This is not surprising, as theseR 2 values lie high on the He 2 entrance barrier (Figure 3.1) and the system prefers a lower energy state. The situation is more interesting forR 1 = 2.5 ˚ A. The Rydberg electron is located on atom 3 forR 2 = 2.8, 2.9, and 3.0 ˚ A for the same reason as above, it doesn’t want to pay the energy price of the entrance barrier at 2.5 ˚ A (the Rydberg electron on the R 1 diatom). In all these three cases, the potential curve has higher energy at R = 2.5 ˚ A than the other distances of R = 2.8, 2.9, and 3.0 ˚ A so the Rydberg electron prefers to be on the farthest atom, atom 3. The case of R 1 = 2.5 ˚ A and R 2 = 2.7 ˚ A is a compromise. The transition is delocalized with a node on the central atom. The Rydberg electron prefers the outer atoms as that electronic arrangement leads to least repulsion and lowest energy. Figure 3.4a shows the most relevant part of the lowest triplet surface of the linear 3-atom Helium system to the creation He 2 , i.e., the entrance barrier region viewed from large R 1 with R 2 < 2 ˚ A. Panel (b) shows a top view of the peak region. The blue sheet dropping sharply toward smallR 1 andR 2 is headed towards the region of strong binding i.e. the minimum of the potential surface. The grid (spacing of 0.1 ˚ A) puts 52 Figure 3.3: 1D, 3-atom lowest adiabat slices: NTOs are shown for 2.5 and 5.0 ˚ A (note dashed lines in top panel). The NTOs show holes and electrons on the left and right as indicated in the top panel; atom labels (1,2,3) are as per Figure 3.2a. As with MOs, red and blue orbital colors only have meaning within an NTO (phases of1). Orbital cutoff is chosen for viewing convenience. 53 Figure 3.4: (a) is a view of energy versus R 1 and R 2 for 3 atoms. (b) shows a top view. The rectangle contains the (R 0 1 ; R 0 2 ) combinations used in the trajectories. Gradients enable one to see how trajectories are launched from different (R 0 1 ; R 0 2 ). The arrows originate at the points where the gradients are evaluated. initial arrangements (R 0 1 , R 0 2 ) in perspective insofar as where they are located on the adiabat. TheR 1 andR 2 distances of the linear cases we examined lie within the black rectangle. The gradients of the surface explain the short time dynamics. For example, the forces within the rectangle are directed mainly alongR 1 . This results in lengthening of R 1 which causes R 2 to shorten. The contour energies shown in panel (b) are in units of cm 1 and show the isoenergetic curves on the surface. Figure 3.5 is the radial distribution function g(r) for liquid He, retrieved from reference 35. 35–37 It dictated the R 1 andR 2 values used throughout. 3.2.1 3-atom linear configurations As mentioned before, we propagated trajectories on the lowest triplet adiabat for many combinations ofR 0 1 andR 0 2 . The simulations provide insight into the exciton localization leading to He 2 formation. Figure 3.6a shows results for several representative cases. The R 2 (t) oscillation corresponds to the vibrational motion of the bound He 2 diatom while 54 Figure 3.5: The smooth curve is the calculated radial distribution function g(r) for helium at 1.21 K. Points show experimental data. Retrieved from reference 35. the third He atom moves away (R 1 (t)). The oscillation period is 100 fs which is quite long relative to the 20 fs period for energies near the bottom of thea 3 + u well. This reflects the highly anharmonic oscillation due to the He 2 being created high in the a 3 + u potential. For the linear 3-atom model, values ofR 0 2 smaller than 2.6 ˚ A always yield He 2 . Results for more initial configurations are given in Appendix A. Consider the trajectory with (R 0 1 ,R 0 2 ) = (3.2 ˚ A, 2.7 ˚ A) in Figure 3.6a. The system passes through the deep well before returning to the barrier peak region. The He 2 moiety gained internal energy due to the initial repulsion by the third atom early in the trajectory. This internal energy is enough to overcome the barrier and thus, dissociation takes place whenR 2 (t) returns to the region of the barrier peak. The same is true for configurations withR 0 1 = 3.2 - 3.5 ˚ A andR 0 2 = 2.7 ˚ A, as well as those withR 0 1 = 3.0 - 3.5 ˚ A andR 0 2 = 2.8 ˚ A. The traces for these dissociative cases are similar to the one for (R 0 1 , R 0 2 ) = (3.2 ˚ A, 2.7 ˚ A) and are included in Appendix A. The above examples capture the gist of He 2 creation for the chosen 3-atom (R 0 1 ,R 0 2 ) configurations. TheR 2 pair always entered the region of strong binding, as evidenced 55 Figure 3.6: Trajectories for (a) linear 3-atom cases (Figure 3.2a), and (b) linear 4- atom cases (Figure 3.2b). Note that R 2 (t) oscillation transfers to R 3 (t) in the first two panels. (c) Snapshots of a trajectory for 3-atoms with R 0 1 = 3.1 ˚ A, R 0 2 = 2.5 ˚ A, and 0 = 90 (Figure 3.2c). Smaller initial distances (R 0 2 ) are characteristic of the rising part of g(r). Larger initial distances (R 0 1 ) are characteristic of the peak region of g(r). by its rapid inward and outward passages through the region of the potential minimum. There were no exceptions. However, theR 2 pair did not always remain intact. In some initial configurations, theR 2 pair gained enough internal energy sufficient to overcome the barrier and dissociate upon it’s return to the barrier peak. Such results invite one to contemplate that a fourth atom might stabilize He 2 by removing some of its internal 56 energy. After all, with 3 atoms there was limited opportunity to remove enough energy from He 2 to prevent dissociation upon its return to the barrier region. 3.2.2 4-atom linear configurations Above observations led us to examine the linear 4-atom configurations next. Both the symmetric (R 0 1 = R 0 3 ) and non-symmetric (R 0 1 6= R 0 3 ) arrangements of atoms (Figure 3.2b) were studied. For the symmetric configurations, the middle two atoms always pass through the region of the potential minimum. The results for all the symmetric 4-atom trajectories are shown in Appendix A. For R 0 2 = 2.7 ˚ A, all R 0 1 values yielded He 2 ; for R 0 2 = 2.8 ˚ A, R 0 1 = 3.1-3.3 ˚ A yielded He 2 ; and for R 0 2 = 2.9 ˚ A, R 0 1 = 3.1 and 3.2 ˚ A yielded He 2 . It is clear that He 2 formation is more efficient when a 4th He atom is available to disperse some of the internal energy of the He 2 diatom. Non-symmetric arrangements are more interesting (Figure 3.6b). Initial distances between He atoms are crucial to whether the He 2 is stabilized or dissociates. Exciton transfer took place in 10 of the 24 initial configurations considered. All linear 3-atom and 4-atom trajectories are given in Appendix A. The above results show the subtle nature of these systems. Creation of the lowest triplet adiabat at t = 0 endows it with initial interatomic distances characteristic of su- perfluid helium. This results in weak attractive and repulsive interactions that launch, and participate in, the ensuing dynamics. How these dynamical processes play out is not obvious a priori. A detailed analysis is necessary. A complete theoretical descrip- tion includes non-adiabatic transitions among multiple adiabats, e.g., using the surface- hopping approach. We have used the adiabatic dynamics approach. In our simulations, the nulcei follow the lowest triplet adiabat and the electrons instantaneously adjust to 57 the nuclear positions. This overestimates the rate of exciton hopping between sites, es- pecially in the following 10-atom cases. Non-adiabatic effects will be included in our future work. 3.2.3 3-atom nonlinear configurations Next, we studied dynamics in nonlinear 3-atom configurations. Trajectories for nonlin- ear 3-atom cases used 0 = 30 , 60 , 90 (Figure 3.2c), and T-shaped initial configura- tions. The results for 30 were close to those for the corresponding 1D 3-atom cases. The orthogonal component of force from the atom at 30 is still small in most cases to make a significant difference in the dynamics compared to the linear counterparts. 0 = 60 results are also close to those for 0 = 30 . ForR 0 2 = 2.7 ˚ A, the tendency for stabilization of He 2 is higher than the corresponding 0 = 30 and the linear con- figurations. For example, for R 0 1 = 3.2 ˚ Aand R 0 2 = 2.7 ˚ A, the He 2 is stabilized in this case whereas it dissociates in the linear and 0 = 30 arrangements. The interplay of weak attractive and repulsive interactions is apparent in these results. A slight change in distance between the atoms (0.1 ˚ A) can lead to either stabilization or dissociation of the He 2 diatom. For 90 , even when the diatom initially localizes the exciton, it is discouraged from bonding because of the angular momentum it receives by being pushed from outside its center-of-mass. ForR 0 2 = 2.5 ˚ A, 0 = 90 , andR 0 1 values of 3.0 - 3.5 ˚ A, He 2 survived in all cases exceptR 0 1 = 3.2 ˚ A. Figure 3.6c shows snapshots of the trajectory forR 0 2 = 2.5 ˚ A, 0 = 90 , andR 0 1 = 3.1 ˚ A. ForR 0 2 = 2.6 ˚ A and 0 = 90 , onlyR 0 1 = 3.5 ˚ A yielded He 2 , and forR 0 2 = 2.7 - 2.9 ˚ A and 0 = 90 , He 2 did not survive in any of the trajectories. For T-shaped starting geometries, exciton localization at the diatom was retained for R 0 2 = 2.4 - 2.6 ˚ A and allR 0 1 values considered. The third atom and He 2 simply moved 58 apart. However, forR 0 2 = 2.7 - 2.9 ˚ A, repulsion from the third atom pushes theR 2 pair of He atoms away and it never moves toward the potential minimum. See Appendix A for trajectories of all non-linear arrangements. 3.2.4 10-atom 2D configurations Following the 3-atom and 4-atom models, we extended the model to 2D 10-atom con- figurations. Nine different initial configurations were considered and their trajectories on the lowest energy adiabat were calculated. Figure 3.7: Snapshots of a trajectory for 10-atoms showing exciton hopping. In each of the nine trajectories, highly vibrationally excited He 2 was formed within 100 fs. This system features a dense manifold of excited states with multiple crossings 59 and seams. It is highly likely for the trajectories to encounter multiple crossings during the 500 fs dynamics. We performed dynamics on the first triplet surface using the Born Oppenheimer approximation. On multiple occasions during the 500 fs dynamics, a pair of ground state He atoms came closer than the He 2 at its outer turning point, Figure 3.7. As a consequence of the Born Oppenheimer approximation, the exciton jumps to the new closer pair of He atoms: He 2 + He 0 -He 0 ! He-He + He 0 2 irrespective of the distance between the two sites. This is an unphysical result as the hopping probability decays exponentially with distance between sites.45 An accurate modeling of the formation of He 2 in the large 10 atom systems should include non-adiabatic transitions which will be included in our future work. It should be noted, however, that the He 2 yield is not affected by exciton hopping. Namely, He 2 is forthcoming regardless of whether the original He 2 transfers its exciton to a nearby He-He pair. The computational results are summarized as follows. In the linear and 0 = 30 3- atom cases, propagation on the lowest triplet adiabat usingR 0 1 andR 0 2 values in accord with g(r) results in a diatom being pushed toward the He 2 potential minimum. It does not always stay bound. Stabilization of the He 2 moiety formed with the initial push is achieved with an atom on the opposite side. Neighboring atoms provide means for dispersing some of the internal energy of the He 2 helping its stabilization. Sideways arrangements ( 0 = 90 and T-shaped) are repulsive. Preliminary results with 10-atom cases indicate facile He 2 production. Non-adiabatic effects are apparent in the 10-atom larger systems and will be included in our future work. 3.3 Computational Details Thea 3 + u andc 3 + g potential curves and 3, 4, and 10-atom adiabats were computed us- ing the equation-of-motion for excitation energies coupled-cluster approach with single 60 and double excitations (EOM-EE-CCSD). 38, 39 We used the doubly-augmented Dun- ning double- basis set (d-aug-cc-pVDZ) for diffuse wave functions of Rydberg states. Vibrational probability densities for a 3 + u were computed using the ‘eig’ function in MATLAB. Natural transition orbitals (NTOs) 40–43 provide insight into correlated ex- cited states by providing a convenient and compact way to visualize electronic transi- tions. Electron and hole orbitals are obtained by singular value decomposition (SVD) of the one-particle transition density matrix: exc (r h ;r e ) = X K K h K (r h ) e K (r e ) (3.1) Orbitals h K (r h ) and e K (r e ) denote hole and electron states, respectively, and the index K denotes the NTO corresponding to the singular value K . A small number of K usually dominate. This facilitates interpretation of excited electronic states. Ab initio molecular dynamics (AIMD) 44–46 was performed using gradients on the lowest triplet adiabat computed “on the fly”. Energy was calculated at each step using EOM-EE- CCSD with the d-aug-cc-pVDZ basis set. Nuclei were propagated using the Velocity Verlet algorithm. Trajectories were launched with zero initial kinetic energy, consistent with superfluid conditions and were propagated for 500 steps (484 fs) using a time-step of 40 a.u. (0.96755 fs). We used an energy convergence threshold of 10 6 and a conver- gence tolerance for the SCF cycle of 10 8 during the AIMD simulations. Calculations were carried out using the Q-Chem electronic structure package, 47, 48 and the natural transition orbitals (NTO’s) were plotted using Jmol. 49 61 3.4 Conclusions We have performed ab initio molecular dynamics in various Helium (1D and 2D) sys- tems with the gradients on the first triplet adiabatic surface computed “on the fly”. The results reveal that the creation of He in large helium droplets is followed by rapid, efficient production of highly vibrationally excited He 2 , despite a large barrier in the corresponding gas phase association reaction (Figure 3.1). The formation of He 2 in he- lium droplets is assisted by the interaction with the neighboring He atoms which remove the internal energy of the He 2 pair upon it’s return to the barrier region and stabilize it. The He 2 forms a bubble that travels to the droplet surface where it participates in the creation of He + 4 . Non-adiabatic effects are important, especially, in the larger 10-atom systems, where exciton hopping is observed over large distances. Our future work will include larger systems, non-adiabatic transitions, and more sophisticated sampling of initial conditions. 62 3.5 Appendix A: Trajectories for 3-atom linear geome- tries, 3-atom non-linear geometries and 4 atoms constrained to a straight line. Figure 3.8: Trajectories for linear 3-atom cases for R 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.6 - 2.8 ˚ A. The R 2 pair does not always remain intact, but it always passes through the region of the potential minimum as noted by the downward thrusts of the red traces. 63 Figure 3.9: Trajectories for 2D 3-atom geometries with = 30 . The initial dis- tances lie in the range R 0 1 = 3.0 - 3.4 ˚ A and R 0 2 = 2.6 - 2.8 ˚ A. Results are similar to those for the corresponding 1D 3-atom cases. 64 Figure 3.10: Trajectories for 2D 3-atom geometries with = 60 . The initial dis- tances lie in the range R 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.6 - 2.8 ˚ A. Results show He 2 formation takes place in more initial configurations than the corresponding 1D 3- atom cases (Figure 3.8) as well as the 2D 30 geometries (Figure 3.9). 65 Figure 3.11: Trajectories for 2D 3-atom geometries with = 90 . The initial dis- tances lie in the rangeR 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.4 - 2.6 ˚ A. The diatom pair, even when it initially localizes the exciton, is discouraged from bonding because of re- pulsion and the angular momentum it receives by being pushed from outside its center-of-mass. 66 Figure 3.11: Trajectories for 2D 3-atom geometries with = 90 . The initial dis- tances lie in the rangeR 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.7 - 2.9 ˚ A. The diatom pair, even when it initially localizes the exciton, is discouraged from bonding because of re- pulsion and the angular momentum it receives by being pushed from outside its center-of-mass. 67 Figure 3.12: Trajectories for 2D 3-atom geometries with T-shaped initial configu- rations. The initial distances lie in the range R 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.4 - 2.6 ˚ A. For R 0 2 = 2.4 - 2.6 ˚ A, initial localization of the exciton at the diatom is retained and the third atom separates from He 2 . 68 Figure 3.12: Trajectories for 2D 3-atom geometries with T-shaped initial configu- rations. The initial distances lie in the range R 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.7 - 2.9 ˚ A. For R 0 2 = 2.7 - 2.9 ˚ A, He 2 did not survive in any of the trajectories. 69 Figure 3.13: Trajectories for symmetric (R 0 1 = R 0 3 ) arrangements of 4-atoms ar- ranged in a linear fashion. The initial distances lie in the range R 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.6 - 3.0 ˚ A. The middle two atoms pass through the region of the potential minimum. For R 0 2 = 2.6 ˚ A, 2.7 ˚ A, all R 0 1 values yielded He 2 ; for R 0 2 = 2.8 ˚ A, R 0 1 = 3.0 - 3.3 ˚ A yielded He 2 ; for R 0 2 = 2.9 ˚ A, R 0 1 = 3.0 - 3.2 ˚ A yielded He 2 ; and for R 0 2 = 3.0 ˚ A, R 0 1 = 3.0 ˚ A and 3.1 ˚ A yielded He 2 . 70 Figure 3.13: Trajectories for symmetric (R 0 1 = R 0 3 ) arrangements of 4-atoms ar- ranged in a linear fashion. The initial distances lie in the range R 0 1 = 3.0 - 3.5 ˚ A and R 0 2 = 2.6 - 3.0 ˚ A. The middle two atoms pass through the region of the potential minimum. For R 0 2 = 2.6 ˚ A, 2.7 ˚ A, all R 0 1 values yielded He 2 ; for R 0 2 = 2.8 ˚ A, R 0 1 = 3.0 - 3.3 ˚ A yielded He 2 ; for R 0 2 = 2.9 ˚ A, R 0 1 = 3.0 - 3.2 ˚ A yielded He 2 ; and for R 0 2 = 3.0 ˚ A, R 0 1 = 3.0 ˚ A and 3.1 ˚ A yielded He 2 . 71 Figure 3.14: Trajectories for non-symmetric arrangements of 4-atoms arranged in a linear fashion. R 0 1 is fixed at 3.4 ˚ A, the other initial distances lie in the range R 0 2 = 2.6 - 2.9 ˚ A and R 0 3 = 2.8 - 3.3 ˚ A. Exciton transfer took place in 10 of the 32 initial configurations. 72 Chapter 3 References [1] J. Peter Toennies and Andrey F. Vilesov, Angewandte Chemie International Edi- tion 43, 2622 (2004). [2] M. Y . Choi, G. E. Douberly, T. M. Falconer, W. K. Lewis, C. M. Lindsay, J. M. Merritt, P. L. Stiles, and R. E. Miller, International Reviews in Physical Chemistry 25, 15 (2006). [3] W. E. Callegari, C.; Ernst, HandbookofHigh-resolutionSpectroscopy. John Wiley & Sons, 2011. [4] Frank Stienkemeier and Kevin K Lehmann, Journal of Physics B: Atomic, Molec- ular and Optical Physics 39, R127 (2006). [5] Josef Tiggesb¨ aumker and Frank Stienkemeier, Phys. Chem. Chem. Phys. 9, 4748 (2007). [6] Charles Bernando Sean M. O. O’Connell Deepak Verma Rico Mayro P. Tanyag, Curtis F. Jones and Andrey F. Vilesov, ColdChemistry: MolecularScatteringand ReactivityNearAbsoluteZero. The Royal Society of Chemistry, 2018. [7] J. Gspann, Surface Science 106, 219 (1981). [8] H. Buchenau, E. L. Knuth, J. Northby, J. P. Toennies, and C. Winkler, The Journal of Chemical Physics 92, 6875 (1990). [9] H. Buchenau, J. P. Toennies, and J. A. Northby, The Journal of Chemical Physics 95, 8134 (1991). [10] Berton E. Callicoatt, Kirk F¨ orde, Lilian F. Jung, Thomas Ruchti, and Kenneth C. Janda, The Journal of Chemical Physics 109, 10195 (1998). [11] Luis F. Gomez, Evgeny Loginov, Russell Sliter, and Andrey F. Vilesov, The Jour- nal of Chemical Physics 135, 154201 (2011). [12] Michael Renzler, Matthias Daxner, Nikolaus Weinberger, Stephan Denifl, Paul Scheier, and Olof Echt, Phys. Chem. Chem. Phys. 16, 22466 (2014). [13] H. Sch¨ obel, P. Bartl, C. Leidlmair, S. Denifl, O. Echt, T. D. M¨ ark, and P. Scheier, The European Physical Journal D 63, 209 (2011). [14] Jeonghee Seong, Kenneth C. Janda, Nadine Halberstadt, and Fernand Spiegel- mann, The Journal of Chemical Physics 109, 10873 (1998). 73 [15] B. von Issendorff, H. Haberland, R. Fr¨ ochtenicht, and J.P. Toennies, Chemical Physics Letters 233, 23 (1995). [16] J. W. Keto, F. J. Soley, M. Stockton, and W. A. Fitzsimmons, Phys. Rev. A 10, 872 (1974). [17] J. W. Keto, F. J. Soley, M. Stockton, and W. A. Fitzsimmons, Phys. Rev. A 10, 887 (1974). [18] J. W. Keto, M. Stockton, and W. A. Fitzsimmons, Phys. Rev. Lett. 28, 792 (1972). [19] J. C. Hill, O. Heybey, and G. K. Walters, Phys. Rev. Lett. 26, 1519 (1971). [20] R. L. Brooks, J. L. Hunt, and D. W. Tokaryk, The Journal of Chemical Physics 91, 7408 (1989). [21] C.-C. Hu, Photodetachmentofmetastableheliummoleculesfromheliumclusters, PhD thesis, University of Rhode Island, 1999. [22] T. Jianh, C. Kim, and J.A. Northby, Surface Review and Letters 03, 377 (1996). [23] Cheolkyu Kim, Photodetachment studies of negatively charged and metastably excitedheliumnanodroplets, PhD thesis, The University of Rhode Island, 1997. [24] S. Yurgenson, C.-C. Hu, C. Kim, and J.A. Northby, The European Physical Journal D - Atomic, Molecular, Optical and Plasma Physics 9, 153 (1999). [25] Peter J. Knowles and John N. Murrell, The Journal of Chemical Physics 102, 9442 (1995). [26] Peter J. Knowles and John N. Murrell, Molecular Physics 87, 827 (1996). [27] Jordan Fine, Deepak Verma, Curtis F. Jones, Curt Wittig, and Andrey F. Vilesov, The Journal of Chemical Physics 148, 044302 (2018). [28] S. L. Fiedler and J. Eloranta, Journal of Low Temperature Physics 174, 269 (2014). [29] J. Eloranta and V . A. Apkarian, The Journal of Chemical Physics 115, 752 (2001). [30] J. Eloranta, N. Schwentner, and V . A. Apkarian, The Journal of Chemical Physics 116, 4039 (2002). [31] C. Focsa, P.F. Bernath, and R. Colin, Journal of Molecular Spectroscopy 191, 209 (1998). [32] Michele Pavanello, Mauricio Cafiero, Sergiy Bubin, and Ludwik Adamowicz, In- ternational Journal of Quantum Chemistry 108, 2291 (2008). 74 [33] David R. Yarkony, The Journal of Chemical Physics 90, 7164 (1989). [34] M. Hartmann, R. E. Miller, J. P. Toennies, and A. Vilesov, Phys. Rev. Lett. 75, 1566 (1995). [35] D. M. Ceperley, Rev. Mod. Phys. 67, 279 (1995). [36] E. C. Svensson, V . F. Sears, A. D. B. Woods, and P. Martel, Phys. Rev. B 21, 3638 (1980). [37] H. N. Robkoff and R. B. Hallock, Phys. Rev. B 24, 159 (1981). [38] Rodney J. Bartlett, Wiley Interdisciplinary Reviews: Computational Molecular Science 2, 126 (2012). [39] Anna I. Krylov, Annual Review of Physical Chemistry 59, 433 (2008), PMID: 18173379. [40] Anatoliy Luzanov, A A. Sukhorukov, and V ´ E. Umanskii, 10, 354 (1976). [41] Stefanie A. Mewes, Felix Plasser, Anna Krylov, and Andreas Dreuw, Journal of Chemical Theory and Computation 14, 710 (2018), PMID: 29323887. [42] Felix Plasser, Michael Wormit, and Andreas Dreuw, The Journal of Chemical Physics 141, 024106 (2014). [43] Richard L. Martin, The Journal of Chemical Physics 118, 4775 (2003). [44] M. C. Payne, M. P. Teter, D. C. Allan, T. A. Arias, and J. D. Joannopoulos, Rev. Mod. Phys. 64, 1045 (1992). [45] R. Car and M. Parrinello, Phys. Rev. Lett. 55, 2471 (1985). [46] Mark E Tuckerman, Journal of Physics: Condensed Matter 14, R1297 (2002). [47] Yihan Shao et al., Molecular Physics 113, 184 (2015). [48] Anna I. Krylov and Peter M.W. Gill, Wiley Interdisciplinary Reviews: Computa- tional Molecular Science 3, 317 (2013). [49] Jmol: an open-source java viewer for chemical structures in 3d. http://www.jmol.org/. 75 Chapter 4: The Trivial Crossing Problem of FSSH 4.1 Introduction Tully’s fewest switches surface hopping (FSSH) 1 is one of the most popular approaches for simulating quantum-classical dynamics in a large variety of systems. 2–5 In this ap- proach, the nuclei are treated classically, while electrons are treated quantum mechani- cally, and transitions (hops) among coupled excited states incorporate feedback between electronic and nuclear subsystems. Its application to a broad range of systems has also brought forth various limitations of FSSH associated with nuclear 6, 7 and quantum deco- herence, 8–12 trivial crossings, 13 and superexchange, 14 among others. 5 Various modifica- tions to the standard FSSH algorithm have been suggested to overcome these limitations and significant progress has been made to enhance the efficiency and accuracy of FSSH when applied to the most challenging chemical and physical phenomena such as inter- nal conversion in multichromophore systems that can experience coherence and phase interference and where energy levels are dense and cross. 13–26 The non-crossing rule states that diatomic molecules cannot have level crossings and their adiabatic states always “avoid” each other. 27 Polyatomic molecules, on the other hand, always have “unavoided” crossings, generally termed as conical intersections. 28–30 76 These conical intersections are avoided along 2 dimensions and cross along the re- maining N-2 dimensions (where N is the total number of nuclear degrees of freedom). However, in extended polyatomic molecules, e.g., dimers, trimers, etc. only neighbor- ing molecules are significantly coupled and the coupling between the non-neighboring molecules is negligible. In such a case, the energy gap between two non-interacting po- tential energy surfaces (PESs) dominated by non-neighboring molecules is very small leading to trivial unavoided crossings. Trivial unavoided crossings are effectively (N-1)- dimensional intersections and the direct passage of a trajectory through the intersection seam is unavoidable. Trivial crossings present a huge numerical challenge in such ex- tended systems as the non-adiabatic coupling is zero except at the exact crossing point where it becomes infinite due to the vanishing energy gap. Dynamics propagated with a finite time-step would generally miss this crossing point and give unphysically long pop- ulation transfer time-scales and unrealistic long-range energy transfer. 31 The treatment of trivial crossings introduces another difficulty i.e. their identification. 13 Several strategies have been proposed to identify trivial crossings between non- interacting potential energy surfaces. 13, 32–34 In the adiabatic representation, tracking the state identities over time allows distinction between real and trivial crossings. This is achieved by calculating the overlap between nearest energy states and enforcing a unit hopping probability if the overlap is greater than a chosen threshold. 13 Other approaches take advantage of the diabatic representation to eliminate the spike in the non-adiabatic coupling. In self-consistent FSSH (SC-FSSH) 15 as well as the direct trajectories with surface hopping (DTSH) 34 method, the wavefunction is propagated in the diabatic basis and nonadiabatic couplings are not involved in the electronic dynamics. An alternative approach is to avoid the analytical calculation of the nonadiabatic coupling matrix el- ements at a finite number of time-steps altogether. The norm-preserving interpolation 77 (NPI) scheme is one such approach in which the time-derivative coupling is calculated from a set of relatively simple analytical expressions derived from the continuous in- terpolation of adiabatic electronic wavefunctions. 35, 36 NPI has been demonstrated to predict very accurate state transition probabilities when applied to a simple model of a weakly avoided crossing. Here, we analyze the electron-vibrational energy relaxation following photoexcita- tion of an extended spirobifluorene molecule using FSSH-based simulations. 37 Spiro linked compounds have mutually perpendicular conjugated-systems connected via a non-conjugated sp 3 hybridized atom. The tetrahedral character of the connecting atom maintains perpendicular arrangement of the rigid planar structures of the spiro linked units and leads to efficiently suppressed interactions between the-systems. Thus, the potential energy surfaces dominated by different polyfluorene units will exhibit triv- ial unavoided crossings. However, vibrational motion at 300 K can induce coupling between the two polyfluorene units of the spirobilfuorene resulting in both trivial cross- ings as well as real crossings between excited states. Therefore, we separated the two polyfluorene segments by 100 nm to ensure only trivial crossings between states lo- calized on different segments. Simulations of the photoinduced dynamics of both non- interacting separated polyfluorenes as well as weakly interacting spirobifluorene are per- formed within the Non-adiabatic Excited State Molecular Dynamics (NA-ESMD) 13, 38 framework using the standard FSSH, 1 and the additional Min-Cost correction. The chapter is organized as follows. The next section provides a brief overview of the numerical NA-ESMD approach and the Min-Cost state tracking algorithm. The electronic transition densities used to follow changes in localization and computational details are also described in this section. In Section 4.3 we present a comparison of the 78 simulated photoinduced dynamics of spirobifluorene and separated polyfluorene units with the two methods. Finally, Section 4.4 summarizes our findings and conclusions. 4.2 Theoretical Methodology 4.2.1 Non-adiabatic Excited State Molecular Dynamics The NA-ESMD framework provides an efficient and sufficiently accurate method for simulating photoinduced dynamics in large molecular systems on picosecond timescales. It has been specifically developed to simulate electronic and vibrational relaxation processes following photoexcitation in extended conjugated molecules in- volving multiple coupled electronic excited states. The general NA-ESMD framework is built utilizing standard FSSH, a mixed quantum classical approach in which an ensem- ble of nuclear trajectories evolve independently on a single PES at any given time and hops between excited states can occur based on the non-adiabatic coupling vector and a stochastic switching procedure. Excited state energies, gradients, and nonadiabatic coupling terms are calculated “on the fly” using the collective electron oscillator (CEO) approach. 39, 40 Configuration Interaction singles (CIS) 41 formalism is used to describe the correlated excited states implemented with the Austin model 1 (AM1) 42 semiempir- ical Hamiltonian. Additional details of the NA-ESMD approach, implementation, and limitations can be found elsewhere. 4, 13, 38, 43 NA-ESMD further employs a Min-Cost adiabatic state reassignment algorithm to identify unavoided crossings by tracking state identities over time. 13 New states at the current time step (i) are assigned in terms of old states at the preceding time step (i-1). The correspondence between new and old states is found based on the maximization of the trace of the square of the overlap matrix with elements defined as 79 s (t;t + t) (r; R(t)) (r; R(t + t)) = X n;m g (t) nm g (t + t) mn (4.1) If a maximum overlap greater than a defined threshold is identified, it is interpreted as a trivial crossing and the adiabatic states are reassigned by interchanging their popu- lations. The NA coupling is canceled, and the unit hopping probability is assigned. 4.2.2 Transition Density Analysis Electronic transition density (TD) matrices provide a useful means to track changes in the electronic wavefunction localization during dynamics independent of the adiabatic state populations. The diagonal elements of the TD matrices ( g ) nm =h (r;R(t))jc y m c n j g (r;R(t))i (4.2) computed within the CEO approach, 39, 40 represent the net change in the electronic den- sity distribution on an atomic orbital caused by an excitation from the ground state to an excited state. Here g (r;R(t)) and (r;R(t)) represent the CIS adiabatic ground and excited state wavefunctions, respectively. ‘n’ and ‘m’ are indices referring to the atomic orbital (AO) basis functions and c y m and c n are creation and annihilation opera- tors. The fraction of TD localized on a molecular segment is obtained by summing the contributions from each atom in the segment (index A) as ( g ) 2 Xunit = X n A ;m A ( g ) 2 n A m A (4.3) The usual normalization condition P n;m ( g ) 2 nm = 1 holds for the CIS approximation. 80 4.2.3 Simulation Details We have modelled the photoinduced dynamics of the spirobifluorene molecule as well as two polyfluorene segments separated by 100 nm whose optimized geometry are shown in Figure 4.1a. Simulations have been performed at 300 K using Min-Cost method and the standard FSSH. First, ground state MD was performed on spirobifluorene using a Langevin thermostat with a friction coefficient of 20 ps 1 and a classical time-step of 0.5 fs. 500 snapshots of spirobifluorene nuclei positions and momenta were collected from the thermally equilibrated ground-state MD trajectory to provide initial conditions for the excited state dynamics. 500 separated geometries were generated from the spirobi- fluorene geometries by moving one of the segments 100 nm along the vector connecting the origin to the central tetrahedral C atom. As the central C atom is shared between the polyfluorene units, one C atom and a total of 4 H atoms are added to satisfy valence. For each system (spirobifluorene and separated segments) and method (Min-Cost, and standard FSSH), 500 independent NA-ESMD trajectories were started from the initial configurations by instantly promoting the system to an excited state that has a fre- quency selected according to a Frank-Condon window defined as g (r;R) = f 2 exp[T 2 (E laser ) 2 ] (4.4) Heref represents the normalized oscillator strength for state, andE laser repre- sents the energy of a laser pulse centered at 350 nm. A Gaussian pulse,f(t) =exp( t 2 2T 2 ), with T 2 = 42.5 fs corresponding to a FWHM (full width at half maximum) of 100 fs has been used. All 500 trajectories were then propagated using Langevin dynamics for 1 ps with an electronic time-step of 0.025 fs and a nuclear time-step of 0.1 fs for the spirob- ifluorene system and 0.02 fs for the separated system. Ten electronic states and their 81 corresponding non-adiabatic coupling vectors were included in the simulations. For the Min-Cost method, the quantum time step was further reduced by a factor of 40 in order to identify and treat the trivial crossings. 4.3 Results and discussion 4.3.1 Absorption Spectra The chemical structure of the spirobifluorene, shown in Figure 4.1a(top), reveals the tetrahedral character of the central C atom connecting the two polyfluorene units main- taining the perpendicular orientation. Figure 4.1b shows the initial localization of the electronic transition densities for the two lowest energy excited states, S 1 and S 2 , of the spirobifluorene system. The states are completely localized on different segments con- firming that the spiro-linkage inhibits electronic coupling between the two segments. However, vibrational fluctuations at 300 K can induce coupling. On the other hand, separating the polyfluorene segments by 100 nm prohibits coupling between them. The structure in Figure 4.1a (bottom) shows the corresponding orientation of isolated polyfluorene units separated by 100 nm. The separation ensures that the excited states S 1 and S 2 are non-interacting and only experience trivial unavoided crossings between them. Thus, analysis of this system allows for a clear and unambiguous comparison of different methods for treatment of trivial crossings. Figure 4.1c shows the simulated absorption spectra at 300 K for both systems. The solid lines correspond to the spirobifluorene and the dashed lines correspond to the sep- arated system. The spectra have a high intensity peak at 420 nm and a lower intensity peak around 360 nm. The contribution of the individual excited states is also shown. The lower energy peak has contributions only from S 1 and S 2 while the peak at 360 nm 82 Figure 4.1: (a) Chemical structure of spirobifluorene (top) and the correspond- ing structure of polyfluorene segments separated by 100 nm (bottom). (b) Natu- ral transition orbitals for the S 0 ! S 1 and S 0 ! S 2 transitions at the optimized ground-state geometry for spirobifluorene. (c) Simulated absorption spectra with contributions of different excited states for spirobifluorene (solid lines) and sepa- rated polyfluorenes (dashed lines) at 300 K. has contributions from all higher energy states considered, i.e. S 3 -S 10 . The spectra for the two systems almost overlap with each other. The slight difference in the S 1 and S 2 absorption between the two systems is a consequence of the weak interaction between the polyfluorene units in spirobifluorene at 300 K induced by vibrational motion. As a consequence of the weakly coupled spirobifluorene units, the system is expected to ex- perience both real (interacting) and trivial (non-interacting) crossings during dynamics. 83 4.3.2 Photoinduced dynamics of non-interacting polyfluorene seg- ments In order to study the effect of trivial crossings on the relaxation rate using standard FSSH as well as the Min-Cost state tracking method, we studied the relaxation dynamics in the separated system. The large distance of 100 nm between the two polyfluorene units pro- hibits coupling between the segments even at 300 K and ensures only trivial unavoided crossings between excited states localized on different segments. The initial excitation generated at 350 nm corresponds to the range of states between S 3 and S 10 . Figure 4.1a shows the evolution of state populations computed as the fraction of trajectories in S 1 , S 2 , S 3 , S 4 and S 5 during dynamics. In standard FSSH, these are adiabatic states, while in the Min-cost method, the populations correspond to diabatic states named according to their initial energy ordering. For clarity, we introduce the notation S a n and S d n to distin- guish between adiabatic (a) and diabatic (d) states wheren refers to their original energy ordering. The electronic population decays to S d 3 on an ultrafast time-scale which then relaxes to the lowest excited state S d 1 within the first 500 fs for Min-cost (see Table 4.1). The NA coupling between states S d 1 and S d 2 , S d 2 and S d 3 as well as S d 1 and S d 3 are shown in Figures 4.3 & 4.5. The NA coupling between states S d 1 and S d 3 is significantly larger thanS d 2 andS d 3 resulting in population transfer directly to state S d 1 from S d 3 . Relaxation from S d 3 to S d 2 is relatively rare, owing to the smaller NA coupling, resulting in S d 2 never accumulating any significant fraction of the electronic population. Meanwhile, the sharp NA coupling peaks between states S d 1 and S d 2 correspond to trivial crossings. In contrast, for standard FSSH, S a 2 retains significant electronic population through- out the entire 1 ps dynamics. Since standard FSSH cannot detect trivial crossings, the 84 Figure 4.2: Simulation results for polyfluorene segments separated by 100 nm. Columns 1 and 2 correspond to (i) standard FSSH, and (ii) Min-Cost methods, re- spectively. (a) Average population on each electronic surface as a function of time obtained from the fraction of trajectories in each state. (b) Time-dependence of the average fraction of TD localized on segment 1 [( gf ) 2 A > 0:5 at t = 1000 fs] and seg- ment 2. (c) Time dependence of the average participation number. (d) Histogram of the number of TD localization changes between the polyfluorene segments during dynamics. 85 entire dynamics proceeds through adiabatic states that change their electronic charac- ter. This means that when the S d 1 and S d 2 have crossed, the electronic character of S a 2 and S a 1 are interchanged resulting in a large NA coupling between S a 2 and S a 3 rather than S a 1 and S a 3 , until the states re-cross. This allows population transfer to S a 2 where it be- comes trapped because the only pathway to S a 1 is through a trivial crossing, which is not detected in standard FSSH. Energy transfer can be followed by tracking the evolution of the transition density localized in each polyfluorene segment. In order to analyze this, we differentiate the polyfluorene chains based on the final TD localization denoted f. In each trajectory, the segment with the larger fraction of TD (( gf ) 2 A > 0:5) at t = 1000 fs is labelled as segment 1 and the other is labelled segment 2. Figure 4.2b shows the time dependence of the TD localized in segment 1 (blue line) and 2 (red line) averaged over the ensemble of trajectories for both standard FSSH and the Min-cost method. In standard FSSH, the final TD is localized on segment 1 (by definition), but at the initial time, the TD is localized in segment 2 in about half of the trajectories. Meaning that an unphysical long range energy transfer occurs from segment 2 to segment 1. The artificial energy transfer is the direct result of neglecting trivial crossings between S d 1 and S d 2 localized on different segments. This failure results in evolution on an adiabatic state whose elec- tronic character is changing due to trivial crossings resulting in an unphysical transfer of electronic TD from one segment to another over a distance of 100 nm. In comparison, the implementation of the Min-Cost state tracking to identify trivial crossings shows practically no change in transition density localization. In that case, the initial excita- tion is primarily localized in segment 1 in all the trajectories and it remains localized in segment 1 for the entire dynamics. That clearly indicates that trivial crossings between 86 excited states are identified and the nuclear trajectory hops between the crossing adia- batic states. In other words, the electronic population stays in the same diabatic state and avoids sudden unphysical changes in the transition density localization. Figure 4.3: (a) Absolute value of the non-adiabatic coupling terms (NACT) be- tween states (S 1 , S 2 ), (S 2 , S 3 ) and (S 1 , S 3 ) averaged over all trajectories, (b) Relative probability of changing TD localization between segments as a function of time, for polyfluorene segments separated by 100 nm. Columns 1 and 2 correspond to (i) standard FSSH, and (ii) Min-Cost methods, respectively. To provide further insight into the changes in transition density localization, the participation number is analyzed. The extent of (de)localization of the excitation can be seen by the participation number, defined as, PN = " X segmentX ( g ) 2 segmentX 2 # 1 (4.5) where X = 1 and 2. By definition, the participation number can range from 1PN 2. A value of PN = 1 signifies complete TD localization within one polyfluorene unit, 87 whereas PN = 2 implies complete delocalization between the two segments. The evolu- tion of participation number during dynamics is plotted in Figure 4.2c. For all methods, the value of PN = 1 provides clear evidence that the transition density is indeed always localized on a single segment, as expected for non-interacting segments separated by 100 nm, and the system does not experience any delocalization. The changes in local- ization can be seen more clearly in Figure 4.2d. Here, the average number of changes in TD localization per trajectory between the polyfluorene segments is shown. In standard FSSH, the average number of changes in TD localization per trajectory is 56 whereas in Min-Cost, there are 0.08 changes in TD localization per trajectory. Figure 4.3b shows the relative probability that a trajectory will experience a change in transition density localization at any given time during dynamics. As can be seen, in standard FSSH, changes in TD take place throughout dynamics with a large probability whereas in Min- Cost, the probability of changes in TD localization is negligible. 4.3.3 Photoinduced dynamics of interacting spirobifluorene Next, the NA-ESMD simulations were repeated for the spirobifluorene molecule shown in Figure 4.1a (top) to analyze the performance of the trivial crossing algorithms in the weak coupling case when both real crossings between interacting states and trivial cross- ings between non-interacting states can exist. The spirobifluorene system excited at 350 nm results in the initial excitation of S 3 -S 10 , consistent with the initial excitation created in the non-interacting separated polyfluorene segments. Figure 4.4a shows the evolution of the population for S 1 , S 2 , S 3 , S 4 and S 5 following photoexcitation of spirobifluo- rene. Similar to the separated system, the electronic population from the high energy excited states first undergoes ultrafast relaxation to the intermediate S d 3 within 100 fs. Subsequently, the population of S d 3 decays on a relatively slow time-scale concomitant 88 Figure 4.4: Simulation results for spirobifluorene. Columns 1 and 2 correspond to (i) standard FSSH, and (ii) Min-Cost method, respectively. (a) Average population on each electronic surface as a function of time obtained from the fraction of tra- jectories in each state. (b) Time-dependence of the average fraction of TD localized on segment 1 [( gf ) 2 A > 0:5 at t = 1000 fs] and segment 2. (c) Time dependence of the average participation number. (d) Histogram of the number of TD localization changes between the polyfluorene segments during dynamics. 89 with the rise of population in S d 1 . Again, S d 2 never sustains any significant electronic population owing to the larger NA coupling between S d 3 and S d 1 compared to S d 3 and S d 2 shown in Figure 4.5a. However, in the case of spirobifluorene, S a 2 remains unpopu- lated in standard FSSH whereas S a 2 retained a significant fraction of population in the non-interacting system modeled by standard FSSH. In this case, the population from S a 2 can relax to S a 1 via quantum hops between interacting states. The absence of popula- tion accumulation in S 2 , even for standard FSSH, indicates electronic coupling between chromophores leading to efficient downhill non-adiabatic transitions. To follow energy transfer in the spirobifluorene, the evolution of the average TD localized in each polyfluorene chain is shown in Figure 4.4b. Again, the segment with the larger final fraction of TD (( gf ) 2 A > 0:5) at t = 1000 fs is labelled segment 1 and the other is labelled segment 2. For all simulations, the initial TD has equal probability of being localized on either segment. During dynamics, the TD localized on segment 2 decreases and is transferred to segment 1. Since this is an ensemble average, it is not clear whether the initial TD is delocalized between the two segments or whether the initial excitation was localized on segment 1 in half of the trajectories and localized on segment 2 in the remaining trajectories. This is clarified by the evolution of the partici- pation number shown in Figure 4.4c. It reveals that initially delocalized excitations are created which become completely localized on one polyfluorene segment within 500 fs, as seen by the decay of participation number from an initial value of PN = 1.25 to PN = 1.03 for both methods. Energy transfer from one segment to the other occurs in spirobifluorene owing to the electronic coupling induced relaxation via quantum hops between interacting localized states. Alternatively, adiabatic relaxation via trivial cross- ings between non-interacting states does not result in a change in TD localization as demonstrated in the separated polyfluorene system. 90 Figure 4.5: (a) Absolute value of the non-adiabatic coupling terms (NACT) be- tween states (S 1 , S 2 ), (S 2 , S 3 ) and (S 1 , S 3 ) averaged over all trajectories, (b) Relative probability of changing TD localization between segments as a function of time, for spirobifluorene. Columns 1 and 2 correspond to (i) standard FSSH, and (ii) Min-Cost methods, respectively. The average number of changes in the TD localization between segments is shown in Figure 4.4d. In standard FSSH, the localization changes an average of 5 times per trajectory. The changes in localization can arise either due to non-adiabatic transitions between localized interacting states or due to missed trivial crossings where the system remains on the same adiabatic state that changes its electronic character. On the other hand, in Min-Cost, changes in localization can only occur via non-adiabatic transitions involving a change in electronic character of the active state but not at trivial crossings because the electronic character is preserved. Here, the average number of changes in TD localization is reduced to 4. Changes in TD localization persist throughout dynamics for both methods without any significant deviation, as shown in Figure 4.5b. Thus, the 91 effect of trivial crossings is masked by the presence of coupling between the excited states in the weak coupling limit. 4.3.4 Relaxation time-scales The relaxation time-scale can be found by tracking the rise in population of S 1 ob- tained by fitting the S 1 curves in Figures 4.2a and 4.4a with the function 1exp t . The relaxation time-scales are presented in Table 4.1 for spirobifluorene and separated polyfluorene units. For both systems, Min-Cost, has faster relaxation rates than standard FSSH. This is a direct consequence of the trivial crossings encountered between states S 1 and S 2 during dynamics. Standard FSSH fails to identify trivial crossings, there- fore, the system remains on the upper adiabatic surface at trivial crossings and reaches the lowest excited state only via non-adiabatic transitions resulting in a longer relax- ation time-scale. Relaxation to S d 1 in the separated system can only occur through trivial crossings. In contrast, the spirobifluorene system has the additional relaxation path- way through interacting states with significant non-adiabatic coupling. Therefore, the spirobifluorene generally exhibits a faster relaxation timescale than the non-interacting separated polyfluorene system, as expected. For Min-Cost, relaxation to S d 1 can occur both via non-interacting trivial crossings as well as non-adiabatic transitions between interacting states. Table 4.1: Rise-time of the average population of state S 1 (S a 1 for FSSH, S d 1 for the Min-Cost method) for separated system and spirobifluorene. Standard FSSH Min-Cost 100 nm separated polyfluorenes 427.9 335.7 Spirobifluorene 271.1 230.7 92 4.4 Conclusions We have performed simulations of the photoinduced relaxation dynamics of a spirob- ifluorene molecule as well as for isolated polyfluorene segments separated by 100 nm using standard FSSH and Min-Cost state tracking method that solves the trivial crossings problem in FSSH. The performance of Min-Cost state tracking algorithm implemented within the NA-ESMD framework has been evaluated based on comparison to standard FSSH. In the separated polyfluorene system, the only pathway of energy relaxation from state S 2 to the lowest excited state S 1 is via trivial crossings. In standard FSSH, the elec- tronic population that relaxes to S a 2 does not find this pathway and remains in S a 2 over the course of the dynamics up to 1000 fs. However, in the corrected method, any population that relaxes to S d 2 is quickly transferred to S d 1 via a trivial crossing resulting in faster relaxation timescales for Min-Cost, compared to standard FSSH. In FSSH, the TD lo- calization changes from one segment to the other an average 56 times per trajectory as it fails to identify trivial crossings and misses to hop. The Min-Cost state tracking method identifies most of the trivial crossings and has negligible TD localization changes be- tween the polyfluorene segments as expected. However, the algorithm does not identify exactly all of the trivial crossings as it is sensitive to the chosen parameters for overlap threshold between adiabatic states to identify a trivial crossing as well as the number of times the quantum time-step is reduced in the trivial crossing region. 93 4.5 Appendix A: Natural Transition Orbitals for S 1 - S 10 excited states of spirobifluorene at 300 K 94 95 96 Chapter 4 References [1] John C. Tully, The Journal of Chemical Physics 93, 1061 (1990). [2] Linjun Wang, Alexey Akimov, and Oleg V . Prezhdo, The Journal of Physical Chemistry Letters 7, 2100 (2016), PMID: 27171314. [3] Svetlana Kilina, Dmitri Kilin, and Sergei Tretiak, Chemical Reviews 115, 5929 (2015), PMID: 25993511. [4] Tammie Nelson, Sebastian Fernandez-Alberti, Adrian E. Roitberg, and Sergei Tre- tiak, Accounts of Chemical Research 47, 1155 (2014), PMID: 24673100. [5] Joseph E. Subotnik, Amber Jain, Brian Landry, Andrew Petit, Wenjun Ouyang, and Nicole Bellonzi, Annual Review of Physical Chemistry 67, 387 (2016), PMID: 27215818. [6] Linjun Wang, Guangjun Nan, Xiaodi Yang, Qian Peng, Qikai Li, and Zhigang Shuai, Chem. Soc. Rev. 39, 423 (2010). [7] Linjun Wang, Alexey V . Akimov, Liping Chen, and Oleg V . Prezhdo, The Journal of Chemical Physics 139, 174109 (2013). [8] Eric R. Bittner and Peter J. Rossky, The Journal of Chemical Physics 103, 8130 (1995). [9] Ahren W. Jasper and Donald G. Truhlar, The Journal of Chemical Physics 127, 194306 (2007). [10] Oleg V . Prezhdo, The Journal of Chemical Physics 111, 8366 (1999). [11] Michael D. Hack and Donald G. Truhlar, The Journal of Chemical Physics 114, 9305 (2001). [12] Michael J. Bedard-Hearn, Ross E. Larsen, and Benjamin J. Schwartz, The Journal of Chemical Physics 123, 234106 (2005). [13] Sebastian Fernandez-Alberti, Adrian E. Roitberg, Tammie Nelson, and Sergei Tre- tiak, The Journal of Chemical Physics 137, 014512 (2012). [14] Linjun Wang, Dhara Trivedi, and Oleg V . Prezhdo, Journal of Chemical Theory and Computation 10, 3598 (2014), PMID: 26588504. 97 [15] Linjun Wang and Oleg V . Prezhdo, The Journal of Physical Chemistry Letters 5, 713 (2014), PMID: 26270842. [16] Linjun Wang, Andrew E. Sifain, and Oleg V . Prezhdo, The Journal of Physical Chemistry Letters 6, 3827 (2015), PMID: 26722878. [17] Linjun Wang, Andrew E. Sifain, and Oleg V . Prezhdo, The Journal of Chemical Physics 143, 191102 (2015). [18] Philip Shushkov, Richard Li, and John C. Tully, The Journal of Chemical Physics 137, 22A549 (2012). [19] Neil Shenvi, Joseph E. Subotnik, and Weitao Yang, The Journal of Chemical Physics 135, 024101 (2011). [20] V Gorshkov, Sergei Tretiak, and Dmitry Mozyrsky, 4, 2144 (2013). [21] Joseph E. Subotnik and Neil Shenvi, The Journal of Chemical Physics 134, 024105 (2011). [22] Alexey V . Akimov, Run Long, and Oleg V . Prezhdo, The Journal of Chemical Physics 140, 194107 (2014). [23] Linjun Wang and David Beljonne, The Journal of Physical Chemistry Letters 4, 1888 (2013), PMID: 26283125. [24] Alexey V . Akimov, Dhara Trivedi, Linjun Wang, and Oleg V . Prezhdo, Journal of the Physical Society of Japan 84, 094002 (2015). [25] Heather M. Jaeger, Sean Fischer, and Oleg V . Prezhdo, The Journal of Chemical Physics 137, 22A545 (2012). [26] Alexey V . Akimov and Oleg V . Prezhdo, Phys. Rev. Lett. 113, 153003 (2014). [27] J. von Neuman and E. Wigner, Physikalische Zeitschrift 30, 467 (1929). [28] David R. Yarkony, Rev. Mod. Phys. 68, 985 (1996). [29] Fernando Bernardi, Massimo Olivucci, and Michael A. Robb, Chem. Soc. Rev. 25, 321 (1996). [30] M. Baer, Beyond Born-Oppenheimer: Electronic Nonadiabatic Coupling Terms andConicalIntersections. Wiley, 2006. [31] Tammie Nelson, Sebastian Fernandez-Alberti, Adrian E. Roitberg, and Sergei Tre- tiak, Chemical Physics Letters 590, 208 (2013). 98 [32] E. Fabiano, T.W. Keal, and W. Thiel, Chemical Physics 349, 334 (2008), Electron Correlation and Molecular Dynamics for Excited States and Photochemistry. [33] Christian Evenhuis and Todd J. Mart´ ınez, The Journal of Chemical Physics 135, 224110 (2011). [34] G. Granucci, M. Persico, and A. Toniolo, The Journal of Chemical Physics 114, 10608 (2001). [35] Garrett A. Meek and Benjamin G. Levine, Chemical Physics 460, 117 (2015), Novel aspects of the Jahn-Teller effect in molecular and solid-state systems. [36] Yinan Shu and Benjamin G. Levine, The Journal of Physical Chemistry C 120, 23246 (2016). [37] Tobat P. I. Saragi, Till Spehr, Achim Siebert, Thomas Fuhrmann-Lieker, and Josef Salbeck, Chemical Reviews 107, 1011 (2007), PMID: 17381160. [38] Tammie Nelson, Sebastian Fernandez-Alberti, Vladimir Chernyak, Adrian E. Roit- berg, and Sergei Tretiak, The Journal of Chemical Physics 136, 054108 (2012). [39] Sergei Tretiak and Shaul Mukamel, Chemical Reviews 102, 3171 (2002), PMID: 12222985. [40] Shaul Mukamel, Sergei Tretiak, Thomas Wagersreiter, and Vladimir Chernyak, Science 277, 781 (1997). [41] D. J. Thouless, Thequantummechanicsofmany-bodysystems[by]D.J.Thouless. Academic Press New York, 2d ed. edition, 1972. [42] Michael J. S. Dewar, Eve G. Zoebisch, Eamonn F. Healy, and James J. P. Stewart, Journal of the American Chemical Society 107, 3902 (1985). [43] Tammie Nelson, Sebastian Fernandez-Alberti, Vladimir Chernyak, Adrian E. Roit- berg, and Sergei Tretiak, The Journal of Physical Chemistry B 115, 5402 (2011), PMID: 21218841. 99 Chapter 5: Ehrenfest Dynamics with Decoherence and Detailed Balance 5.1 Introduction Quantum nonadiabatic (NA) processes, such as electron and proton transport, 1 exci- ton relaxation, 2, 3 dissociation, 4, 5 and recombination, 6, 7 and energy transfer, 8–11 are widespread in nature. 12–16 A purely quantum description of large systems, typically involved in such processes, is too complex and time consuming and thus, approximate semiclassical treatments are unavoidable. Such phenomena can generally be modelled by representing dynamics of a few key quantum particles coupled to classical degrees of freedom. Mixed quantum-classical dynamics (MQCD) methods, in which electrons are treated at the quantum level and nuclei are treated classically, have become the most appropriate approach to model these processes. Fewest switches surface hopping (FSSH), 17 is one of the most popular NA molecular dynamics (MD) techniques. In this approach, the dynamics of classical particles evolv- ing on a potential energy surface is described by Newtonian equations, while quantum particles are propagated according to the time-dependent Schr odinger equation. NA transitions (hops) among coupled electronic states incorporate feedback between elec- tronic and nuclear subsystems. The FSSH probabilities are designed to minimize the 100 number of state switches. Once a surface hop is assigned, the nuclear velocities are adjusted in the direction of the NA coupling vector so as to conserve the total energy of the system. If the nuclear kinetic energy is less than the energy required for the quan- tum transition, the hop is rejected. Such velocity rescaling and hop rejection produce detailed balance between the upward and downward transitions, approximately repro- ducing the Boltzmann distribution of the excited state populations. 18, 19 An ensemble of surface hopping trajectories are simulated independently, and the population of each quantum state is given by the fraction of trajectories that are assigned to that state. In the absence of velocity rescaling, FSSH probabilities satisfy internal consistency, i.e. the fraction of trajectories on a given state at any time equals the quantum population of that state obtained by the time-dependent Schr odinger equation. Individual trajectories always propagate on a pure state (eigenstate) of the system. Owing to these features, the FSSH method has been successful in many applications, making it immensely popular for modelling NA processes. Its conceptual simplicity and ease of implementation have helped FSSH retain its popularity. The popularity of FSSH has also brought forward its limitations, as demonstrated by the numerous types of corrections and modifications to the standard FSSH algorithm. A detailed discussion of such methods can be found in this review. 20 Since FSSH trajec- tories hops between potential energy surfaces, it is not straightforward to apply FSSH to systems with a continuum of states, and a discretization of the continuum is required for FSSH dynamics. 21 In large scale systems, such as extended polyatomic molecules, molecular crystals or assemblies of quantum dots, there may exist electronic states lo- calized in different parts of a system and coupled very weakly. 22 Such non-interacting states may have similar energy and cross during dynamics, leading to trivial crossings in which the NA couplings behave as delta functions of time. Very small time-steps are 101 required for identification of such crossings, which is extremely time-consuming. 23 Fur- ther, FSSH is known to give more accurate results in the adiabatic representation rather than diabatic representation, because the velocity adjustment process differs strongly in the two representations. The representation dependence of the results may be viewed as an undesirable property. 24 All MQCD methods are intrinsically limited by the classi- cal mechanical description of the nuclear motion. Many nuclear quantum effects, such as tunneling and zero-point vibrational energy that are important for low temperature dynamics, cannot be easily described by classical trajectories. 25, 26 Last but not least, loss of quantum coherence in the electronic subsystem induced by coupling to nuclear motions has to be addressed. 27–34 Ehrenfest dynamics is another primary MQCD method. Here, classical particles evolve on a mean potential energy surface which is averaged over all quantum states and weighted by the corresponding populations. Ehrenfest dynamics is fully determin- istic and hence, more time-efficient than FSSH which is based on a stochastic algorithm and requires averaging over multiple realizations of the random process to exhibit good statistics. Additional sampling may be needed for rare events. 30 Although the Ehren- fest technique is simple and easy to implement, the mean-field approximation is only adequate when either the quantum and classical subsystems are weakly coupled to each other, or when the classical subsystem reacts similarly to all quantum states included in the active space. Another serious limitation of the Ehrenfest approach is its inability to collapse the system’s wavefunction on a pure state. After passing through a region of strong electronic coupling, the trajectory evolves with a wavefunction that is a combina- tion of nearby adiabatic eigenfunctions even in regions of negligible coupling between the states. 31 However, this deficiency of the Ehrenfest method can be resolved if the 102 results of the Ehrenfest dynamics are interpreted as the average values of the proper- ties and their distributions. Moreover, both the Ehrenfest MD and the FSSH techniques are fully quantum coherent, they do not account for decoherence. 32 Many decoherence correction schemes for FSSH as well as Ehrenfest dynamics have been proposed. 33, 34 Additionally, the standard Ehrenfest method does not satisfy detailed balance between transitions downward and upward in energy, and therefore, it cannot generate thermody- namic equilibrium between quantum and classical subsystems. In particular, the Ehren- fest method does not produce the Boltzmann equilibrium populations of the quantum states at a finite temperature. 3, 35 A proper theory that accounts for decoherence and de- tailed balance is of fundamental importance to accurately simulate processes happening in complex condensed matter systems. In this work, we propose a modified Ehrenfest method with corrections to account for both decoherence and detailed balance. Decoherence is introduced via the coher- ence penalty functional (CPF) method developed by Akimov et al. 34 The CPF method dynamically penalizes development of coherences during evolution of quantum degrees of freedom via an additional term in the classically mapped Hamiltonian, thus preserv- ing the overall Hamiltonian structure of the equations of motion. The detailed bal- ance correction is applied to the coupling matrix of the Hamiltonian, as proposed by Bastidaetal. 36 for simulating vibrational energy transfer between quantum systems and classical baths. We have implemented the two corrections to the Ehrenfest technique within the PYXAID software package 35, 37 and performed NA dynamics in a hybrid molecule/metal-tip system to analyze the performance. 103 5.2 Theory: Formulation of the Ehrenfest method with decoherence and detailed balance corrections We begin with a derivation of the Ehrenfest technique. Then we present the decoherence and detailed balance corrections to the Ehrenfest technique, and formulate equations of motions that include both corrections. After that, we describe the optical response theory formalism to compute the decoherence rates, and provide details of the simulation of charge injection, relaxation and recombination in a molecule/metal-tip system studied experimentally 38, 39 by scanning tunneling microscopy (STM). 5.2.1 The Ehrenfest method Consider the time-dependent (TD) Schr odinger equation (SE): i~ @ (~ r (3n) ; ~ R (3N) ;t) @t =H(~ r (3n) ; ~ R (3N) ;t) (~ r (3n) ; ~ R (3N) ;t) (5.1) where~ r (3n) and ~ R (3N) are coordinates of n electrons and N nuclei, respectively. For simplicity of notation, the electronic and nuclear degrees of freedom will be denoted by r andR, respectively. The non-relativistic molecular Hamiltonian is defined as H(r;R;t) = N X J=1 1 2M J r 2 J n X j=1 1 2 r 2 j + X J<K Z J Z K ~ R J ~ R K + X j<k 1 j~ r j ~ r k j X J;j Z J ~ R J ~ r j (5.2) 104 whereM J andZ J are the mass and the charge of theJ-th nucleus respectively. Com- bining all the terms except the nuclear kinetic energy term defines the electronic Hamil- tonian operator. H(r;R;t) = N X J=1 1 2M J r 2 J +H el (r;R;t) (5.3) Next, consider the Born-Oppenheimer separation ansatz for the wavefunction (r;R;t). The total wavefunction (r;R;t) is described by a linear combination of products of electronic wavefunctions for different adiabatic states, i (r;R(t)), and nuclear wave functions, i (t;R(t)), associated with the electronic states: (r;R;t) = X i i (t;R(t)) i (r;R(t)) (5.4) The adiabatic Born-Oppenheimer representation is based on the fact that the mass of the nuclei is much larger than the mass of the electrons, and the electrons move rapidly in the nuclear field. Using the representation of the total wavefunction given in Eq. 5.4, inserting it into the TD-SE, Eq. 5.1, and projecting it onto the corresponding elec- tronic states, one obtains the dynamics of nuclear wavepackets correlated with different electronic states: i~ @ i (t;R(t)) @t = X j ( T nucl + ~ E j (t;R) ij ~ 2 ~ d (1) ij M ~ r~ 2 ~ d (2) ij 2M ) i (t;R(t)) (5.5) where ~ d (n) ij are the nth order non-adiabatic couplings between electronic statesi andj, defined as ~ d (1) ij = h i (r;R(t))j ~ r j (r;R(t))i (5.6a) ~ d (2) ij = h i (r;R(t))jr 2 j (r;R(t))i (5.6b) 105 The quantities ~ E i (t;R)) are the adiabatic potential energy surfaces, corresponding to the electronic states i (r;R(t)). They are the eigenvalues of the stationary electronic SE, parametrically dependent on the nuclear coordinatesR(t). H el (r;t;R) i (r;t;R(t)) = ~ E i (t;R) i (r;t;R(t)) (5.7) The diagonal parts of the term~ 2 ~ d (2) ij 2M in Eq. 5.5 can be considered as a correction to the adiabatic energies: E i (R) = ~ E i (R)~ 2 ~ d (2) ij 2M (5.8) The off-diagonal terms representing the second order NA coupling between electronic states are often neglected under the assumptions that adiabatic wavefunctions are slowly varying with respect to nuclear coordinates. 40–42 However, these terms may become im- portant in some cases, and one must be careful when using this ansatz. 43 After applying these assumptions, Eq. 5.5 simplifies to: i~ @ i (t;R(t)) @t = X j ( (T nucl +E j (t;R)) ij ~ 2 ~ d (1) ij M ~ r ) j (t;R(t)) (5.9) The next approximation involves treating the nuclear degrees of freedom classically, using the following standard quantum-classical correspondence rules. T nucl ~ 2 2M r 2 ! P 2 2M (5.10) ~ 2 ~ d (1) ij M r!i~ ~ d (1) ij ~ P M (5.11) 106 The solution of the TD-SE is expressed on a finite set of adiabatic basis functions defined as the eigenfunctions of the stationary SE, Eq. 5.7,f i g,i = 0;:::;N b 1, whereN b is the basis size: (r;R;t) = N b 1 X i=0 c i (t) i (r;R(t)) (5.12) The notation i (r;R(t)) indicates functional dependence of the electronic wave func- tions, i , on the electronic coordinates,r, and its parametric dependence on the nuclear coordinates,R. It should be noted that with this definition, there exists no explicit func- tional dependence of (r;R;t) on the nuclear positionsR. Using the ansatz, Eq. 5.12 together with Eqs. 5.10 and 5.11, one converts the TD- SE, Eq. 5.9, into the equations of motion for the amplitudes,c i (t): i~ @c i @t = N b 1 X j=0 E i ij i~ ~ d (1) ij ~ P M ! c j (5.13) Solution of Eq. 5.13 describes evolution of the electronic degrees of freedom. The first term represents evolution of the coefficient of an eigenfunction of the Hamiltonian, and the second term couples the electronic states with each other via the nuclear kinetic energy and the non-adiabatic coupling. 107 In the Ehrenfest method, the nuclei move on the average electronic potential energy surface, subject to the mean-field force: 44 ~ F mf =h (r;R;t)jr R H el (r;t;R)j (r;R;t)i = * X i c i (t) i (r;R(t))jr R H el (r;t;R)j X j c j (t) j (r;R(t)) + = * X i c i c i i jr R H el j i + * X i;j j6=i c i c j i jr R H el j j + = X i jc i j 2 ~ F i + X i;j c i c j (E i E j )d (1) ij (5.14) where ~ F i = i j @H el @R j i (5.15) is the Hellmann-Feynman force corresponding to the adiabatic electronic state i . The mean-field force is computed on the basis of the contribution of each adiabatic state to the time-dependent electronic wave function. The first term couples the populations of the electronic statesjc i j 2 to the nuclear trajectory, while the second term with the coefficients c i c j includes interferences between these states. The system of coupled differential equations 5.13 and 5.14 constitute the Ehrenfest method. The trajectories generated in the Ehrenfest dynamics are fully deterministic and continuous. 5.2.2 Coherence penalty functional As noted above, the Ehrenfest dynamics does not account for decoherence that would be induced within the electronic subsystem if the nuclei were treated quantum mechan- ically. Following the work of Alexey Akimov, Run Long and Oleg Prezhdo, 34 we add the decoherence correction to the dynamics via the CPF method as follows. 108 In the variablesq i =Re(c i ) andp i =Im(c i ), the semiclassical TD-SE, i~ @c i @t = X j E i (R(t)) ij i~ P M d (1) ij c j (5.16) becomes equivalent to the following Hamiltonian equations of motion _ p i = @H @q i (5.17a) _ q i = @H @p i (5.17b) with the Hamiltonian given by H = X i E i 2~ (q 2 i +p 2 i ) P M X i;j d (1) ij p i q j (5.18) The following augmented Hamiltonian to account for decoherence effects has been pro- posed: ~ H =H + X i;j i6=j ij (q 2 i +p 2 i )(q 2 j +p 2 j ) (5.19) where the term (q 2 i +p 2 i )(q 2 j +p 2 j ) =jc i c j j 2 is simply a square of the magnitude of coherence between the pair of states i and j. This term penalizes coherence develop- ment between pairs of states. The constant ij determines the magnitude of the penalty for each pair of states. It is interpreted as the decoherence rate and is calculated as de- scribed in section 5.2.4. Each term (q 2 i +p 2 i )(q 2 j +p 2 j ) creates an additional potential that biases system’s dynamics to choose its evolution pathways such that they minimize the developed coherences (provided ij > 0). The sum in the right-hand side of Eq. 5.19, is called the CPF. 109 When the system accumulates coherence between the pair of statesi andj, the point in the phase space of the dynamical system is located uphill the bias potential. The gra- dient of the effective energy surface acts in the opposite direction and moves the phase space point away from the regions of large coherences. In this way, the system avoids developing large coherences. Such avoidance of the regions with large coherences in turn affects the average populations of the electronic states as functions of time, by changing the rates of transitions between them. Further details and illustrations of the CPF method can be found here. 34 5.2.3 Detailed balance and thermodynamic equilibrium Next, the detailed balance correction is applied based on the work of Adolfo Bastida, Alberto Requena and Beatriz Miguel. 36 They show that the coupling matrix elements between the basis states modified with the quantum correction factor 2 i~ _ Rd qc kj =i~ _ Rd kj 2 1 +e ~! kj KT !1 2 (5.20) results in transition rates that obey the detailed balance condition. Here, d (1) ij is repre- sented byd ij . In order to maintain the symmetry of the coupling matrix, keeping the overall Hamil- tonian Hermitian, the following modification is made, d sqc jk =d sqc kj = k d qc jk j d qc kj ; k>j (5.21a) d sqc jj = d jj ; j;k = 1;:::;N (5.21b) where k =jc k j is the absolute value of the wavefunction coefficient. 36 110 Overall, the following scaling for the off-diagonal Hamiltonian matrix elements was proposed: ~ H ij =H ij j s 2 1 +exp( E k B T ) H ji i s 2 1 +exp( E k B T ) (5.22) where E =E i E j The Hamiltonian considered in the original approach was real-valued. 36 In contrast, the NA coupling is imaginary. In such case, hermiticity requires that ~ H ij = ~ H ji . In order to maintain the NA Hamiltonian Hermitian, we modified the original scaling to the following form: H ij ! ~ H ij : Im[ ~ H ij ] = Im[H ij ] j s 2 1 +exp( E k B T ) i s 2 1 +exp( E k B T ) (5.23a) Re[ ~ H ij ] = Re[H ij ] (5.23b) Such scaling preserves the Hermitian character of the Hamiltonian and ensures that its off-diagonal elements vanish when state populations reach the Boltzmann distribution. The final form of the Hamiltonian, ~ H ij for the proposed Ehrenfest method with the decoherence and detailed balance corrections is given by Im[ ~ H ij ] = Im[H ij ] j s 2 1 +exp( E k B T ) i s 2 1 +exp( E k B T ) (5.24a) Re[ ~ H ij ] = Re[H ij ] + X i;j i6=j ij (q 2 i +p 2 i )(q 2 j +p 2 j ) (5.24b) The method, called Ehrenfest-DDB (decoherence detailed balance), has been imple- mented in the PYXAID software package. 36, 37 111 5.2.4 Pure-dephasing/decoherence time We estimate the decoherence time as the pure-dephasing time for pairs of initial and final states using the optical response theory formalism. 36, 45 Fluctuations in the energy gap between the electronic states due to nuclear motions are characterized by the auto- correlation function (ACF). C(t) =hE(t)E(0)i T (5.25) The brackets indicate canonical averaging. Fourier transformation of the ACF character- izes the phonon modes that couple to the electronic subsystem. The resulting spectrum is known as the influence spectrum, or spectral density. The pure-dephasing function is computed using the second-order cumulant expansion of the optical response func- tion. 36, 37 D cumu (t) =exp(g(t)) (5.26) whereg(t) is g(t) = t Z 0 d 1 1 Z 0 d 2 C( 2 ) (5.27) The pure-dephasing functions are fitted by a Gaussian, exp[0:5 t 2 ] to extract the pure-dephasing time constant, . The fitting results are presented in Table 5.1. Table 5.1 also contains the non-adiabatic couplings between pairs of adjacent electronic states. The decoherence rate, ij , that enters Eq. 5.24b is computed as the inverse of the pure- dephasing time,. 112 5.2.5 Simulation details In order to test the proposed method, we have modelled an electric-field induced charge and energy transfer in the system composed of a silver tip represented by a Ag 20 tetra- hedral cluster and a Zn(II)-etioporphyrin I molecule, as shown in Figure 5.1a. This test-system has been chosen for its complex electronic structure, including energetically well-separated molecular states and densely-packed levels of the metal tip. The sys- tem choice has been motivated by the experiments showing current-induced molecular luminescence. 38, 39 We employed the plane-wave density functional theory approach, as implemented in the QUANTUM ESPRESSO package, 46 to describe the electronic properties of the sys- tem. The Perdew-Burke-Ernzerhof (PBE) 47 exchange-correlation functional with norm- conserving pseudopotential and scalar-relativistic correction was used. The kinetic en- ergy cutoff was set to 60 Ry for geometry optimization and molecular dynamics and to 800 Ry for charge density calculations. The Grimme’s DFT-D2 dispersion correction was included. 48, 49 The molecule-tip distance was set to 6 ˚ A. Following the geometry optimization, a 1 ps long thermalization dynamics was performed at 300 K using the Andersen thermostat. Ground state MD was then performed to generate a 1 ps micro- canonical trajectory with a classical time-step of 1 fs to provide initial conditions for the excited state dynamics. The non-adiabatic dynamics were performed at 300 K using the Ehrenfest-DDB method as well as decoherence induced surface hopping (DISH). 33 We have performed the dynamics with DISH for comparison of the results of the two methods. DISH is a popular surface hopping method for NAMD that includes decoherence effects and satisfies detailed balance. Both methods have been implemented in the PYXAID pack- age 37, 47 and use the classical path approximation. Motivated by the fact that the nuclear 113 dynamics in the system under investigation involves vibrational motions near the equi- librium geometry, which are sampled sufficiently accurately during the 1 ps MD trajec- tory, we repeated the NA Hamiltonian multiple 1000 times in order to study long-time dynamics, up to 1 ns. This approximation cannot sample slow motions of the molecule with respect to the tip, however, it represents well faster nuclear motions that contribute most to the NA coupling. For each studied physical process (as described in Figure 5.1b) and method (Ehrenfest-DDB and DISH), 100 initial conditions were sampled from the 1ps MD tra- jectory and used as initial conditions for NAMD. The NA dynamics were initiated by instantly promoting the system to a particular excited state. Additionally for the DISH scheme, 100 realizations of the stochastic surface hopping process were sampled for each initial condition. The populations of different electronic states in Ehrenfest-DDB were computed as squares of the coefficients of expansion of the time-dependent wave- function in the adiabatic basis, Eq. 5.12. In DISH the state populations were obtained as fractions of trajectories in each adiabatic state. The atomic structure and orbital visu- alizations have been prepared with the XCRYSDEN software. 50 5.3 Results: Charge injection, cooling and recombina- tion in a metallic-tip/porphyrin system The optimized structure of the system under investigation 38, 39 is shown in Figure 5.1a, while Figure 5.1b summarizes the different electronic processes studied with the two methods, Ehrenfest-DDB and DISH. 33 First, the electron injection from the tip to the molecule is analyzed in a step-wise manner, with transition S4!S3 (initial phase of the electron transfer) and the complete electron injection process (transitions from S4 114 Figure 5.1: (a) Top view on the tip-porphyrin system. (b) The studied nonadiabatic processes and charge densities of the corresponding states. through S2) studied separately. Then, the combined electron injection and cooling pro- cesses (transitions from S4 through S1) are investigated. Finally, the full dynamics (S4 through S0) and charge recombination within the molecule (transition S1!S0) are mod- elled with the two methods. Table 5.1 reports the pure-dephasing times and the average absolute NA coupling for all pairs of states. Besides summarizing the key electronic pro- cesses, Figure 5.1b also reveals the nature of relevant electronic states. The molecular highest molecular orbital (HOMO) representing state S0 is localized on the metal cen- ter and complexing nitrogen atoms of the porphyrin molecule. The lowest unoccupied molecular orbital (LUMO) and LUMO+1, representing states S1 and S2, are composed of the carbon orbitals of the porphyrin. LUMO+2 (S3) is delocalized over the tip and the porphyrin, demonstrating strong top-molecule interaction in this state. LUMO+3 (S4) is fully localized on the silver tip. Figure 5.2a shows phonon-induced fluctuations in the energies of all electronic states during 500 fs of the MD trajectory, while Figure 5.2b zooms onto the excited states 115 Table 5.1: Pure-dephasing time (fs) (lower triangle) and time-averaged absolute value of nonadiabatic coupling (meV) (upper triangle) S0 S1 S2 S3 S4 S0 - 0.91 0.73 0.35 0.17 S1 5.2 - 34.99 7.22 2.25 S2 4.2 38.1 - 20.90 3.85 S3 4.4 31.5 84.6 - 3.63 S4 4.4 23.9 61.9 36.1 - and shows the energy gaps between nearest neighbor states. The gap between states S1 and S2 fluctuates most among the excited states, often reaching below 0.025 eV , which corresponds to k B T at room temperature. The NA coupling is also largest for this pair of states, Figure 5.2c. This is both because the energy gap between these states reaches small values, Figure 5.2b, and since the states are localized on and couple to high frequency modes of the porphyrin molecule, Figures 5.1b and 5.2e. Because most electronic transitions couple to one or two phonon modes, Figure 5.2e, with the excep- tion of the S2!S1 transition, the ACF functions of the energy gap fluctuations decay slowly, Figure 5.2d. The initial ACF value, corresponding to the average gap fluctuation squared, Eq. 5.25, is by far the largest for the S1!S0 transition. The large initial ACF value rationalizes 37 why the pure-dephasing function, Eq. 5.26, decays fastest for the S1!S0 transition. Figures 5.3-5.8 show evolution of excited state populations for various scenarios differing in the initially excited state and number of states involved. Figures 5.3-5.8 are 116 Figure 5.2: (a) Phonon-driven evolution of the electronic state energies, S0-S4. (b) Energy gaps between excited states, S1-S4. (c) Absolute value of nonadiabatic cou- pling. (d) Un-normalized autocorrelation functions. (e) Influence spectra. (f) Pure- dephasing functions for transitions between adjacent electronic states. The energy gap between states S1 and S2 drops below k B T=0.025eV at room temperature (part b). Transition between this pair of states exhibits the largest nonadiabatic coupling (part c) and couples to high frequency vibrations (part e), arising from the por- phyrin molecule that supports these states. The energy gap fluctuation is largest for the S0-S1 state pair (part d) which exhibits fastest dephasing (part f). 117 constructed the same way. Namely, the Ehrenfest-DDB data are shown in the upper row, and the DISH results are given in the lower row. The panels in the left column present the data up to 1 ns, while the panels in the right column show the same results for the first 100 ps of the evolution. Figure 5.3: Dynamics of electronic state populations displayed up to 1 ns (left pan- els) and focused on the first 100 ps (right panels) calculated within the Ehrenfest- DDB method (upper panels) and the DISH method (lower panels) for the partial electron injection process (S4!S3). Both methods produce similar timescales, with Ehrenfest-DBB marginally slower than DISH. Whereas the DISH method shows almost complete relaxation, the Ehrenfest-DDB scheme stops at about 6% short of complete relaxation. The dynamics for the electron injection process was started with the tip/porphyrin system excited to state S4. Figure 5.3 focuses on states S3 and S4. The population relaxation is very similar for Ehrenfest-DDB and DISH, with Ehrenfest-DDB showing 118 slightly slower dynamics. This result indicates that Ehrenfest-DDB can provide reason- able relaxation times. At long times, approaching the thermodynamic limit, the relax- ation slows down sharply for Ehrenfest-DDB, and the populations reach finite steady state values. In comparison, the DISH relaxation continues until all electronic popula- tion relaxes from the initial state S4 to state S3. The following analysis rationalizes the observed result. Note that detailed balance is achieved in the modified Ehrenfest method 36 by multiplying the coupling matrix element by a correction that decays to zero when the system approaches Boltzmann equilibrium. i.e., the closer the system is to equilibrium, the smaller the modified coupling, and the slower the dynamics. Specifically, consider a 2-level system with energies of states k andj equal toE and 0 respectively. For simplification, assume E k B T = 1. At equilibrium, the Boltzmann quantum populations of the two states are 2 k = 1 1 +e (5.28a) 2 j = e 1 +e (5.28b) wheree is Euler’s number. Let the populations 2 k and 2 j differ from their equilibrium values byx, such that k = r 1 1 +e +x (5.29a) j = r e 1 +e x (5.29b) 119 Plugging the values from Eqs. 5.29a and 5.29b into Eq. 5.24a we obtain Im[ ~ H jk ] =Im[H jk ] p 2e 1 +e p 1 +x(1 +e) s 1x 1 +e e !! (5.30) To estimate the value of the coupling between the two states whenx approaches 0, we Taylor expand the two functions ofx atx = 0 in the above equation and obtain p 1 +x(1 +e) s 1x 1 +e e ! = (1 +e) 2 2e x +O(x 3 ) (5.31) Therefore, Im[ ~ H jk ]/Im[H jk ]jcxj (5.32) wherec is a constant. Eq. 5.32 demonstrates that the detailed balance correction makes the coupling go to zero, as the system approaches equilibrium. For this reason, the original Ehren- fest modified for detailed balance method, 36 as well as Ehrenfest-DDB, require infinite time to achieve true Boltzmann distribution of populations. In practice, one usually requires fairly short-time dynamics in order to obtain the characteristic timescale. For example, the 100 ps data in Figure 5.3 is sufficient to perform an exponential fitting and deduce the quantum transition time. Although surface hopping methods such as DISH 38 can achieve the thermodynamic equilibrium more closely than the current ver- sion of Ehrenfest-DDB, they are more computationally expensive, since they involve stochastic sampling that has to be repeated multiple times for each initial condition. In comparison, Ehrenfest-DDB needs only one trajectory for each initial condition. Fur- ther, Ehrenfest-DDB requires only propagation of the Schr odinger equation (with the modified Hamiltonian), while DISH and related surface hopping methods require both 120 propagation of the Schr odinger equation as well as additional steps, making surface hopping more algorithmically complex. Figure 5.4: Dynamics of electronic state populations displayed up to 1 ns (left pan- els) and focused on the first 100 ps (right panels) calculated within the Ehrenfest- DDB method (upper panels) and the DISH scheme (lower panels) for the com- plete electron injection process (S4!S3!S2). Although the initial relaxation rate is marginally slower with Ehrenfest-DDB as compared to DISH, electronic popu- lations equilibrate at a similar state with both methods. Following the above two state dynamics, a slightly more complex system was con- sidered. The dynamics for the complete charge injection process, Figure 5.1b, including states S4, S3 and S2 was started with the system excited to state S4. Figure 5.4 shows the change in the electronic population of the three states versus time. The Ehrenfest- DDB and DISH schemes arrive at similar distributions of populations of the three states 121 at long time, left 5.4. On the shorter time-scale, right column of Figure 5.4, the dy- namics is slightly slower for Ehrenfest-DDB than DISH. Overall, the behavior of the two methods in the current example involving states S2-S4 is similar to the previous example involving states S3-S4, Figure 5.3. Figure 5.5: Dynamics of electronic state populations displayed up to 1 ns (left pan- els) and focused on the first 100 ps (right panels) calculated within the Ehrenfest- DDB scheme (upper panels) and the DISH scheme (lower panels) for the combined injection and cooling processes (S4!S1). The initial dynamics are similar and slightly slower with Ehrenfest-DDB compared to DISH. At the long-time limit the populations of states S1 and S2 reach different values. The different behavior of the two methods can be attributed to the facts that the energy gap between these states drops below kBT (Figure 5.2b) and that the nonadiabatic coupling is large (Figure 5.2c). 122 As the next step, we added the state S1 to the above set of states and studied the dynamics of the electron injection followed by the cooling process, including states S1, S2, S3 and S4. The system was excited to state S4 at the initial time for both Ehrenfest- DDB and DISH. The evolution of the state populations is shown in the upper row of Figure 5.5 for Ehrenfest-DDB and in the lower row of Figure 5.5 for DISH. Ehrenfest- DDB produces different long-time limit for populations of states S1 and S2 compared to DISH. The difference can be attributed to the energy gap and NA coupling values for the S1-S2 state pair, Figure 5.2b,c. In particular, the energy gap frequently drops below 0.025 eV , which corresponds tok B T at room temperature. The fluctuation of the energy gap from less thank B T to nearly 10 timesk B T makes the concept of equilibrium populations ill-defined. Further, the NA coupling between states S1 and S2 reaches over 20k B T , strongly mixing these states, and once again, making equilibrium populations ill-defined. The early time dynamics is slower in Ehrenfest-DDB than in DISH, similarly to the previous examples. Figure 5.6 shows the results of simulations involving states S1 and S2 only, with the initial population in S2. In DISH, the populations of S1 and S2 are similar, whether or not the higher energy states S3 and S4 are included. The difference in the S1 and S2 populations between Figures 5.5 and 5.6 is more pronounced in Ehrenfest-DDB, in- dicating that population equilibration is a more complex process showing non-trivial dependence on fluctuations in the energy gap, NA coupling and number of states. The Ehrenfest-DDB equations of motion are non-linear, with non-linearities arising in both detailed balance and decoherence terms. Therefore, the dependence of the ensuing dy- namics on the parameters of the time-dependent Hamiltonian can be rather involved and is hard to analyze analytically, beyond simple arguments such as those presented in Eqs. 5.28a-5.32. 123 Figure 5.6: Dynamics of electronic state populations displayed up to 1 ns (left pan- els) and focused on the first 100 ps (right panels) calculated within the Ehrenfest- DDB scheme (upper panels) and the DISH scheme (lower panels) just for the pair of states S1 and S2 that exhibit different long-time limits during the S4!S1 dynam- ics shown in Figure 5.5. The differences arise because the S1-S2 energy gap drops below k B T (Figure 5.2b) and the nonadiabatic coupling is large (Figure 5.2c). Next, we performed dynamics involving all states from S0 through S4, describing the complete set of processes from the excitation, to the electron injection, to the cooling, and to the recombination. As before, Figures 5.3-5.5, all of the electronic population was assigned to state S4 at the initial time. The results are shown in Figure 5.7. The major difference from the previous simulations is that the ground state S0 is strongly separated in energy from the excited states by about 1 eV , Figure 5.2a. Further, the NA 124 Figure 5.7: Dynamics of electronic state populations displayed up to 1 ns (left pan- els) and focused on the first 100 ps (right panels) calculated within the Ehrenfest- DDB scheme (upper panels) and the DISH scheme (lower panels) for the complete set of processes, including charge injection, cooling and recombination (S4!S0). The initial dynamics are similar and slightly slower in Ehrenfest-DDB compared to DISH. At long times, the Ehrenfest-DDB dynamics slows down considerably, with the larges differences observed for the populations of states S0, S1 and S2. See Figures 5.5 and 5.6 for additional insights. coupling between S0 and all excited states is at least an order of magnitude smaller (0.2- 0.9 meV) than the NA coupling between the excited states (2-35 meV), Table 5.1, and pure-dephasing is an order of magnitude shorter (4-5 fs vs. 25-85 fs). The large energy gap, and small NA couplings and pure-dephasing times lead to very slow population of the ground state S0. Even the 1 ns simulations are not sufficient to reach equilibrium populations in both Ehrenfest-DDB and DISH. The 100 ps data, right panels in Figure 125 5.7, are similar for the two methods, with the exception of S1-S2 equilibration, discussed above, Figures 5.5 and 5.6. The ground state S0 population reaches about 20% within 100 ps in both Ehrenfest-DDB and DISH, indicating that both methods predict similar charge recombination times. At the long time, left panels of Figure 5.7, the Ehrenfest- DDB dynamics slows down considerably, even though the population of the S0 state continues to grow, as required by thermal equilibration. Figure 5.8: Dynamics of electronic state populations displayed up to 1 ns (left pan- els) and focused on the first 100 ps (right panels) calculated within the Ehrenfest- DDB scheme (upper panels) and the DISH scheme (lower panels) for the recom- bination process (S1!S0). Both DISH and Ehrenfest-DDB require more than 1 ns to reach equilibrium, owing to the large energy gap (Figure 5.2a) and small nonadiabatic coupling (Figure 5.2c and Table 5.1) between the states. Ehrenfest- DDB approaches equilibrium faster with two states (S1!S0) than with five states (S4!S0), Figure 5.7. 126 Figure 5.8 focuses on charge recombination due to non-radiative relaxation from the first excited state S1 to the ground state S0. The short time dynamics is similar between Ehrenfest-DDB and DISH. Judging by the 100 ps data, right panels in Figure 5.8, Ehrenfest-DDB even shows faster relaxation than DISH. The population of the ground state reaches 40% in Ehrenfest-DDB and slightly less than 30% in DISH. The ground state population reaches 50% at about the same time with the two methods, while the longer time decay is significantly slower in Ehrenfest-DDB, left panels of Figure 5.8. The current work demonstrates that both detailed balance and decoherence can be introduced into the Ehrenfest method, making is applicable for simulation of large re- alistic systems and long timescales. At the same time, further studies are required to investigate how the detailed balance modification of the Ehrenfest method affects the quantum transition rates and approach to equilibrium, with and without inclusion of de- coherence. There exist several quantum correction factors that account for detailed bal- ance. 33 Bastida et al. 36 employedQ(!) = 2 1+exp(~!) . Other choices includeQ(!) = ~! 1+exp(~!) and Q(!) = ~!, as well as more complicated versions. Each quantum correction will produce a different scaling of the coupling matrix elements, leading to different rates of quantum transitions. Further, Kleinekath ofer and co-workers 51 pointed out that the original approach does not reproduce the high-temperature limit, and pro- posed an additional normalization. They divided the scaled coupling by the difference i j , which restores the high-temperature limit of equally populated states, but causes divergence as the system approaches this limit and i j ! 0. A more detailed analy- sis and further tests are needed to order to establish which type of correction works best for a particular class of problems. 127 5.4 Conclusions We have presented the Ehrenfest-DDB method for NAMD simulations that accounts for decoherence effects and detailed balance between transitions upward and downward in energy. Both features are essential for modeling excited state dynamics in large sys- tems on long time-scales. The methodology is based on the mixed quantum-classical Ehrenfest approach that couples quantum mechanical expectation values to classical variables. This simplest quantum-classical approximation is rooted in the Ehrenfest the- orem and can be derived in different ways. 21, 36, 45, 51 Decoherence effects are missing in the Ehrenfest method because the \bath" is treated classically, and the evolution of the quantum subsystem remains unitary even when coupled to an environment. The Ehrenfest approach cannot describe quantum-classical thermalization because it is a mean-field method and lacks correlations of classical trajectories with different quan- tum states. Decoherence effects are incorporated into the Ehrenfest approach via an additional term that penalizes development of coherences during evolution of the quan- tum degrees of freedom. In the classically mapped Hamiltonian, the coherence penalty term adds an energy penalty to the total energy. 34 The decoherence timescale is esti- mated as the pure-dephasing time of the optical response theory. 52, 53 Detailed balance is achieved by scaling the off-diagonal matrix elements of the Hamiltonian with a quantum correction factor known in semiclassical approximations to quantum time-correlation functions. 54 The detailed balance correction is a function of temperature and energy gap between a given pair of states. Both modifications modify the Schr odinger equa- tion making it non-linear. At the same time, the key feature of the Ehrenfest method, which that describes quantum-classical dynamics fully deterministically and with wave- functions, is preserved. In contrast, other approaches such as surface hopping or master 128 equations, either introduce stochastic processes (hops) that require extensive sampling, or deal with N 2 dimensional density matrices rather than N dimensional wavefunctions. Thus, the Ehrenfest-DDB method is computationally more efficient than the alternative techniques. We illustrated the Ehrenfest-DDB method by performing calculations on a hybrid molecule-tip system studied experimentally. The system undergoes several important processes that are typical of modern applications of nanoscale materials, including charge injection, relaxation and recombination. The Ehrenfest-DDB results are com- pared to calculations performed using DISH, which is a popular surface hopping tech- nique with decoherence and detailed balance effects. The Ehrenfest-DDB approach gives timescales that are similar and slightly longer than those of DISH. In the long-time limit, the Ehrenfest-DDB dynamics slow down considerably, due to the main feature of the detailed balance correction that drives the off-diagonal matrix elements of the Hamil- tonian to zero as the system approaches Boltzmann equilibrium. Nevertheless, correct timescales can be obtained with relatively short simulations. The current version of Ehrenfest-DDB uses one of the several possible quantum correction factors responsible for detailed balance, 55 following the original publication. 36 The performance of other detailed balance corrections 51, 56 deserve further studies that are planned in near future. The Ehrenfest-DDB method, being an extension of the Ehrenfest dynamics, is easy to implement and computationally efficient. Avoiding the need to generate ensembles of stochastic trajectories, it is capable of describing thermal relaxation and loss of quan- tum coherence, making it appropriate for calculations on large condensed matter and nanoscale systems. 129 Chapter 5 References [1] D. G. Evans, A. Nitzan, and M. A. Ratner, The Journal of Chemical Physics 108, 6387 (1998). [2] John Morelli and Sharon Hammes-Schiffer, Chemical Physics Letters 269, 161 (1997). [3] Jian-Yun Fang and Sharon Hammes-Schiffer, The Journal of Chemical Physics 110, 11166 (1999). [4] Anirban Hazra, Alexander V . Soudackov, and Sharon Hammes-Schiffer, The Jour- nal of Physical Chemistry Letters 2, 36 (2011), PMID: 26295211. [5] Puja Goyal and Sharon Hammes-Schiffer, The Journal of Physical Chemistry Let- ters 6, 3515 (2015), PMID: 26275870. [6] I. Presiado, Y . Erez, R. Gepshtein, and D. Huppert, The Journal of Physical Chem- istry C 114, 3634 (2010). [7] Jing Huang, Likai Du, Jun Wang, and Zhenggang Lan, The Journal of Physical Chemistry C 119, 7578 (2015). [8] Xing Gao, Qian Peng, Yingli Niu, Dong Wang, and Zhigang Shuai, Phys. Chem. Chem. Phys. 14, 14207 (2012). [9] Mirjana Eckert-Maksi´ c and Ivana Antol, The Journal of Physical Chemistry A 113, 12582 (2009), PMID: 19603774. [10] Run Long, Niall J. English, and Oleg V . Prezhdo, The Journal of Physical Chem- istry Letters 5, 2941 (2014), PMID: 26278240. [11] Haiming Zhu, Ye Yang, Kim Hyeon-Deuk, Marco Califano, Nianhui Song, Youwei Wang, Wenqing Zhang, Oleg V . Prezhdo, and Tianquan Lian, Nano Letters 14, 1263 (2014), PMID: 24359156. [12] Jin Liu and Oleg V . Prezhdo, The Journal of Physical Chemistry Letters 6, 4463 (2015), PMID: 26505613. [13] H. Marciniak, M. Fiebig, M. Huth, S. Schiefer, B. Nickel, F. Selmaier, and S. Lochbrunner, Phys. Rev. Lett. 99, 176402 (2007). [14] Henning Marciniak, Igor Pugliesi, Bert Nickel, and Stefan Lochbrunner, Phys. Rev. B 79, 235318 (2009). 130 [15] Bruno Ehrler, Mark W. B. Wilson, Akshay Rao, Richard H. Friend, and Neil C. Greenham, Nano Letters 12, 1053 (2012), PMID: 22257168. [16] Akshay Rao, Mark W. B. Wilson, Justin M. Hodgkiss, Sebastian Albert-Seifried, Heinz B¨ assler, and Richard H. Friend, Journal of the American Chemical Society 132, 12698 (2010), PMID: 20735067. [17] Spiridoula Matsika and Pascal Krause, Annual Review of Physical Chemistry 62, 621 (2011), PMID: 21219147. [18] Sharon Hammes-Schiffer and Alexei A. Stuchebrukhov, Chemical Reviews 110, 6939 (2010), PMID: 21049940. [19] Alexey V . Akimov, Amanda J. Neukirch, and Oleg V . Prezhdo, Chemical Reviews 113, 4496 (2013), PMID: 23627277. [20] Chengbo Guan, Ning Wu, and Yang Zhao, The Journal of Chemical Physics 138, 115102 (2013). [21] Lars V . Sch¨ afer, Gerrit Groenhof, Martial Boggio-Pasqua, Michael A. Robb, and Helmut Grubm¨ uller, PLOS Computational Biology 4, 1 (2008). [22] John C. Tully, The Journal of Chemical Physics 93, 1061 (1990). [23] Priya V . Parandekar and John C. Tully, The Journal of Chemical Physics 122, 094102 (2005). [24] Linjun Wang, Alexey Akimov, and Oleg V . Prezhdo, The Journal of Physical Chemistry Letters 7, 2100 (2016), PMID: 27171314. [25] Neil Shenvi, Sharani Roy, and John C. Tully, The Journal of Chemical Physics 130, 174107 (2009). [26] Tammie Nelson, Sebastian Fernandez-Alberti, Adrian E. Roitberg, and Sergei Tre- tiak, Chemical Physics Letters 590, 208 (2013). [27] Eric R. Bittner and Peter J. Rossky, The Journal of Chemical Physics 103, 8130 (1995). [28] Benjamin J. Schwartz, Eric R. Bittner, Oleg V . Prezhdo, and Peter J. Rossky, The Journal of Chemical Physics 104, 5942 (1996). [29] Oleg V . Prezhdo, The Journal of Chemical Physics 111, 8366 (1999). [30] Amber Jain, Ethan Alguire, and Joseph E. Subotnik, Journal of Chemical Theory and Computation 12, 5256 (2016), PMID: 27715036. 131 [31] Shu Chun Cheng, Chaoyuan Zhu, Kuo Kan Liang, Sheng Hsien Lin, and Don- ald G. Truhlar, The Journal of Chemical Physics 129, 024112 (2008). [32] Ross E. Larsen, Michael J. Bedard-Hearn, and Benjamin J. Schwartz, The Journal of Physical Chemistry B 110, 20055 (2006), PMID: 17020394. [33] Heather M. Jaeger, Sean Fischer, and Oleg V . Prezhdo, The Journal of Chemical Physics 137, 22A545 (2012). [34] Alexey V . Akimov, Run Long, and Oleg V . Prezhdo, The Journal of Chemical Physics 140, 194107 (2014). [35] Priya V . Parandekar and John C. Tully, Journal of Chemical Theory and Compu- tation 2, 229 (2006), PMID: 26626509. [36] Adolfo Bastida, Carlos Cruz, Jos´ e Z´ u˜ niga, Alberto Requena, and Beatriz Miguel, 417, 53 (2006). [37] Alexey V . Akimov and Oleg V . Prezhdo, Journal of Chemical Theory and Compu- tation 9, 4959 (2013), PMID: 26583414. [38] X. H. Qiu, G. V . Nazin, and W. Ho, Science 299, 542 (2003). [39] Joanna Jankowska and Oleg V . Prezhdo, The Journal of Physical Chemistry Letters 9, 3591 (2018), PMID: 29897769. [40] Eyal Neria and Abraham Nitzan, The Journal of Chemical Physics 99, 1109 (1993). [41] Pavel Jungwirth and R. Benny Gerber, The Journal of Chemical Physics 104, 5803 (1996). [42] Todd J. Martinez, M. Ben-Nun, and R. D. Levine, The Journal of Physical Chem- istry 100, 7884 (1996). [43] Alexey V . Akimov and Oleg V . Prezhdo, Journal of Chemical Theory and Compu- tation 10, 789 (2014), PMID: 26580053. [44] Folkmar A. Bornemann, Peter Nettesheim, and Christof Sch¨ utte, The Journal of Chemical Physics 105, 1074 (1996). [45] Sergei Egorov and J.L. Skinner, 293, 469 (1998). [46] Paolo Giannozzi et al., Journal of Physics: Condensed Matter 21, 395502 (2009). [47] John P. Perdew, Kieron Burke, and Matthias Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996). 132 [48] Grimme Stefan, Journal of Computational Chemistry 27, 1787. [49] Barone Vincenzo, Casarin Maurizio, Forrer Daniel, Pavone Michele, Sambi Mauro, and Vittadini Andrea, Journal of Computational Chemistry 30, 934. [50] A Kokalj, Journal of molecular graphics & modelling 17, 176 (1999). [51] Mortaza Aghtar, J¨ org Liebers, Johan Str¨ umpfer, Klaus Schulten, and Ulrich Kleinekath¨ ofer, The Journal of Chemical Physics 136, 214101 (2012). [52] N. F. Mott, Mathematical Proceedings of the Cambridge Philosophical Society 27, 553 (1931). [53] John B. Delos, Walter R. Thorson, and Stephen K. Knudson, Phys. Rev. A 6, 709 (1972). [54] Klaus Hepp, Comm. Math. Phys. 35, 265 (1974). [55] Oleg V . Prezhdo and Vladimir V . Kisil, Phys. Rev. A 56, 162 (1997). [56] S. (Shaul) Mukamel, Principles of nonlinear optical spectroscopy. Oxford : Ox- ford University Press, 1999. 133 Chapter 6: Future Work In the preceding chapters, we have seen the applications of adiabatic and non- adiabatic molecular dynamics methods to study electronic dynamical processes in sev- eral nanoscale and condensed phase systems. The chemical process to be described dictates the method used for the theoretical analysis as different methods/algorithms include differing levels of the necessary physics for any given problem. The work pre- sented in Chapter 2 utilizes the self-consistent charge density functional tight binding (SCC-DFTB) method due to the large number of atoms in the problem. 1, 2 However, a major challenge in the description of such large systems is their limited flexibility in terms of choice of model chemistry. To perform time-domain atomistic simulations in such large systems with a reasonable computational cost, DFTB is an excellent choice. However, DFTB, being an approximate method has an intrinsic level of accuracy which may or may not be sufficient for the description of the problem at hand. One must be cautious and DFTB results should be tested before application to new types of structures and reactions against higher level methods. 3 Also, it may not be possible to obtain the desired level of accuracy for all the properties of interest simultaneously. A different set of parameters may be required for computing different properties and such testing needs to be done before application to different chemical environments. 134 In Chapter 3, we performed ab initio molecular dynamics simulations to study the formation of He 2 in helium droplets. We employed the equation-of-motion coupled- cluster method with single and double substitutions for excitation energies (EOM-EE- CCSD), and the d-aug-cc-pVDZ basis to describe the diffuse wave functions of Rydberg states. The dynamics were performed on the first triplet adiabatic surface of the Helium systems studied. This results in exciton hopping between very distant sites especially in the larger 10-atom cases as a consequence of the Born Oppenheimer approximation. A complete theoretical description should include non-adiabatic transitions among mul- tiple adiabats. One way to include non-adiabatic effects is using the surface hopping approach. This is especially important in the 10-atoms systems where there is a dense manifold of closely lying electronic states and the rate of exciton hopping between sites is grossly overestimated in the adiabatic framework. In Chapter 4, we discussed the Trivial Crossing Problem of FSSH 4 and it’s effect on the electronic dynamics of a photoinduced spirobifluorene. We also presented the Min- Cost state tracking algorithm 5–7 which successfully treats the trivial crossing problem of FSSH. A few other strategies to deal with the trivial crossing problem of FSSH have been proposed. Notable among them are the self-consistent fewest switches surface hopping (SC-FSSH) 8 by Linjun Wang and Oleg Prezhdo and the direct trajectories with surface hopping (DTSH) 9 by G. Granucci and co-workers. Both the latter corrections are based on dynamics in the diabatic basis since the nonadiabatic couplings are not involved in the electronic dynamics when a diabatic basis is used. A comparison of all three strategies with FSSH as well as with each other will shed more light on the trivial crossing issue of FSSH and provide a better understanding of how each corrected algorithm treats such crossings and which type of correction works best for dynamics in systems with trivial crossings. 135 In chapter 5, we proposed a new semiclassical method based on Ehrenfest method that accounts for decoherence and detailed balance and showed it’s successful applica- tion to a hybrid molecule-tip system. This work demonstrates that both detailed balance and decoherence can be introduced into the Ehrenfest method, making it applicable for simulation of large realistic systems and long timescales. However, further studies are required to investigate how the detailed balance modification of the Ehrenfest method affects the quantum transition rates and approach to equilibrium, with and without inclu- sion of decoherence. As both the detailed balance and the decoherence corrections in- troduce non-linear terms to the Schr odinger equation, such studies will provide a clearer picture of the effects of both the correction terms separately as well as together to the overall dynamics. It should be noted that there exist several quantum correction factors that account for detailed balance. 10 In our implementation, we followed Bastida et al. 11 who employed Q(!) = 2 1 +e ~! (6.1) Other QCFs which satisfy the detailed balance condition include Q(!) = ~! 1e ~! (6.2) and Q(!) =e ~!=2 (6.3) as well as more sophisticated versions. Each quantum correction will produce a different scaling of the coupling matrix elements, leading to different rates of quantum transitions. Thus, a comprehensive study of these different QCFs is essential to understand their 136 behavior at small and large ! as well as compare their accuracy for several different problems. Further, Kleinekath ofer and co-workers 12 pointed out that the symmetry corrected couplings for detailed balance as used by Bastida et al. and us does not reproduce the high-temperature limit. They proposed a modified scheme with an additional normal- ization d sqc jk =d sqc kj = k d qc jk j d qc kj k j ; E j >E k and j 6= k d sqc jk =d sqc kj = (d qc jk +d qc kj )=2; E j >E k and j = k d sqc jj =d jj (6.4) They divided the scaled coupling by the difference j k , which restores the high- temperature limit of equally populated states, but causes divergence as the system ap- proaches this limit and j k ! 0. To avoid the singularity, they set the coupling value for equal coefficients at a finite value equal to the average of the upward and downward quantum corrections. A more detailed analysis and further tests are needed in order to establish which type of correction works best for a particular class of problems. 137 Chapter 6 References [1] G. Seifert, The Journal of Physical Chemistry A 111, 5609 (2007), PMID: 17439198. [2] M. Elstner, D. Porezag, G. Jungnickel, J. Elsner, M. Haugk, Th. Frauenheim, S. Suhai, and G. Seifert, Phys. Rev. B 58, 7260 (1998). [3] M. Elstner, Theoretical Chemistry Accounts 116, 316 (2006). [4] John C. Tully, The Journal of Chemical Physics 93, 1061 (1990). [5] Tammie Nelson, Sebastian Fernandez-Alberti, Vladimir Chernyak, Adrian E. Roit- berg, and Sergei Tretiak, The Journal of Chemical Physics 136, 054108 (2012). [6] Tammie Nelson, Sebastian Fernandez-Alberti, Vladimir Chernyak, Adrian E. Roit- berg, and Sergei Tretiak, The Journal of Physical Chemistry B 115, 5402 (2011), PMID: 21218841. [7] Tammie Nelson, Sebastian Fernandez-Alberti, Adrian E. Roitberg, and Sergei Tre- tiak, Accounts of Chemical Research 47, 1155 (2014), PMID: 24673100. [8] Linjun Wang and Oleg V . Prezhdo, The Journal of Physical Chemistry Letters 5, 713 (2014), PMID: 26270842. [9] G. Granucci, M. Persico, and A. Toniolo, The Journal of Chemical Physics 114, 10608 (2001). [10] Sergei Egorov and J.L. Skinner, 293, 469 (1998). [11] Adolfo Bastida, Carlos Cruz, Jos´ e Z´ u˜ niga, Alberto Requena, and Beatriz Miguel, 417, 53 (2006). [12] Mortaza Aghtar, J¨ org Liebers, Johan Str¨ umpfer, Klaus Schulten, and Ulrich Kleinekath¨ ofer, The Journal of Chemical Physics 136, 214101 (2012). 138

Abstract (if available)

Abstract Adiabatic as well as non-adiabatic molecular dynamics methods are essential tools to simulate many important physical processes. Molecular dynamics simulations provide mechanistic information that often cannot be detected by experiment alone and also serve as a predictive tool to stimulate new research. However, sophisticated modeling of ultrafast processes, especially those involving excited states, is intellectually challenging and computationally demanding. The work presented here provides several diverse examples of challenging gas phase and condensed phase ultrafast photoinduced processes. For each system, we offer a detailed account of the computational modeling scheme employed for an accurate description of the chemical and physical properties of interest. At the same time, we highlight the limitations of the particular dynamics methods used. ❧ In Chapter 1, we start with an overview of the adiabatic and the non-adiabatic dynamics. This is followed by a description of the primary non-adiabatic dynamics methods, the Ehrenfest dynamics and the Fewest Switches Surface Hopping method. We also describe the major limitations of the two approaches, especially the neglect of decoherence effects among others. ❧ This is followed by an example of the application of adiabatic molecular dynamics in Chapter 2 to study the atomistic origin of room temperature coherence between the light-hole and heavy-hole excitations in two-dimensional CdSe and CdSe/CdZnS core/shell nanoplatelets, synthesized and studied recently by optical spectroscopy. The molecular dynamics calculations allowed us to mimic closely the experimental systems and measurements, and to elucidate the role of defects that are inevitably present in experiment. Chapter 3 provides another application of quantum-classical molecular dynamics to study the conversion of He(2³S) to He₂(a₃Σᵤ⁺) in large helium nanodroplets. He ⃰₂ is a likely participant in the production of He₄⁺ on surfaces of large helium droplets which shows an anomalously intense peak in mass spectra of large helium droplets. Ab initio molecular dynamics with highly accurate and robust equation-of-motion coupled-cluster (EOM-CC) theory, showed rapid and efficient conversion of He ⃰ to He ⃰₂ in systems composed of small numbers of Helium atoms. ❧ The Born-Oppenheimer dynamics is followed by the non-adiabatic dynamics with electronic transitions in chapters 4 and 5. Chapter 4 analyzes the fewest switches surface hopping (FSSH) algorithm and brings attention to one of its limitations, the trivial crossing problem. It provides a comparison of the standard fewest switches surface hopping (FSSH) algorithm with the Min-Cost adiabatic state identity tracking algorithm, which treats the trivial crossing problem of FSSH. The performance of the two methods is investigated for the internal conversion and energy transfer dynamics of a weakly coupled spiro-linked polyfluorene (spirobifluorene) system. The results reveal that FSSH fails to identify trivial crossings, therefore, the system remains on the upper adiabatic surface at such crossings and reaches the lowest excited state only via non-adiabatic transitions resulting in a longer relaxation time-scale as compared to the Min-Cost method. ❧ Chapter 5 addresses the other primary mixed quantum classical dynamics method, the Ehrenfest method. Ehrenfest method does not produce the Boltzmann equilibrium populations of the quantum states at a finite temperature. Additionally, both the Ehrenfest MD and the FSSH techniques are fully quantum coherent, they do not account for decoherence. In this chapter, we present a semiclassical approach for non-adiabatic molecular dynamics based on the Ehrenfest method with corrections for decoherence effects and detailed balance. The decoherence effects are described via a coherence penalty functional that drives dynamics away from regions in Hilbert space characterized by large values of coherences. Detailed balance is described by modification of the off-diagonal matrix elements with a quantum correction factor used in semiclassical approximations to quantum time-correlation functions. Both decoherence and detailed balance corrections introduce non-linear terms to the Schrodinger equation. At the same time, the simplicity of the Ehrenfest method, which involves fully deterministic dynamics and requires a single trajectory for each initial condition, is preserved. In contrast, surface hopping techniques contain stochastic elements and require averaging over multiple realizations of stochastic process for each initial condition. The method is applied to a hybrid nanoscale system consisting of a fluorophore molecule and an STM tip studied experimentally. By investigating current-induced charge injection, cooling, and recombination, we demonstrate that the Ehrenfest-decoherence-detailed-balance (Ehrenfest-DDB) method produces timescales that are similar to those obtained with the decoherence induced surface hopping (DISH), which is a popular nonadiabatic molecular dynamics technique used with many condensed matter and nanoscale systems. The Ehrenfest-DDB technique provides efficient means to study quantum dynamics in large, realistic systems. ❧ Finally, in chapter 6, we outline future directions and discuss further investigations of the developed computational scheme to establish the range of problems it can be applied to with useful accuracy.

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

Photoexcitation and nonradiative relaxation in molecular systems: methodology, optical properties, and dynamics

PDF

Theoretical modeling of nanoscale systems: applications and method development

PDF

Non-adiabatic molecular dynamic simulations of charge carrier dynamics in two-dimentional transition metal dichalcogenides

PDF

Towards robust dynamical decoupling and high fidelity adiabatic quantum computation

PDF

Development and application of robust many-body methods for strongly correlated systems: from spin-forbidden chemistry to single-molecule magnets

PDF

Controlling electronic properties of two-dimensional quantum materials: simulation at the nexus of the classical and quantum computing eras

PDF

Out-of-equilibrium dynamics of inhomogeneous quantum systems

PDF

Spectroscopic signatures and dynamic consequences of multiple interacting states in molecular systems

PDF

Electronic structure analysis of challenging open-shell systems: from excited states of copper oxide anions to exchange-coupling in binuclear copper complexes

PDF

Quantum molecular dynamics and machine learning for non-equilibrium processes in quantum materials

PDF

Electronically excited and ionized states in condensed phase: theory and applications

PDF

Molecular orientation from sum frequency generation spectroscopy: case study of water and methyl vibrations

PDF

Photodissociation dynamics of atmospherically relevant small molecules in molecular beams

PDF

Molecular-scale studies of mechanical phenomena at the interface between two solid surfaces: from high performance friction to superlubricity and flash heating

PDF

Reactive and quantum molecular dynamics study of materials: from oxidation and amorphization to interfacial charge transfer

PDF

Dynamics of water in nanotubes: liquid below freezing point and ice-like near boiling point

PDF

Molecular simulations of water and monovalent ion dynamics in the electroporation of phospholipid bilayers

PDF

Multiscale and multiresolution approach to characterization and modeling of porous media: From pore to field scale

PDF

Kinetic Monte Carlo simulations for electron transport in redox proteins: from single cytochromes to redox networks

Asset Metadata

Creator Nijjar, Parmeet Kaur (author)

Core Title Adiabatic and non-adiabatic molecular dynamics in nanoscale systems: theory and applications

School College of Letters, Arts and Sciences

Degree Doctor of Philosophy

Degree Program Chemistry (Chemical Physics)

Publication Date 02/14/2019

Defense Date 01/18/2019

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag adiabatic dynamics,Boltzmann equilibrium,decoherence,Ehrenfest method,fewest switches surface hopping,non-adiabatic dynamics,OAI-PMH Harvest

Format application/pdf (imt)

Language English

Contributor Electronically uploaded by the author (provenance)

Advisor Prezhdo, Oleg (committee chair), Dawlaty, Jahan (committee member), El-Naggar, Moh (committee member), Wittig, Curt (committee member)

Creator Email parmeet.iitkgp@gmail.com,pnijjar@usc.edu

Permanent Link (DOI) https://doi.org/10.25549/usctheses-c89-121057

Unique identifier UC11675102

Identifier etd-NijjarParm-7068.pdf (filename),usctheses-c89-121057 (legacy record id)

Legacy Identifier etd-NijjarParm-7068.pdf

Dmrecord 121057

Document Type Dissertation

Format application/pdf (imt)

Rights Nijjar, Parmeet Kaur

Type texts

Source University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Access Conditions The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the a...

Repository Name University of Southern California Digital Library

Repository Location USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA