6728
J. Phys. Chem. B 2001, 105, 6728-6737
Influence of Solvation Environment on Excited State Avoided Crossings and Photodissociation Dynamics† N. Yu, C. J. Margulis,‡ and D. F. Coker*,§ Department of Chemistry, Boston UniVersity, 590 Commonwealth AVenue, Boston, Massachusetts 02215 ReceiVed: March 8, 2001; In Final Form: May 9, 2001
Photodissociation dynamics of I2 in rare gas liquids and solid matrixes have been the focus of much experimental and theoretical research. Experimental studies of both A and B state photoexcitation leading to dissociation have been thoroughly investigated by Apkarian and Schwentner, showing very dissimilar behavior between solid and liquid phases. B state excitation of I2 in liquid rare gases results in fast dissociation in a few hundred femtoseconds, while the B state lifetime in the solid is measured to be on the order of a few picoseconds. The current paper addresses various important issues related to the implementation of the nonadiabatic dynamics methods in the situation where quasi-nonavoided crossings take place and play a crucial role in determining the dynamics. Reliable treatment of quasi-nonavoided crossings in nonadiabatic MD calculations enables the accurate computation of pump-probe signals that reproduce experimental results. Analysis of trajectory data provides a detailed understanding of the interactions between local solvation environment symmetries and electronic coupling symmetries that are fundamentally important in governing nonadiabatic relaxation in condensed phases. The paper also explores the fundamental issue of the existence of approximate “selection rules” for nonadiabatic electronic relaxation arising from symmetries of electronic state coupling and the local symmetries present in ordered or even disordered condensed phase environments.
1. Introduction In the gas phase the lifetime of I2 prepared in its electronic B state by photoexcitation is many orders of magnitude longer than in the condensed phase. This is due to the fact that the B state is bound and very small couplings exist between it and the different electronic states that cross this state at various bond lengths. (Figure 1 displays models of these crossing surfaces which are some of the potential curves used in our calculations.) These are strictly speaking nonavoided crossings; i.e., they involve pairs of states which, when labeled by symmetry, have energies that can become degenerate and thus cross one another at certain nuclear configurations. In the condensed phase, however, interactions between the solvent and each of the iodine atoms breaks the symmetry and induces couplings among these otherwise uncoupled states. One would imagine that at very low solvent density these interactions would be weak, and that the symmetry eigenstates of the gas phase solute would be a very good approximation to the actual eigenstates of the system. Increasing the density of the solvent can only create stronger interactions with the solute and give rise to appreciable offdiagonal elements in the Hamiltonian matrix written in the gas phase eigenstate basis set. An exception to this is when the density is increased under the constraint of preserving some of the symmetries of the gas phase molecular electronic states in the molecule plus condensed phase system. Although this is an idealization, this situation is very much like what happens in solid rare gas matrixes doped with I2 molecules. The environment around the I2 at very low temper†
Part of the special issue “Bruce Berne Festschrift”. Present address: Department of Chemistry, Columbia University, New York, NY 10027. § Schlumberger Visiting Professor of Theoretical Chemistry, Department of Chemistry, Cambridge University, Lensfield Road, Cambridge CB2 1EW, U.K. ‡
Figure 1. Isolated molecule excited state potential energy (y-axis in cm-1) curves as functions of bond length (x-axis in Ångstroms). The B-state curve which is the lowest energy state from the (J ) 3/2, J′ ) 1/ ) spin-orbit excited state manifold of I is crossed by various states 2 2 from the (3/2, 3/2) spin-orbit ground state manifold. In this paper we explore the influence of solvation environment on B-state crossings with the following states: Closest crossing 1Πu, next furthest 3Π2g, and the next furthest 3Πg.
atures is highly symmetrical, and the electronic states that cross the B state as the I2 molecule vibrates are quasi-degenerate with this state. Quasi-degenerate states are responsible for quasinonavoided crossings and these in turn result in the longer lifetime of the B state in solids as compared to the relatively short excited electronic state lifetimes observed in liquids.1 In experiments these quasi-nonavoided crossings thus manifest themselves in the observation that the B state lifetime in solid rare gases is more than an order of magnitude longer than in typical liquid environments. The short lifetime in the liquid of
10.1021/jp0108925 CCC: $20.00 © 2001 American Chemical Society Published on Web 06/23/2001
J. Phys. Chem. B, Vol. 105, No. 28, 2001 6729
Photodissociation Dynamics of I2 course arises due to the highly asymmetrical local environment expected in a fluid. The lifetime of an I2 molecule excited to its B state in the solid is, however, many orders of magnitude shorter than in the gas phase due to the fact that as the excited state molecule vibrates in the highly symmetrical environment of a rare gas solid, it will slowly transfer energy to the surroundings or convert some of its vibrational energy into its own center-of-mass translational motion and move in an asymmetrical way relative to its surroundings, thus making the local environment it experiences less symmetrical and so lifting quasi-degeneracies between its electronic states and resulting in real nonadiabatic couplings between these states. The aim of this paper is to explore how the evolving local symmetries and asymmetries of a condensed phase environment can influence the coupling between molecular electronic states of different symmetries and thus effect the rate of nonadiabatic electronic relaxation of excited molecules in these dynamical environments.1,2 We will thus employ a surface hopping approach that can reliably treat systems with complex patterns of avoided and quasi-nonavoided crossings that arise when molecules with dense manifolds of excited electronic states are solvated in condensed phases. Quasi-nonavoided crossings occur in calculations when two or more electronic energy states appear to switch their symmetry in less than one electronic time step no matter how small this electronic step is chosen (even for steps of the order of 10-5 fs). In other words, the coupling between states is very small and localized in space in such a way that it is of the order of the error in the computation of the energy eigenvalues of the system. When such a situation arises, extreme care must be taken to ensure that the nonadiabatic algorithm follows states of appropriate symmetry rather than making nonphysical transitions. In the method section of this paper we address this issue and devise a detection scheme to cope with these quasinonavoided crossings so that quasi-classical surface hopping dynamics will be able to follow the correct symmetry of the molecular electronic wave function during these instantaneous state switchings. The reliable treatment of quasi-nonavoided crossings and quasi-degeneracy is in principle important not only in highly symmetrical environments, such as in solids, as discussed above. As molecular dissociation takes place, the system evolves from an arrangement of avoided and quasi-nonavoided molecular state crossings, to a different pattern of quasi-nonavoided crossing states characteristic of environmentally perturbed atoms. Even in disordered liquid environments there can thus be problems with state degeneracy that arise, for example, when the two iodine radicals become sufficiently separated that they move as separate atoms rather than experience bonding interactions. Under these long bond length circumstances the different total angular momentum orientation states on each radical are very close in energy and there will be a very different pattern of quasi-nonavoided crossings compared to those of the short bond length molecular states. These perturbed atom state crossings represent the random reorientation of the total angular momenta of the hole in each radical’s p-shell. The accurate treatment of these long bond length quasi-nonavoided crossings is probably not that critical for describing germinate recombination dynamics, as all the crossing surfaces at these long bond lengths give very similar forces. The accurate treatment of molecular state quasi-degeneracy at short bond lengths with the very different forces associated with different states will, however, be crucial to understanding the differences in early time nonadiabatic
molecular electronic relaxation in solids and liquids, and these types of processes will be the focus of this paper. Previous calculations from this group3,4 provided an apparently accurate description of excited state predissociation dynamics in liquids, but due to inadequate treatment of the quasinonavoided crossings present in the early time low-temperature solid phase dynamics, these earlier studies gave estimates for the solid phase lifetime that were more than an order of magnitude too fast. In this paper we show how accurate treatment of quasi-nonavoided crossings give lifetimes in both the solid and liquid phases that are in excellent agreement with experiments. On a more fundamental note we explore the findings of these calculations and their connection with the experiments in terms of the symmetries of the coupling between electronic angular momentum states and the local solvent pair structure found even in a liquid. The findings of these calculations are thus understood in terms of approximate “selection rules” for condensed phase nonadiabatic electronic relaxation that arise due to the interplay between electronic state symmetries and local solvent structure. The paper is organized as follows: In section 2 we briefly outline the standard surface hopping approach, which uses the instantaneous adiabatic electronic state representation. Next we describe how this approach must be adapted to treat condensed phase systems in which quasi-nonavoided crossings play a crucial role in electronically nonadiabatic relaxation. We also outline how we use these trajectories to compute experimental pump-probe signals. Section 3 details the results of applying these modified methods to explore nonadiabatic electronic relaxation of photoexcited I2 in liquid and solid rare gases. We show that the important features in the experimental results are reproduced by these calculations, and we discuss our finding that nonadiabatic relaxation in this model condensed phase system seems to take place through specific channels determined by interplay of electronic state symmetry and local environmental symmetry. We show why these approximate nonadiabatic “selection rules” are even active in the liquid. Concluding remarks are presented in section 4. 2. Methods 2.1. Nonadiabatic Molecular Dynamics. Despite the considerable recent advances in semiclassical5-11 and quasiclassical12-29 methods for treating excited state nonadiabatic dynamics, methods based on surface hopping ideas are by far the most practical approach currently available for incorporating the possibility of nonradiative electronic transitions into condensed phase molecular dynamics calculations. In this section we describe the problems that arise from quasi-nonavoided crossings and show how they can be solved within the framework of the surface hopping (SH) approach. Similar problems arising from state degeneracy are expected to plague any approach to nonadiabatic dynamics that uses the instantaneous adiabatic representation. We assume that the reader is familiar with the fundamental aspects of the SH method.12-17 In SH one solves the equations of motion of classical mechanics over potential energy surfaces that result from diagonalizing an electronic Hamiltonian at the current nuclear positions. Simultaneously, the time dependent Schrodinger equation is solved for the electronic degrees of freedom moving under the influence of the time varying field generated by the classically evolving nuclei. This quantum equation is propagated solely as an auxiliary equation with the purpose of deciding over which electronic potential energy surface the classical dynamics is supposed to be carried out. A
6730 J. Phys. Chem. B, Vol. 105, No. 28, 2001
Yu et al.
swarm of these quasi-classical trajectories is evolved and its dynamics is supposed to represent the evolution of the nuclear density under the noninteracting trajectory approximation.18 This nuclear density can branch at avoided crossing points and hence display the different output channels in a nonadiabatic event. The SH algorithm models transition probabilities to the different electronic states as first-order kinetic rates; hence gif the probability of changing the classical dynamics from the current surface “i” to new surface “f” in an electronic time step δ is computed using Tully’s fewest switches model of the process14 according to the following result
(
)
Ffi〈φf|∂/∂t|φi〉 gif ) -2δRe Θ[Re(Ffi〈φf|∂/∂t|φi〉)] Fii
(1)
In this expression Θ is the Heaviside step function, and φk are the instantaneous adiabatic eigenstates in terms of which we write the time dependent mixed state electronic wave function ψ(t), the solution of the quantum auxiliary equation discussed above, as ψ(t) ) ∑jaj(t)φj. The time dependent expansion coefficients here are used to define the electronic density matrix elements Ffi(t) ) a/f (t) ai(t) appearing in the above expression for the transition probability. These probabilities of transition are fed into a stochastic procedure that decides if a switch of the label of the occupied state will occur or not (see refs 14 and 15 for details). If a switch in the identity of the occupied state occurs during a nuclear time step composed of many (n) small electronic steps, i.e., ∆ ) nδ, the impulsive nuclear forces active as a result of the transition are computed using the following expression for the space and time localized Pechukas force.17,30
Ffi ) -Re
{
}
(Ef - Ei) D ˆ fi D ˆ fi‚R˙ ∆
(2)
In this expression D ˆ fi ) 〈φf|∇R|φi〉 are the nonadiabatic or derivative coupling matrix elements and only the directions of these vectors in the 3N dimensional coordinate space of the N particle system are important for determining the nuclear forces during nonadiabatic transitions. We propagate the expansion coefficients ak(t) ) 〈φk(t)|ψ(t)〉, which determine the density matrix elements used to compute the hopping probabilities given in eq 1, according to the following short time approximate propagation scheme
〈φj(t+δ)|ψ(t+δ)〉 )
∑i 〈φj(t+δ)|φi(t)〉e-i(E (t+δ)+E (t))δ/2p〈φi(t)|ψ(t)〉 j
i
(3)
and we employ the following estimator for the ∂/∂t matrix elements
〈φf|∂/∂t|φi〉 ) 〈φf|∇R|φi〉‚R˙ ) -
〈φf|∇RH|φi〉 ‚R˙ (Ef - Ei)
(4)
to compute this final ingredient of these hopping probabilities. The important feature of these results is that it is unnecessary to follow the basis set phase as must often be done when other estimators of these quantities are employed to compute hopping probabilities.31.32 In our semiempirical diatomics-in-molecules (DIM) electronic structure calculations,3,33,34 we write the electronic eigenfunctions in terms of a complete orthonormal set of polyatomic basis functions Φj, which are assumed to vary slowly
with nuclear configuration. Thus, the eigenstate basis functions are written as
φn )
∑j ΓnjΦj
(5)
and Γ are the DIM eigenvectors. Using the arguments presented in ref 3, we thus compute the various quantities needed for the calculations discussed above as
〈φf|∂/∂t|φi〉(t+kδ) ) Γ/fm(t+kδ) Γin(t+kδ)
∇R〈Φm|H|Φn〉(t)‚R˙ (t) ∑ ∑ m n (E (t+kδ) - E (t+kδ)) f
(6)
i
Here the ∂/∂t matrix elements at a point k electronic time steps within the current nuclear time step, which began at time t, are computed assuming both the nuclear velocities and Hamiltonian matrix elements in the polyatomic function basis set are slowly varying on the electronic time scale. The electronic coefficient propagation can be done using the rapidly varying DIM eigenvectors defined above by employing the following expression
〈φj(t+δ)|ψ(t+δ)〉 )
∑i ∑n Γ/jn(t+δ) Γin(t)e-i(E (t+δ)+E (t))δ/2p〈φi(t)|ψ(t)〉 j
i
(7)
Finally, we compute the nonadiabatic coupling vectors D ˆ fi ) 〈φf|∇R|φi〉, which are used in eq 2 by finite difference only when a hop between surfaces is accepted. 2.2. Quasi-Nonavoided Crossings. The nonadiabatic dynamics of molecules in fluctuating, partially ordered environments is complicated due to transient symmetries in the environment that can result in transient symmetries of the electronic states. Electronic state symmetries can be used to great advantage in isolated molecule calculations, as they mean that the Hamiltonian can often be block-diagonalized. Thus, rather than diagonalizing the full electronic Hamiltonian, we need only diagonalize each individual block of coupled states. Suppose, for example, that the symmetry of an N state electronic problem is such that there are two blocks of states that are uncoupled from one another, one containing L states and the other M states say. Thus, we need only diagonalize the separate L × L and M × M blocks rather than the full N × N matrix in order to obtain the eigenvalues and eigenvectors of the full system. For an excited molecule in a fluctuating, quasi-ordered environment like a crystal at finite temperature, for example, a higher symmetry situation like this may exist for a brief period and as energy from the excited molecule is deposited into the surroundings, or the molecule moves in an asymmetrical way relative to its environment, the symmetries may be washed out, leading to a situation of lower symmetry, perhaps more like that of a molecule in a disordered liquid environment. Technically, it is easier (though obviously more expensive) just to ignore these transient symmetries and blindly diagonalize the full electronic Hamiltonian each time step. Such an approach, however, will fail when used in conjunction with the surface hopping method described above unless special care is exercised. To see this, suppose the symmetry of the system in some configuration is such that the Hamiltonian can be broken into two blocks as described above. States within the same block are coupled and so they must exhibit avoided crossings. States from the different blocks, however, are completely uncoupled, so they can have quasi-degenerate energies and thus show nonavoided or quasi-nonavoided crossings. In this case the
Photodissociation Dynamics of I2
J. Phys. Chem. B, Vol. 105, No. 28, 2001 6731
energy denominators in eq 4 or 6 will vanish, making these estimators of the ∂/∂t matrix elements diverge, giving an illdefined result for the transition probability in eq 1. The hopping algorithm must thus be carefully designed to do what is physically realistic at these quasi-nonavoided crossing points. Obviously, if these states are truly from different instantaneous symmetry blocks, and thus with no couplings between them, then there should be no nonadiabatic mixing and hence no hopping between them. These divergent transition probabilities predicted by eq 1 in this situation are meaningless, and the hops implied by them should be ignored. We thus need an approach to identify when these quasi-nonavoided crossings are occurring and once identified the algorithm should continue with same symmetry state occupied rather than make a nonphysical hop to a state of different symmetry. Outside such quasi-nonavoided crossing regions, the regular surface hopping approach should yield transition probabilities that are reliable provided the approximations of surface hopping are satisfied. The method we employ to identify these quasi-nonavoided crossing situations involves monitoring the overlap of the eigenstate occupied at electronic time t + nδ with the same energy ordered state at t + (n + 1)δ. Thus, if during an electronic time step δ this overlap falls less than some tolerance, i.e.
〈φi(t+nδ)|φi(t+(n+1)δ)〉 < tol
(8)
then the occupied state must have crossed some other state (the new state i labeled according to energy) during δ. Such a crossing means that the two states are essentially uncoupled and no transitions between them should in reality occur. At this point we thus turn off the state hopping so as to avoid any spurious values of gif that result because this quantity cannot be reliably estimated under such quasi-nonavoided crossing circumstances and we search through the new states to find the one with the strongest overlap with the state occupied at the beginning of the nuclear time step, i.e., t. This strongest overlapping state is made the newly occupied state, thus preserving state symmetry in this situation where no transition should occur. In our calculations we varied the tolerance parameter over the range 0.01 < tol < 0.1 and obtained converged results. We thus used tol ) 0.05 in all calculations reported here. It is important to note that the prescription given in eq 11 does not rely on identifying symmetric environments, so it ought to work automatically even when an accidental degeneracy is encountered. Situations of accidental degeneracy will be generic in disordered environments and our method should in principle be able to handle them. An approach such as that outlined above, however, will encounter problems if the start of one of the electronic time steps coincides closely with the “point” of degeneracy. When the eigenvalues and eigenfunctions are determined numerically with finite precision, these crossings no longer occur at points; rather, there will be crossing regions corresponding to crossing time intervals. At short internuclear distances when the forces are large and the state energies vary rapidly, such problems will be considerably less important, and our algorithm should step over these quasi-nonavoided crossings, employ the overlap criterion to realize that state branching should not occur, and proceed on the same symmetry state. At long bond lengths, however, many surfaces converge and the couplings between states become small so our approach will begin to confuse the occupied state identity under these circumstances. In this paper, however, we will use the overlap detection scheme described above to explore the difference
between pre-dissociation dynamics of I2 in liquid and solid solvation environments. As we shall see in the next section, this predissociation dynamics involves motion through regions in which the states cross one another at sufficiently steep slopes so that the time-stepping algorithm rarely lands in a small degenerate region; rather, it steps over such small “points” which it successfully detects using the overlap criterion. We stop the runs once the bond extension is beyond Rmax ) 4.7 Å to avoid entering completely dissociated regions where the states can become quasi-degenerate over extended regions and detection becomes problematic for the current approach. As mentioned in the Introduction, in these dissociated regions the quasi-degenerate states actually involve reorientation of the total angular momentum of the separated radicals by their local solvent environment. We are currently developing techniques that can reliably treat the passage of the system through such extended regions of quasi-degeneracy, and the application of these methods to describe the longer time nonadiabatic dynamics that leads to germinate recombination will be the subject of a future publication.35 2.3. Semiempirical Electronic Structure. The final key ingredient of the dynamics calculations reported here is the electronic structure approach we employ to generate the excited state potential surfaces and couplings. We refer the reader to the series of papers3,33,34 in which this approach was developed for the details and also to the similar DIM approach of reference.36 Alternative Hamiltonians describing the solvent induced coupling between the B state and various other electronic states of iodine in condensed environments have been proposed.37,38 Apkarian and co-workers1 have used these various formulations to analyze their experimental results and suggest that our approach offers an accurate representation of these couplings. 2.4. Calculation of Pump-Probe Signals. In the next section we report pump-probe signals calculated from our ensemble of nonadiabatic MD trajectories. These signals we computed using the following expression4
σ(t) ∼
∫-∞∞dt′exp[-(t - t′)2/2δ2]S(t′)
(9)
assuming laser pulses whose autocorrelation function is a Gaussian in time with width δ ) 100 fs1, and
S(t) ) 1 N
∑ ∑f |µfjad(Rk(t))|2 exp[-{∆Vfj (Rk(t)) - hνprobe}2/2∆2]
N k)1
k
k
(10) where we sum over the ensemble of surface hopping trajectories k that occupy the jkth adiabatic electronic eigenstate at time t and over all possible final adiabatic electronic eigenstates f. Here ∆ is the energy spread of the probe laser pulse (we used 0.05 eV in these calculations1), νprobe is the probe laser frequency (we used νprobe ) 580 nm in this work), and ∆Vfjk ) Vf(Rk(t)) - Vjk(Rk(t)) is the energy gap between the adiabatic electronic states. Finally, in this expression the dipole matrix elements µad fjk between these states are assumed to be the same constant value for all pairs of states that are separated from the occupied state by about hνprobe. 3. Results and Discussion 3.1. Overview of Experimental Findings. The methods described in the previous section have been employed to explore
6732 J. Phys. Chem. B, Vol. 105, No. 28, 2001
Figure 2. Chart of experimental pump-probe signals from ref 1. The upper panel shows these results for excitation of I2 in liquid CCl4 with λpump ) λprobe ) 580 nm. The lower panel shows the experimental results for excitation of I2 in solid krypton with λpump ) λprobe ) 556 nm.
the effect of solvent environmental symmetry on excited state nonadiabatic relaxation dynamics. The system we have chosen to study, I2 in liquid and solid rare gases, is motivated by the availability of exceptionally detailed microscopic experimental information from many experimental studies that have been conducted on this system over the years39-46 that make it the paradigm for much of our understanding of molecular relaxation in condensed phases. In addition, the availability of detailed, apparently reliable, electronic Hamiltonian information, and the experience in developing flexible multistate potential models describing excited states of iodine-containing molecules in condensed phase environments3,33,34 make this problem the perfect testing ground for developing this understanding of the effect of environmental symmetry on nonadiabatic electronic relaxation. The experiments that originally inspired our calculations come from the laboratory of Apkarian and co-workers,1,47,48 and the recent work of Schwenter and his group2 is also closely related. Typical experimental pump-probe results, here reproduced from the work of Apkarian’s group,1 are presented in Figure 2. In the liquid (top panel) these results show at most two or three features associated with excited B state vibrations lasting for on the order of 1 ps before rapid nonadiabatic electronic relaxation depletes B state population and these vibrational features decay. In the solid (lower panel), however, the experiments indicate coherent B state vibrational features persisting for on the order of 10 or more vibrational periods out to times of more than 3. ps. Thus, the experiments reveal
Yu et al. that the electronic relaxation dynamics is nearly an order of magnitude slower in the solid than in the liquid. 3.2. Simulation Details. We performed simulations of the ensuing dynamics after B state photoexcitation of I2 in solid and liquid rare gases. To avoid finite size effects, we used a simulation box containing 500 atoms in these calculations. Our results have been averaged over only a small ensemble of 16 trajectories, all of which show qualitatively similar features so, though we do not expect these results to be completely converged especially at the longer times, the qualitative behavior of our results is expected to generally be reliable. The solid argon samples used in these studies were prepared from an initial fcc lattice, at a density F* ) 1.072 78, from which two nearest neighbor argon atoms were removed and the I2 was inserted into this defect site. This configuration was then equilibrated at 40 K for 10 ps, following which a set of independent initial configurations (each separated by at least 1 ps) with energy gaps resonant with the chosen pump laser wavelength (in these calculations we used λpump ) 533 nm for both the liquid and solid) were sampled and used as initial phase space points for propagation on the B state after photoexcitation.4,49 Our studies of liquid state photodissociation dynamics were performed in liquid xenon prepared in much the same way as described above, except that the initial density of the lattice from which the fluid samples were obtained was about 30% less (F* ) 0.7) and the samples were equilibrated at 300 K. Initial independent resonant conditions were generated following the same procedure. 3.3. Pump-Probe Signals and Nonadiabatic Trajectories. In Figure 3 we present the results of our calculations of the time dependent pump-probe signals in liquid xenon (upper panel) and solid argon (lower panel). Given the small ensemble size used in these calculations, the agreement with the experimental results presented in Figure 2 is quite remarkable. Comparing the experimental and calculated pump-probe signals in liquids (top panels in Figures 2 and 3, respectively), for example, we see that the frequency of the vibrational features, their relative amplitudes, and decay rate are very nearly quantitatively reproduced by our calculations. The differences in the experimental and calculated liquid results are probably real as they are consistent with differences in solvent LennardJones radii. The experimental results were measured for liquid CCl4 solvent, while our calculations are in liquid xenon at about the same solvent density. The larger void space around the I2 molecule in the larger Lennard-Jones radius xenon solvent gives rise to a slightly longer vibrational period in the calculated results than observed in the experiments in the smaller LennardJones radius CCl4 solvent, where the voids around the I2 are tighter. The agreement between the experimental and calculated pump-probe signals in the solid (bottom panels of Figures 2 and 3, respectively) is also qualitatively very good. The calculations reproduce very nearly to the same oscillation frequency and decay rate as observed in the experiments with 7 or 8 clear vibrational features present before the signals decay by about 2500 fs. Differences associated with the relative solvent radius are expected to be smaller here as the radii of krypton (experimental solvent) and argon (solvent used in our solid calculations) are more similar. The main differences observed between the experimental and calculated solid phase pump-probe signals are, however, associated with other real differences in the calculation and experimental setup. The experimental results presented here actually involve probe transitions between the occupied mo-
Photodissociation Dynamics of I2
J. Phys. Chem. B, Vol. 105, No. 28, 2001 6733
Figure 3. Pump-probe signals calculated from our nonadiabatic surface hopping trajectories. The upper panel shows these results for excitation of I2 in liquid xenon and the lower panel shows our results for I2 excitation in solid argon. In both cases with λpump ) 533 nm and λprobe ) 580 nm.
Figure 4. I2 bond length histories for nonadiabatic ensemble of surface hopping trajectories. Upper panel shows bond length histories after photoexcitation with λpump ) 533 nm in liquid xenon (F* ) 0.7, T ) 300 K). Lower panel shows bond length histories after photoexcitation with λpump ) 533 nm in solid argon (F* ) 1.07278, T ) 40 K).
lecular B state and states in the ion pair manifold that are accessed with a probe wavelength 556 nm. In our calculations, however, we use the same 580 nm probe wavelength as used in the liquid studies, which excites from the occupied B state to states in the I*I and I*I* spin-orbit excited state manifold. The probe absorption windows for the experimental B to ion pair state probing employed in the solid measurements occur at longer bond lengths around 4 Å. The 580 nm probe between the B state and the spin-orbit excited states used in the liquid experiments and both our solid and liquid calculations, however, involve probe windows near 3 Å. Consequently, the experimental solid signal shows delayed absorption to the ion pair states as the initially excited B state packet starts at short bond lengths near 2.7 Å and the packet must evolve before it enters the 556 nm probe window. The 580 nm probe signals, on the other hand, show immediate intensity since the initial packet pumped to the B state is born in the 580 nm probe window. In Figure 4 we compare the I-I bond lengths as functions of time for our ensemble of trajectories in liquid xenon (upper panel) and solid argon (lower panel). These trajectories are terminated when this bond length extends beyond 4.7 Å for the first time, as discussed at the end of section 2 since the states are weakly coupled and closely spaced in energy beyond this point and state identities cannot be reliably tracked using the approach employed here. From this figure it is clear that these calculations are capturing the significant differences in predissociation dynamics between solid and liquid solvation environments observed in the experiments.
In the liquid we see trajectories dissociating from the initially excited B state very rapidly. The basic behavior observed in these calculated liquid trajectories involves the initially excited packet starting from the steep repulsive wall of the B state surface. Nearly 40% of these liquid trajectories then undergo rapid dissociation on the very first bond extension in less than 100 fs. The liquid state trajectories that remain undissociated in their B states, however, reach the outer turning points of their first bond extension, exchange energy with the solvent environment, and then pass back to slightly longer bond length inner turning points. During this first bond extension, the ensemble undergoes significant dispersion due to the liquid disorder encountered in the first collision and the ensemble vibrations that rapidly become decoherent beyond this point. On the second and third bond extensions of the liquid ensemble trajectories, we see diminishing waves of dissociation taking place between 500 and 1000 fs and few trajectories make it to a fourth bond extension in the liquid. A typical plot of the state energies as functions of time for these rapidly dissociating trajectories in liquid xenon is shown in the upper panel of Figure 5. This trajectory dissociates through the dominant 3Π2g exit channel state (see Figure 1). This 3Π2g exit state accounts for nearly 70% of the dissociation on B state excitation in our liquid calculations. Other minor dissociation channels observed in our calculations include 3Πg, which accounts for about 20% of dissociation, and 1Πu, which is responsible for less than 10% of dissociation in these liquid state simulations.
6734 J. Phys. Chem. B, Vol. 105, No. 28, 2001
Yu et al.
Figure 6. Solid curve shows the decay of probability of finding a bound excited I2 molecule as a function of time in our liquid xenon calculations, and the dashed curve shows same results for our solid argon study.
Figure 5. Representative bond length and energy level histories for I2 photodissociation dynamics in liquid xenon (upper panel) and solid argon (lower panel). State points are the same as those shown in Figure 4. The lowest curve in each panel is the I2 bond length for the given trajectory in Ångstroms, while the other curves are the state energies in electronvolts. We have added 13 eV to these energy curves to get them on an appropriate scale. The instantaneously occupied state in each panel is the thick curve marked by symbols.
An important aspect to note in these typical energy level plots (see Figure 5) for both the liquid and solid trajectories is that the energy of the trajectories generated with this pump wavelength is typically sufficiently large that all the possible crossing dissociative states just mentioned, and in fact several others at longer bond lengths, are repeatedly crossed during the ensuing photoexcited dynamics. Thus, the fact that we see dissociation occurring preferentially through a particular electronic state, the 3Π2g state, is not due to this being the only state accessible; rather many different states can be reached by the dynamics at these energies, so it is apparent that some approximate selection rule must be operative in these calculations, preferentially selecting out the 3Π2g state as the dominant dissociation channel. The situation observed in our low-temperature solid argon calculations (see lower panels of Figures 4 and 5) is in marked contrast to that found in the liquid xenon calculations described above, and these differences are consistent with the experimental differences between the liquid and solid phases summarized in Figure 2. From the lower panel in Figure 4 we see that the bond length histories of the solid argon ensemble do not start to show significant dissociation dynamics until at least the third or fourth bond extension and many trajectories undergo on the order of 10 vibrational periods before dissociating, some continue to vibrate in the excited B state for as long as 3 ps after initial
excitation. Further, due to the highly structured local environment in the dense crystal, the vibrations of the trajectories in the ensemble remain coherent with one another for more than 1 ps, with relatively little dispersion of the ensemble out to this time. The state energy histories for typical solid phase trajectories (see for example the results presented in the lower panel of Figure 5) show that basically the same channels are responsible for this longer time predissociation dynamics as those that were active in liquid phase dissociation and that the relative importance of the different channels in the solid is similar to what was observed in the liquid, with the dominant dissociation channel again being the 3Π2g state in these studies. Another clear difference apparent when comparing the liquid and solid bond length histories presented in Figure 4 is that not only is there a much larger dispersion of the initial I2 bond extensions in the liquid xenon trajectories, but also the amplitudes of these extensions are generally significantly larger in the liquid. This results from the reasonably sizable differences in density in the liquid and solid state points considered in these calculations. The effective density difference is actually a little larger than 30% due to way we remove solvent atoms in order to make room for the I2 solute; in xenon we are removing larger solvent atoms, making a bigger cavity than in argon, and that cavity is being filled by the same size I2 molecule. The surprising aspect of this observation is that even though the solvent atoms are on average further from the I2 in the liquid, they manage to perturb the electronic states of the molecule significantly more than in the higher density solid, resulting in the much faster liquid state dissociation dynamics. Clearly, the couplings between I2 electronic states induced by interactions with the solvent that are responsible for the dissociation dynamics are not solely determined by the proximity of the solvent atoms; rather they must be strongly dependent on anisotropy of the local solvent geometry. In Figure 6 we show how the fraction of undissociated trajectories that were initially excited to the B state decays as a function of time in our calculations. In the liquid (solid curve) trajectories start to dissociate more or less immediately. In the solid (dashed curve), however, we see that dissociation is delayed for more than 1.5 ps. When the solid ensemble eventually does start to decay, however, it does so at roughly the same rate as we see for the liquid. This delayed dissociation in the solid relative to the liquid suggests that the initial local crystalline order prevents I2 dissociation by suppressing the nonadiabatic transitions out of the initially excited B state.
Photodissociation Dynamics of I2 However, once sufficient disorder in the local crystal environment around the excited solute molecule has been established by moving and reorienting the I2 molecule relative to its surrounding, the electronically excited molecules start to undergo nonadiabatic relaxation at roughly the same rate as they do in the liquid, resulting in the observed decay. 3.4. Electronic Couplings in Condensed Environments. Two related aspects of our observations; the fact that dissociation in the liquid is much faster than in the solid and that some dissociation channels are favored over others, can be understood by considering the anisotropy of the off-diagonal Hamiltonian matrix elements, calculated here with DIM, which couple the I2 gas phase states as a result of interactions with the solvent. As we have discussed in detail elsewhere,3,33,34 the diatomicsin-molecules approach we use here to compute the excited state electronic structure is designed specifically to treat the coupling between different electronic angular momentum states of the iodine radicals as they interact with the locally anisotropic solvent environment. Apkarian and co-workers1 first proposed that anisotropy in the coupling between electronic states correlated with the local solvent structure is actually what is responsible for the significant differences they observed comparing their solid and liquid pump-probe signals. As we shall see below, these anisotropies in electronic coupling, together with the local pair correlation between solvent atoms in the immediate vicinity of the I2 solute are responsible for determining what are essentially condensed phase nonadiabatic transition “selection rules” governing which states will be the dominant nonadiabatic relaxation pathways and so determine the importance of different dissociation channels in condensed phases. To develop an understanding of the fundamental geometrical factors that are responsible for these observations, we now analyze the off-diagonal couplings between the B state (the state in the (J1 ) 3/2, J2 ) 1/2) manifold excited in these experiments) and the states from the (3/2, 3/2) ground state manifold that cross the B state, which are responsible for the different possible nonadiabatic electronic relaxation channels out of this initially excited state. In the gas phase these couplings vanish in our model so here we study their behavior in the I2-single rare gas atom cluster system. In Figure 7 we explore the dependence of these coupling matrix elements as functions of the angle between the I2 bond axis and the vector pointing from the I2 center-ofmass to the solvent atom. Within the DIM approach, just as the couplings displayed in this figure come from adding anisotropic contributions from the interactions of each iodine center with the solvent atom, the off-diagonal elements of the electronic Hamiltonian of the full many-body condensed phase system come from adding signed contributions (these values are generally complex in the angular momentum basis set employed here) just like those displayed in these figures for each of the surrounding atoms in the solvation environment. The general nonadditivity of the final adiabatic potentials occurs as these surfaces are obtained by diagonalizing the full Hamiltonian matrix whose elements are calculated in this approach by simply summing up the contributions like those displayed in this figure from all the solvent atoms in the environment around the I2 molecule in the instantaneous configuration of the system. The first thing to note in all these figures is that solvent atoms above and below the I2 bond, i.e., in the regions around 90° and 270°, make little contribution to the couplings between the B state and any of these different possible dissociation channel states. We see that atoms near the ends of the molecule around 0° and 180°, however, make the most important contributions to the nonadiabatic couplings. The next important thing to realize
J. Phys. Chem. B, Vol. 105, No. 28, 2001 6735
Figure 7. Plots of off-diagonal DIM Hamiltonian matrix elements for the I2-argon atom cluster system as functions of angle between the I2 bond vector and the center-of-mass to argon vector with distances fixed at RI-I ) 3.0 Å, and Rcom-Ar ) 4.2 Å. Solid curves are real parts, and dashed curves are imaginary parts of the matrix elements. The upper panel gives the coupling between the B state and the 1Πu state, the middle panel gives the coupling between the B state and the 3Π2g state, and the bottom panel gives the coupling between the B state and the 3Π state. g
is that the fundamental structural motif in an atomic solid or liquid solvation environment is a pair of atoms. From the top and bottom panels in Figure 7 it is clear that the coupling between the B and 1Πu or the B and 3Πg states, respectively, are antisymmetric functions of angle. Thus, a pair of solvent atoms near one end of the molecule will, on average, give contributions to these couplings which cancel, one atom giving a positive and the other a negative contribution, making electronic relaxation through these channels less likely. From
6736 J. Phys. Chem. B, Vol. 105, No. 28, 2001 the middle panel in Figure 7, however, it is clear that the coupling between the B state and 3Π2g is a symmetric function of angle; thus a pair of atoms near the end of a molecule will give contributions to the coupling with this state, which add up constructively, and unless they are matched up precisely with another pair of atoms at the other end of the molecule, there will be a finite nonadiabatic coupling between the B state and this particular state, which will make the 3Π2g the most likely nonadiabatic relaxation channel out of the B state. Thus, even though the magnitudes of the couplings to the 1Πu and 3Πg are somewhat larger than for the 3Π2g, these local symmetries are responsible for making this state the dominant dissociation channel after B state excitation, as we observe in our calculations. This observation is consistent with the experimental findings,1 which suggest that dissociation takes place through states that cross the B state at bond lengths in the range 3.05 Å < Rc < 3.8 Å, so the experiments suggest that the 1Πu state is an unlikely dissociation channel since this state crosses the B state at shorter bond lengths outside this range. At this point the experiments cannot distinguish between the other possible dissociation channels. From the discussion presented above, the origin of the significant differences observed in dissociation dynamics between the solid and liquid in both our calculations and the experiments is now clear: In the highly ordered environment of the solid, contributions to the coupling between the B and 3Π states from pairs of atoms at one end of the molecule are 2g canceled by pairs of atoms at the other end, giving a fragile balance of these electronic interference effects. The iodine atoms must collide with these end atomic pairs several times to disrupt the local crystalline solvent symmetry between the two ends of the solute molecule. Once this symmetry is broken, interactions with the nearest pair to one end of the molecule will give the finite coupling between the B and 3Π2g states and dissociation can begin to take place. In the liquid, on the other hand, the local disorder will guarantee the nonequivalence of the instantaneous solvation environments at either end of the molecule right from the moment of initial excitation resulting in immediate decay of the signal. 4. Conclusions This paper highlights the important connection between dynamics, local environmental symmetry, and electronic structure in determining nonadiabatic relaxation. The approach we detail in section 2 for recognizing quasi-degeneracy and correctly dealing with these situations when they occur in condensed phase applications is crucial for obtaining results that very accurately reproduce the trends observed in the experiments. With this approach we are able to reliably calculate detailed pump-probe signals that are sensitive to nonadiabatic relaxation and that show strong solvation environment effects. A further challenge is the development of an approach for treating extended regions of quasi-degeneracy such as those encountered when the radicals resulting from the dissociation processes considered here reach larger separations in solution. Accurate descriptions of the electronic dynamics through such regions will be crucial for treating germinate recombination dynamics and understanding how the solvent influences the pathways back to the molecular system. Our finding that the excited state nonadiabatic relaxation dynamics which leads to molecular dissociation in solution apparently exhibits very specific pathways governed by molecular electronic and local environmental symmetries is intrigu-
Yu et al. ing, and the suggestion above that the recombination dynamics is governed by similar geometrical restrictions is worthy of further exploration. It might be very interesting to explore these effects, on dissociation in particular, in studies of excited state photodissociation dynamics in cluster environments in which the local solvent asymmetry can be precisely controlled by varying the cluster size. We expect that small clusters with one or two solvent atoms may show asymmetrical structures, while larger clusters with complete solvation shells offer very symmetric environments. If the nonadiabatic selection rules discussed in this paper are active, the predissociation dynamics in these different cluster environments will occur via different channels that should be easily identified in calculations, but observing the differences experimentally will be challenging, but potentially very rewarding. Acknowledgment. We wish to thank V. A. Apkarian, M. Bargheer, and J. Faeder for insightful discussions and gratefully acknowledge the financial support of the National Science Foundation (Grant No. CHE-9978320) and the Petroleum Research Fund administered by the American Chemical Society (Grant No. 34927-AC-6), and a generous allocation of supercomputer time from the Boston University Center for Computational Science. D.F.C. gratefully acknowledges the support of a Schlumberger Visiting Professorship in the Theoretical Sector of the Chemistry Department at Cambridge University. References and Notes (1) Zadoyan, R.; Stirling, M.; Ovchinnikov, M.; Apkarian, V. A. J. Chem. Phys. 1997, 107, 8446. (2) Bargheer, M.; Dietrich, P.; Donovan, G.; Schwentner, N. J. Chem. Phys. 1999, 111, 8556. (3) Batista, V. S.; Coker, D. F. J. Chem. Phys. 1996, 105, 4033. (4) Batista, V. S.; Coker, D. F. J. Chem. Phys. 1997, 106, 6923. (5) Sun, X.; Miller, W. H. J. Chem. Phys. 1997, 106, 6346. (6) Stock, G.; Thoss, M. Phys. ReV. Lett. 1997, 78, 578. (7) Thoss, M.; Stock, G. Phys. ReV. A. (8) Stock, G.; Muller, U. J. Chem. Phys. 1988, 108, 7516. (9) Stock, G.; Thoss, M. Phys. ReV. A 1999, 59, 64. (10) Bonella, S.; Coker, D. F. J. Chem. Phys. 2001, 114, 7778. (11) Bonella, S.; Coker, D. F. Chem. Phys. 2001, 268, 323. (12) Tully, J. C. In Classical and quantum dynamics in condensed phase simulations; Ciccotti, G., Berne, B., Coker, D., Eds.; World Scientific: Dordrecht, 1998; p 489. (13) Tully, J. C. In Dynamics on Molecular Collisions, Part B; Miller, W. H., Ed.; Plenum: New York, 1976; p 217. (14) Tully, J. C. J. Chem. Phys. 1990, 93 (2), 1061-1071. (15) Coker, D. F. In Computer Simulation in Chemical Physics; Allen, M. P., Tildesley, D. J., Eds.; Kluwer Academic Publishers: Dordrecht, The Netherlands, 1993; p 315-377. (16) Hammes-Schiffer, S.; Tully, J. C. J. Chem. Phys. 1994, 101, 4657. (17) Coker, D. F.; Xiao, L. J. Chem. Phys. 1995, 102, 496. (18) Burant, J. C.; Tully, J. C. J. Chem. Phys. 2000, 112, 6097. (19) Kohen, D. Stillinger, F. H.; Tully, J. C. J. Chem. Phys. 1998, 109, 4713. (20) Sholl, D. S.; Tully, J. C. J. Chem. Phys. 1998, 109, 7702. (21) Ben-Nun, M.; Martinez, T. J. J. Chem. Phys. 1998, 108, 7244. (22) Martinez, T. J.; Ben-Nun, M.; Ashkenazi, G. J. Chem. Phys. 1996, 104, 2847. (23) Prezhdo, O. V.; Rossky, P. J. J. Chem. Phys. 1997, 107, 825. (24) Topaler, M. S.; Allison, T. C.; Schwenke, D. W.; Truhlar, D. G. J. Chem. Phys. 1998, 109, 3321. (25) Bittner, E.; Rossky, P. J. Chem. Phys. 1997, 107, 8611. (26) Martens, C. C.; Fang, J.-Y. J. Chem. Phys. 1997, 106, 4918. (27) Prezhdo, O. V.; Kisil, V. V. Phys. ReV. A 1997, 56, 162. (28) Kapral, R.; Ciccotti, G. J. Chem. Phys. 1999, 110, 8919. (29) Nielsen, S.; Kapral, R.; Ciccotti, G. J. Chem. Phys. 2000, 112, 6543. (30) Pechukas, P. Phys. ReV. 1969, 181, 174. (31) Mei, H. S.; Coker, D. F. J. Chem. Phys. 1996, 104, 4755. (32) Krylov, A. I.; Gerber, R. B.; Coalson, R. D. J. Chem. Phys. 1996, 105, 4626. (33) Margulis, C. J.; Coker, D. F. J. Chem. Phys. 2000, 113, 6113. (34) Margulis, C. J.; Coker, D. F. J. Chem. Phys. 2001, 114, 6744.
Photodissociation Dynamics of I2 (35) Yu, N.; Coker, D. F. Nonadiabatic germinate recombination dynamics in condensed phases. In preparation for submission to J. Chem. Phys. 2001. (36) Buchachenko, A. A.; Stepanov, N. F. J. Chem. Phys. 1996, 104, 9913. (37) Ben-Num, M.; Levine, R. D.; Fleming, G. R. J. Chem. Phys. 1996, 105. (38) Roncero, O.; Halberstadt, N.; Beswick, J. A. J. Chem. Phys. 1996, 104, 7554. (39) Chuang, T. J.; Hoffman, G. W.; Eisenthal, K. B. Chem. Phys. Lett. 1974, 25, 201. (40) Harris, A. L.; Brown, J. K.; Harris, C. B. Annu. ReV. Chem. 1988, 39, 341. (41) Smith, D. E.; Harris, C. B. J. Chem. Phys. 1987, 87, 2709.
J. Phys. Chem. B, Vol. 105, No. 28, 2001 6737 (42) Kang, T. J.; Yu, J.; Berg, M. J. Chem. Phys. 1993, 94, 2413. (43) Scherer, N. F.; Jonas, D. M.; Fleming, G. R. J. Chem. Phys. 1993, 99, 153. (44) Zadoyan, R.; Sterling, M.; Apkarian, V. A. J. Chem. Soc., Faraday Trans. 1996, 92 (11), 1821. (45) Lienau, C.; Williamson, J. C.; Zewail, A. H. Chem. Phys. Lett. 1993, 213, 289. (46) Liu, Q.; Wan, C.; Zewail, A. H. J. Chem. Phys. 1996, 105, 5294. (47) Zadoyan, R.; Li, Z.; Martens, C. C.; Apkarian, V. A. J. Chem. Phys. 1994, 101, 6648. (48) Li, Z.; Zadoyan, R.; Apkarian, V. A.; Martens, C. C. J. Phys. Chem. 1995, 99, 7453. (49) Stirling, M.; Zadoyan, R.; Apkarian, V. A. J. Chem. Phys. 1996, 104, 6497.
