"~': "
'
r'J
21 July 1995
CHEMICAL PHYSICS LETTERS ELSEVIER
Chemical PhysicsLetters 241 (1995) 215-222
A multireference time-dependent coupled cluster study of the intramolecular vibrational relaxation process G. Sree Latha, M. Durga Prasad School of Chemistry, University of Hyderabad, Hyderabad500 134, India
Received 8 February 1994; in final form 26 January 1995
Abstract
A multireference time-dependent coupled cluster (MRTDCCM) study of the intramolecular vibrational relaxation process in a model hydrocarbon is presented. Two Ansatze, the ordinary exp(S) and the normally ordered exp(S) are used to map the model space component to the exact wavefunction. The survival probability of the initial state and the energy in the CH mode are calculated. It is found that the ordinary MRTDCCM provides a better description than the normally ordered approach for the short-time dynamics. 1. Introduction
The dynamics of a localised vibrational excitation in a molecule when collisional and radiative interaction with the environment are negligible can be defined as intramolecular vibrational relaxation (IVR) [1]. It provides insight into many experimental studies on unimolecular reactions and spectroscopic observations. Classical and quantum mechanical perspectives of such dynamics have attracted the attention of several authors [1-15]. An important conclusion that has emerged from these studies is that IVR is essentially a consequence of the presence of overlapping non-linear resonances between several degrees of freedom. Various theories of IVR developed up to 1990 have been reviewed by Uzer and Miller [2]. The quantum dynamical study of IVR requires numerical solution of the time-dependent SchrSdinger equation (TDSE), iO~b/Ot = H~b.
(1)
The brute-force method of solving the TDSE is to expand the time-dependent wavefunction as a linear superposition of the basis vectors of an appropriate Hilbert space and integrating the corresponding Hamiltonian matrix [4,8,9]. The computational resources required to solve such a configuration interaction (CI) approach increase exponentially with the number of degrees of freedom in the system. Although significant progress has been made recently in handling large basis sets for such calculations [4], this basis set bottleneck provides the motivation for developing alternative methods that can provide accurate approximations for the description of IVR. In this context, several authors have discussed the utility of dynamical basis sets for carrying out such quantum dynamical calculations based either on the semi-classical Gaussian wave packet propagation (GWP) technique [16-19] or the time-dependent self-consistent field (TDSCF) approximation [20-24]. More recently Topaler and Makri [5] developed a procedure based on a path integral formulation for the propagator when the coupling between the active
0009-2614/95/$09.50 © 1995 Elsevier Science B.V. All rights reserved SSDI 0009-2614(95)00631-1
216
G. Sree Latha, M. Durga P r a s a d / Chemical Physics Letters 241 (1990) 2 1 5 - 2 2 2
(initially excited) mode and the bath modes is linear in the bath coordinates. In this Letter, we present the results of the application of the multireference time-dependent coupled-cluster method (MRTDCCM) [25-28] for the description of IVR. Our goal is to see whether the MRTDCCM approach can provide an accurate description of IVR at a low-order approximation. In the MRTDCCM approach [28] the time evolution operator is split into two parts. The first part describes the evolution of the projection of the wave packet in a model space consisting of states that are strongly coupled to the initial state. The second part which is analogous to the Moiler wave operator in stationary state theories [29], maps the model space component to the full wavefunction. The dynamics within the model space are treated exactly and the approximations are confined to the construction of the wave operator which accounts for the weaker interaction between the model space and the complementary virtual space. An exponential Ansatz is posited for the wave operator. As a consequence, the generators of the wave operator are additively separable. This allows one to decompose the dynamics in terms of the subsystem dynamics as in the TDSCF approach. As is well known [21] TDSCF breaks down when the two- or higher-body correlations become significant. The MRTDCCM is geared to incorporate such many-body effects easily. The essential details of TDCCM are reviewed in Section 2 and its application to IVR in a model hydrocarbon chain is presented in Section 3.
2. Theory 2.1. T h e s y s t e m
In this section we introduce our system and review its features. We study the IVR in a simple hydrocarbon chain model introduced earlier by Hutchinson et al. [8,15] which consists of 11 carbon atoms with a hydrogen atom attached to one of the terminal carbons. The CC bonds are taken to be harmonic and the CH oscillator is taken to be a Morse oscillator. The initial conditions are defined such that the CH oscillator is in one of its excited eigenstates at t - - 0 . Hutchinson et al. found that
only four of the highest frequency normal modes of the CC chain receive energy in the dynamical flow due to non-linear resonances present in the system [8]. So only these four modes are retained. In terms of these coordinates, the Hamiltonian is given by H = H ° + V,
(2a)
H ° = p2H/2/XCH + O[1 - exp( - otqcn)] 2 4
+½ E (p~ + O:q~).
(2b)
j=l 4
V = - 1 / m c ~_, L l j l p c H p j .
(2C)
j=l
The relevant parameters were taken from Ref. [8]. We have chosen to study this system whose features are well known [8] because it is small enough to obtain converged basis set results for comparison. This enables us to judge the validity of the approximations derived from the MRTDCCM approach. We next specify the basis set. The eigenfunctions of H ° are the oscillator product states I m)l nl)l n2)l n3)l n4), where m is the quantum number of the Morse oscillator and n i is the quantum number of the ith harmonic oscillator. We will present the results for the CH vibrational quantum numbers m = 4 and 5 at t = 0 . We have taken M = max(m) = 4 and 5, for these two initial conditions. Hutchinson et al. [8] found that each chain mode receives only a small amount of energy, so we take N = m a x ( n / ) = 1 for each normal mode harmonic oscillator. The resultant basis sets are N t = 80 for the dynamics of the 14, 0, 0, 0, 0) state and N t = 96 for the 15, 0, 0, 0, 0) state. The MRTDCCM is a Fock space theory. It is thus convenient to rewrite the Hamiltonian in the second quantized form in terms of the single particle basis functions used to construct the Hamiltonian. This gives n°=
E
'~a %' : ana '~ , E%
(2d)
ct,n a
ot,[~,m a m,o,na,ng
× a'~+a a+a'~ a t3
--ra a-- ra0 -- n a~nl3 •
(2e)
G. Sree Latha, M. Durga Prasad/ ChemicalPhysics Letters 241 (1995) 215-222 J ~
217
is taken to be an element of a small subset of strongly interacting states {4'0} of the full Hilbert space. This subspace is called the model space. The initial state evolves in time under the influence of the time evolution operator U that satisfies the TDSE,
~ J
~b= U4'0,
(3a)
iU = HU.
(3b)
Within the framework of the time-dependent effective Hamiltonian theory [28,30], U is partitioned Fig. 1. (a) Diagrammatic representation of the Hamiltonian. (b) diagrammatic representation of the S2 and S3 operators. Double arrows represent the Morse states of the CH oscillator and the single arrows represent the particle and hole lines of the harmonic bath modes.
Here a, fl denote the vibrational modes and m, n denote the basis functions of these modes, a n (a~) are the usual creation and annihilation operators of the nth basis function of the ath mode. Fig. la depicts the terms diagrammatically. Two features of the Hamiltonian are worth noting. First, since the perturbation term contains only momentum coupling operators, the TDSCF would predict only trivial dynamics for this system in which only the overall phase of the wavefunction would change. Second, as noted by Hutchinson et al., several near degeneracies exist among the zeroth-order states of the system due to the CH bond anharmonicity. For example, 15, 0, 0, 0, 0) is nearly degenerate with the 14, 1, 0, 0, 0) and t4, 0, 1, 0, 0) states. The occurrence of these degeneracies causes irreversible IVR from the initially prepared states. These states are coupled to the initial states by two-body operators. Hutchinson et al. have shown that these neardegenerate states can be put into different tiers and only states in the neighbouring tiers (which differ in the occupation numbers of two modes) are coupled by the perturbation (2c). Thus an approximation to the time evolution operator containing only two-body operators should provide a good description.
as
U = UexUM.
(4)
where UM is the model space evolution operator that describes the dynamics of the component of the wavefunction in the model space, ¢ = Pq, = uM 4,0.
(5)
tg +
2.2. The M R T D C C M Ansatz
We now turn to the parametrization of the timedependent wavefunction in the MRTDCCM framework. In this approach, the wave packet 4,0 at t = 0
Here P is the projection operator onto the model space. A time-dependent effective Hamiltonian Hal is postulated that generates UM, i 0 M = HeffgM•
(6)
Note that Haf and UM are defined within the model space only. Substituting Eq. (4) into (3b), premultiplying by Ud ~ and comparing the P P component of the resulting equation with Eq. (6) we find neff = P(U~IHU¢x -- iU~x 1/)ex) P.
(7)
The condition that 4, evolves within the model space is satisfied if and only if [30] Q(UexmHUex - i U ~ ~/.)ex)e = 0.
(8)
Here Q = 1 - P is the projection operator onto the virtual space. The time-dependent coupled-cluster method [28] posits an exponential Ansatz for the time-dependent analog of the Moiler wave operator [29] that maps the model space projection 4, onto the full wavefunction O, ~b= U~x4,,
(9a)
Uex = exp(S).
(9b)
Here the cluster operator S induces transitions from the model space to the virtual space. An important property of the cluster operator is that asymptotically
218
G. Sree Latha, M. Durga Prasad / Chemical Physics Letters 241 (1990) 215-222
it is an additively separable quantity. If the system consists of two non-interacting fragments, the total cluster operator becomes the sum of the cluster operators of the fragments. This follows from the exponential Ansatz (9b) and the additively separable property of the Hamiltonian that generates the time evolution operator. As a consequence, the cluster operator can be parameterized in terms of the various subsystem cluster amplitudes. It is then possible to define a sequence of approximations in which subsystem cluster amplitudes corresponding to larger subsystems are included in a systematic manner. We note in passing that our Ansatz (9b) for U~x is somewhat different from the choice of Guha and Mukherjee [28] who prefer a normally ordered exponential Ansatz for it, U~x = N [exp(S)].
(9c)
The number of terms in Eq. (7) and (8) are fewer if a normally ordered Ansatz is used for U~x. However, we found that the ordinary exponential Ansatz is more accurate than the normally ordered Ansatz at least for short time dynamics as discussed below. Throughout our calculations we have taken the set of states {I m, 0, 0, 0, 0), 0 ~
-7- 7-7" 7"
approach or even a single-reference TDCCM approach would need all the n-body excitation operators to account for the influence of quasi-degenerate states. Second, this model space provides a one-toone correspondence with the Morse oscillator states. Thus it can be used to construct an effective Hamiltonian for the active mode alone. Third, this choice provides a stringent test for the applicability of the effective Hamiltonian approach to the quantum dynamics when the wave packet makes significant excursions outside the model space. With this choice of model space, the cluster operator S consists of operators that excite the bath modes and simultaneously cause scattering among the Morse states, S = S 2 + S3 + . . . .
(10a)
$2 = ~.,(mpi IS21 nhi)a 1m+apianah,, i+ i i
(10b)
i,Pi
etc. Here m, n indicate the Morse states, Pi and h i indicate the virtual and occupied states of the ith harmonic oscillator. Some of these are depicted diagrammatically in Fig. lb where standard diagrammatic notation is used [33,34]. Thus, the active states which are partially occupied in the model space are indicated by double arrowheads, particle states that are always unoccupied in the model space functions are marked by a single arrowhead pointing to the left and the remaining hole states are shown to be moving to the right. Because of the subsystem embedding condition [26,28,31], the one-body excitation operators that cause excitations of the bath modes alone are identically zero. Derivation of the working equations for the cluster matrix elements (Eq. (8)) is conveniently carried out in diagrammatic notation. The diagrams for iS 2 are presented in Fig. 2. Note that when Ansatz (9b) is used, the wave packet receives contributions from higher rank excitations even if the cluster operator is truncated at a low order. To see this, consider the approximation S = S 2. The wave operator is given by U~x = exp(S) =
l + ~ ( m p i l S s [n h i ) a
1m+ anap~ah~ 1 i+ i
+ ~ ~ (mpil S2I khi)(kpjl $21 nhj) k Fig. 2. T h e d i a g r a m m a t i c representation o f the equation o f iS~. The notation is the s a m e as in Fig. 1.
X
1+ 1 i + i
J+
J +
am a n a p i a h i a p i a h i
•.. ,
(11)
G. Sree Latha, M. Durga Prasad / Chemical Physics Letters 241 (1995) 215-222
Acting on the function r n h 1 h 2 • • • ) the third term in this operator generates a three-body excitation in which two of the harmonic modes are excited. In this manner all possible n-body excitations are generated. This is in contrast to the Ansatz (9c) which does not generate such higher many-body excitations due to the normal ordering imposed on it which excludes the possibility of S - S contractions [34]. The influence of these higher body excitations appears in the working equations through diagrams in which explicit S - S contractions appear. The last diagram in Fig. 2 is an example of this category. Note that a three-body excited state is an intermediate state with respect to the model space function on which this diagram acts. This is an example of active state propagator renormalization because it is generated by contracting an effective one-body operator with the rightmost S-operator. The other three diagrams in which S - S contractions appear are examples of vertex renormalization in that they involve effective two- and three-body operators. The effect of many-body excitations beyond the two-body level are accounted for in the MRTDCCM through these terms in the working equations.
a 0.8
0.6 .o 0
0 °4
L D.
O.2
~
3o.o
~o.o
Time
We first present the results of the dynamical calculations for survival probability of the initial state. The survival probability when the initial wave packet is 14, 0, 0, 0, 0) is presented in Fig. 3a. The MRTDCCM calculations were carried out under the S = S 2 and S = S 2 + S 3 approximations. A second calculation was carded out with the normally ordered exponential wave operator of Ref. [28] to see whether the Ansatz (9b) is able to simulate the effect of the many-body excitations correctly. The working equations for this calculation are obtained by dropping the terms containing S - S contractions in the ordinary exp(S) theory. For comparison, we have also carried out basis set exact calculation by expanding the wave packet as
~1= E Cmnln2n3n41mnln2n3n4)
(12)
and integrating the SchriSdinger equation for the C coefficients. The summations over m, n, ... were truncated as mentioned above. As can be seen from
70.0 (fs)
go.o
0.80
f"i /
0.60
°
°o i,.
0.20
.
\ 30.0
3. Results and discussion
219
/I
J/
' "-./2 J
60.0 Time (fs)
90.0
Fig. 3. The survival probability of the initial state when (a) 14, O, O, O, O) is the initial state, and (b) 15, O, O, O, O) is the initial state. The continuous line represents the basis set exact results, the dashed line the normally ordered MRTDCCM calculation at the S 2 level, the dash and dotted line the ordinary MRTDCCM at the S2 level and the dash and double dotted line the MRTDCCM at the S3 level of approximation.
Fig. 3a the MRTDCCM at S = S 2 level provides a good description of the irreversible decay of the initial state, and the inclusion of the three-body operator improves the result. In contrast, the normally ordered Ansatz shows artificial recurrences as early as four vibrational periods. As is well known, the use of the normally ordered exp(S) Ansatz with the two-body approximation is equivalent to carrying out dynamical calculations in the limited basis consisting of the model space and all the singly excited states generated from it, for the one-valence problem
220
G. Sree Latha, M. Durga Prasad / Chemical Physics Letters 241 (1990) 215-222
[32]. This limited basis set is the origin of these recurrences. In ordinary exp(S) theory, on the other hand, two, three and four-body excitations also influence the dynamics due to S - S contractions, though the coefficients of such configurations are parameterized as products of the lower rank excitation operators. As a consequence, it does not show any artificial recurrences. The survival probability corresponding to the initial state 15, 0, 0, 0, 0) is presented in Fig. 3b. The trends in this case are similar to those in Fig. 3a. However, after about five vibrational periods the survival probability starts increasing when the ordinary exp(S) Ansatz is used. We traced this artificial rise to the violation of the norm in the MRTDCCM approach. The MRTDCCM uses a similarity transformation to map the model space component of the wavefunction to the exact wavefunction. As noted in Ref. [30], the similarity transformation based methods are prone to norm violations due to intruder states that develop complex eigenvalues when the approximations invoked violate the Hermiticity of the Hamiltonian. Briefly, the function ~b is strictly confined to the model space only when Eq. (8) is exactly satisfied. In an approximate calculation Eq. (8) is not satisfied exactly. However, Heff is still calculated using Eq. (7) assuming Eq. (8) is satisfied. This is equivalent to propagating the wave packet with a modified Hamiltonian, HA=H-Ue~RUex
1
(13)
where R = Q(UZxlHU~x - iU~x1/.))P.
(14)
The second term in Eq. (13) is, in general, nonHermitian. Consequently, H A can develop complex eigenvalues at some stage of the calculation. When the wave packet develops a significant overlap with the eigenvectors corresponding to such complex eigenvalues, its norm grows beyond one. The energy in the Morse oscillator with the initial state 14, 0, 0, 0, 0) is plotted in Fig. 4a. As may be seen the trends are similar to Fig. 3 in that the ordinary exp(S) theory performs better than the normally ordered exp(S) theory and the inclusion of the three-body operator improves the agreement with the exact result in the correct direction extending the range over which the ordinary MRTDCCM is accu-
a 9.00
~D
L
8 .O0.
7.Of
6 .Off
b
!
I
t0.0
20.0 30.0 40.0 Time ( f s )
I
I
I
50.0
11.5
_ 10f x,, ,
'
20.0 Time
/ i
'
40.0 (fs)
5~1.0
Fig. 4. The energy of the CH mode (a) with the initial state 14, 0, 0, 0, 0) and (b) with the initial state 15, 0, 0, 0, 0). The notations is the same as in Fig. 3.
rate. However, after about four vibrational periods the energy content in the Morse oscillator starts increasing, again due to norm violations. Similar trends are seen for the initial 15, 0, 0, 0, 0) state also (Fig. 4b). In summary, we have applied the multireference generalization of the time-dependent coupled-cluster method to study the dynamics of the intramolecular energy transfer process. For the initial conditions that we have studied, only a small part of the wave packet resides within the model space by about three to four vibrational periods and there are no significant recurrences up to about ten vibrational periods. This rapid and irreversible decay is due to the presence of several overlapping non-linear resonances. All these resonances involve only two modes at a time i.e. the CH mode and one of the bath modes. With respect to the model space, however, they
G. Sree Latha, M. Durga Prasad/ Chemical Physics Letters 241 (1995) 215-222
appear as multi-mode resonances. Thus a CI or single-reference TDCCM-based approach requires the inclusion of multi-mode excitation operators in Uex for converged results. The MRTDCCM approach with ordinary exp(S) Ansatz is able to provide an adequate description of IVR at the two-body level even in this case since it simulates the three-, four-,.., body excitations as a sum of products of the S 2 operator. The major advantage of the method lies here. In general, the number of two-body operators increases only quadratically with the number of degrees of freedom while the overall basis set grows exponentially. Thus the MRTDCCM can be expected to be computationally much less intensive than a corresponding basis set calculation without sacrificing accuracy. Even if the three-body excitation operators are to be included the computational effort would scale only cubically with the number of degrees of freedom. However, one must choose the model space with some care as our calculations show, since if the wave packet makes large excursions outside the model space its performance may be marred by intruder states with complex eigenvalues. This may not be a serious problem, however. In a recent calculation Wyatt and co-workers [4] noted that the dynamics of CH chromophore in benzene are dominated by about four modes. Similar trends were noted by Iung and Leforestier in CD3H. If the model space is chosen to include all possible occupancies of these strongly interacting modes, an MRTDCCM approach should provide an adequate description at the S = S 2 level. While computations of this type are needed before a final judgement is passed on the MRTDCCM approach, we believe that our results indicate that it is a promising approach. The normally ordered exp(S) Ansatz, on the other hand, does not suffer from the intruder state problem. However, since it is equivalent to a limited basis set theory, it is clearly not converged in the basis and would require the inclusion of higher rank operators just as in CI.
Acknowledgement Financial support from DST is gratefully acknowledged, GSL acknowledges a sustaining fellowship from UGC.
221
References [1] M. Quack, in: Energy storage and redistribution in molecules, ed. J. Hinze (Plenum Press, New York, 1983) p. 493. [2] T. Uzer and W.H. Miller, Phys. Rept. 199 (1991) 75. [3] R. Marquardt and M. Quack, J. Chem. Phys. 95 (1991) 4854. [4] R.E. Wyatt, C. lung and C. Leforestier, J. Chem. Phys. 97 (1992) 3458, 3477; R.E. Wyatt and C. lung, J. Chem. Phys. 98 (1993) 3577, 5191, 6758. [5] M. Topaler and N. Makri, J. Chem. Phys. 97 (1992) 9001. [6] A.A. Struchebrukhov and R.A. Marcus, J. Chem. Phys. 98 (1993) 8443. [7] E.L. Sibert III, J.T. Hynes and W.P. Reinhardt J. Chem. Phys. 77 (1982) 3595. [8] J.S. Hutchinson, J.T. Hynes and W.P. Reinhardt Chem. Phys. Letters 108 (1984) 353. [9] C. lung and C. Leforestier, J. Chem. Phys. 90 (1989) 3198; Chem. Phys. Letters 155 (1991) 369. [10] D.H. Lu and W.L. Hase, J. Chem. Phys. 89 (1988) 6723. [11] S.M. Lederman, S.J. Klippenstein and R.A. Marcus, Chem. Phys. Letters 146 (1988) 7. [12] A. Garcia-Ayllon, J. Santamaria and G.S. Ezra, J. Chem. Phys. 89 (1988) 801. [13] G.A. Voth, R.A. Marcus and A.H. Zewail, J. Chem. Phys. 86 (1987) 6000. [14] J. Segall, R.N. Zare, H.R. Dubal, M. Lewrenz and M. Quack, J. Chem. Phys. 86 (1987) 634. [15] J.S. Hutchinson, W.P. Reinhardt and J.T. Hynes, J. Chem. Phys. 79 (1983) 4247. [16] E.J. Holler, J. Chem. Phys. 62 (1975) 1544. [17] J. Kucar and H.D. Meyer, J. Chem. Phys. 90 (1989) 5566. [18] R.D. Coalson and M. Karplus, Chem. Phys. Letters 90 (1982) 301; S.Y. Lee and E.J. Heller, J. Chem. Phys. 76 (1982) 3035. [19] K.G. Kay, J. Chem. Phys. 91 (1989) 170. [20] H.D. Meyer, U. Manthe and L.S. Cederbaum, Chem. Phys. Letters 165 (1990) 73. [21] R.B. Gerber and M.A. Ratner, Advan. Chem. Phys. 70 (1988) 97; M.A. Ratner and R.B. Gerber, J. Phys. Chem. 90 (1986) 20. [22] R.D. Coalson, Chem. Phys. Letters 165 (1990) 443; M. Messina and R.D. Coalson, J. Chem. Phys. 92 (1990) 5297, 5712. [23] M. Durga Prasad, Chem. Phys. Letters 194 (1992) 27. [24] R.H. Bisseling, R. Kosloff, R.B. Gerber, M.A. Ratner, L. Gibson and C. Cerjan, J. Chem. Phys. 87 (1987) 2760. [25] K.L. Sebastian, Phys. Rev. B 32 (1985) 6976. [26] M. Durga Prasad, J. Chem. Phys. 88 (1988) 7005; Lecture Notes in Chemistry 50 (1989) 321; G. Madhavi Sastry and M. Durga Prasad, Theoret. Chim. Acta 89 (1994) 193. [27] J. Arponen, Ann. Phys. 151 (1983) 311. [28] S. Guha and D. Mukherjee, Chem. Phys. Letters 186 (1991) 84. [29] C. Bloch, Nucl. Phys. 6 (1958) 329;
222
G. Sree Latha, M. Durga Prasad / Chemical Physics Letters 241 (1990) 215-222
C. Bloch and J. Horowitz, Nucl. Phys. 8 (1958) 91; J. Descloizeu, Nucl. Phys. 20 (1960) 321. [30] G. Sree Latha and M. Durga Prasad, Theoret. Chim. Acta 86 (1993) 511. [31] D. Mukherjee, Intern. J. Quantum Chem. $20 (1986) 409.
[32] D. Sinha, S.K. Mukhopadhyay, R. Chaudhury and D. Mukherjee, Chem. Phys. Letters 154 (1989) 544. [33] B.H. Brandow, Rev. Mod. Phys. 39 (1967) 771. [34] D. Mukherjee and S. Pal, Advan. Quantum. Chem. 20 (1988) 291.