Tatiana Korona
email: tatiana.korona@chem.uw.edu.pl

Coupled cluster singles and doubles polarization propagator accurate through the third order of Møller-Plesset theory

Keywords: coupled cluster theory, molecular properties, polarization propagator, linear response function Coupled cluster theory, molecular properties, polarization propagator, linear response function

Coupled cluster singles and doubles polarization propagator accurate through the third order of Møller-Plesset theory

Tatiana Korona

Introduction

Coupled cluster (CC) theory [1] has become one of the most popular tools for a description of dynamic electron correlation in molecules. For instance, the CCSD(T) method (CC theory limited to single and double excitations with a noniterative inclusion of triple excitations) [2] has been denoted as a 'golden standard' of quantum chemistry for small-size systems. A notable feature of CC theory is its sizeextensivity [3,4], which is also preserved in truncated CC models. A success of CC theory in describing dynamic electron correlation for the energy creates expectations for a similarly good description of molecular properties. Contrary to the correlation energy, some of these properties, like dipole moments and polarizabilities, can be directly compared to experimental values. It is well known that the development of CC methods for describing the molecular properties requires a lot of additional effort in comparison to the calculation of the CC correlation energy. The main difficulty arises from the nonvariational and projective character of the CC equations [5] (for instance, the CC energy is not manifestly real, if a truncated cluster operator T is used). One can deal with the signalized problems in two inequivalent ways, which can be in general denoted as derivative and expectationvalue methods, respectively. It is noteworthy that both approaches give the same results in the limit of full CC theory.

The first method, which has become a standard nowadays, has been pioneered by Monkhorst [6] and Dalgaard and Monkhorst [7] and is based on an analytic differentiation of the CC energy followed by an elimination of the derivative amplitudes [8][START_REF] Bartlett | Analytical evaluation of gradients in coupled cluster and many body perturbation theory[END_REF][START_REF] Salter | [END_REF][11] or alternatively on the use of the CC energy Lagrangian [12,13] [giving the approach known as gradient or linear-response (LR) CC theory]. This method has been extended to dynamic properties by the use of a time-dependent pseudoenergy [6,7], and is therefore denoted frequently as time-dependent (TD) CC theory of molecular properties. The development of the TD-CC approach has been greatly facilitated by an introduction of the time-dependent CC pseudoenergy Lagrangian [14][15][16]. In the last 15 years the implementations of TD-CC2 [17], TD-CCSD [16], orbital-relaxed TD-CCSD [18], TD-CC3 [19], and even TD-CCSDT [20] and arbitrary order [21] frequency-dependent linear-response properties (also called polarization propagators [22]) and their numerous applications (see e.g. Refs. [23][24][25][26][27][28][29][30]) have been reported in the literature. It should be noted that although the orbital-relaxed TD-CCSD approach has been recently developed [18], the standard nonrelaxed version of TD-CC theory should be rather used, since the relaxation of orbitals results in an additional set of poles (generated by the coupled Hartree-Fock equations) for real frequencies [16].

An approach alternative to TD-CC makes a direct use of the basic quantum mechanics expression for the average value of some operator (or its generalizations for higher-order properties) with the exact wave function replaced by its normalized CC counterpart, and will be called shortly XCC theory in this paper. The development of this method has been hindered for a long time by complications stemming from the presence of the deexcitation operator T † (an adjoint of the usual cluster operator T of CC theory), resulting in an infinite series of products of T † and T (see e.g. Ref. [31]). The solution to this problem -proposed some time ago by Jeziorski and Moszynski [32] -relies on the use of a new excitation operator S (vide infra). In the same publication general formulae for the first-order expectation-value CC properties have been proposed. Recently, this approach has been extended to dynamic second-order properties [33]. The corresponding polarization propagator for CC theory limited to single and double excitations (CCSD) has been derived and implemented in Ref. [34]. It should be emphasized that the XCC approach does not require a costly calculation of the left-hand solution of CC equations (the so-called Λ vector [8]), so the XCC first-order properties can be calculated at virtually no cost, once the T operator is available [35]. In addition, XCC theory makes it possible to introduce cumulants of the reduced two-electron density matrices [36] and it should be used in all cases where properly defined CC cumulants are required, like in the expressions for the exchange energy contributions [37][38][39] of symmetry-adapted perturbation theory (SAPT) of intermolecular interactions [40].

The exact polarization propagator exhibits various symmetry properties, like Hermitian (or time-reversal for purely real or purely imaginary perturbations) and operator exchange symmetries. The time-reversal symmetry is not fulfilled by the TD-CC polarization propagator in its original formulation [14], while the XCC propagator satisfies by definition the Hermitian and time-reversal symmetries, but the truncated forms of this propagator show small (high-order) violations of the operator exchange symmetry [33].

Both TD-CC and XCC formulae for polarization propagators include some terms containing connected triple excitations, which contribute already in the third order of Møller-Plesset (MP) theory. As a result, the TD-CCSD and XCCSD methods for the second-order properties are rigorously correct through the second MP order only [16,33,41]. The aim of this paper is to derive a low-cost formula accounting for the missing MP3 contributions by employing the operator techniques, which enable the elimination of long triple energy denominators, and to propose an extension of this formula for the nonperturbative case.

TD-CC and XCC propagators

General expressions for TD-CC and XCC propagators have been developed elsewhere [14,33]. In this section only the final formulae in a form suitable for the further analysis will be reviewed. For a convenience a scalar product notation: X = Φ|XΦ and X|Y = XΦ|Y Φ (where Φ is the reference determinant) will be used. Operators X and Y acting as perturbing operators for second-order properties are assumed to be Hermitian one-electron operators throughout this paper.

TD-CC polarization propagator

A general formula for the TD-CC propagator for the Hamiltonian H, perturbing operators X, Y , and for ω being a frequency of the perturbing field, is given by the following equation [14] X

; Y ω = 1+Λ|[Y , Ω(ω)] + 1+Λ|[X, Ξ(-ω)] + 1+Λ|[[H, Ω(ω)], Ξ(-ω)] , (1)
where Z denotes a similarity-transformed operator

Z = e -T Ze T . (2)
An excitation operator Λ, which appears already in the formulation of the firstorder LR-CC properties [8,[START_REF] Bartlett | Analytical evaluation of gradients in coupled cluster and many body perturbation theory[END_REF]12,13], is defined as a solution of the following equation

1 + Λ|[H, τ ] = 0, (3)
where τ can be any excitation operator within a considered excitation range [e.g. in the present study excitations are limited to singles, doubles and triples (CCSDT)].

A first-order operator Ω(ω) describes the response of a molecule to the perturbing field X oscillating with a frequency ω (note that Ω(ω) corresponds to t X (-ω) in Ref. [16]). The operator Ω(ω) can be obtained by solving the following linear equation

τ |X + τ |[H, Ω(ω)] + ω τ |Ω(ω) = 0. ( 4
)
An analogous equation for Ξ(-ω) can be retrieved from Eq. ( 4) by replacing X by Y and ω by -ω.

The fact that Eq. (1) does not fulfill the time-reversal symmetry has been realizednoticed already in Ref. [14]. In order to force the correct symmetry usually a symmetrized form of Eq. ( 1) is usually utilized, as proposed in Ref. [15]. This solution has a side-effect that some high-order contributions (in terms of the MP fluctuation operator W = H -F , where F is the Fock operator) appear with a wrong coefficient of 1 2 [33], what however does not lead to any known practical problems. If the excitations are restricted to singles, doubles and triples and the nested commutator expansion

e -T Ze T = Z + [Z, T ] + 1 2! [[Z, T ], T ] + . . . , (5)
is utilized, then a rather lengthy formula for the TD-CCSDT propagator results from Eq. ( 1). Since the present study is focused on the analysis of the low-order parts of the TD-CCSDT propagator, only the contributions relevant to the further derivations will be considered. These contributions are listed explicitly in Table 1.

The same table also reports the lowest order of MP theory, to which a given term contributes in the propagator expression. Finally, the table shows whether listed terms are included already on the TD-CCSD level. It follows from an examination of Table 1 that there is one contribution which is entirely missing, while two more terms are included only partially if triple excitations are neglected. The latter terms miss the triples contributions, which enter indirectly through the single and double excitation operators. Therefore, the TD-CCSD polarization propagator is rigorously correct only through the second order of MP theory [16].

XCC polarization propagator

An alternative approach for the coupled cluster second-order properties has been presented recently by Moszynski et al. [33]. This method is based on a general definition of the polarization propagator [22],

X; Y ω = - ∞ n=1 Ψ 0 |Y Ψ n Ψ n |XΨ 0 E n -E 0 + ω - ∞ n=1 Ψ 0 |XΨ n Ψ n |Y Ψ 0 E n -E 0 -ω (6)
where E k and Ψ k are eigenvalues and eigenfunctions for the kth electronic state, respectively. It is noteworthy that the second term of the r.h.s. of Eq. ( 6) can be obtained from the first one either by taking its complex conjugate followed by replacing of ω by -ω * (the operation called generalized complex conjugate -g.c.c. in Ref. [33], which is related to the time-reversal symmetry), or by exchanging X and Y operators and replacing of ω by -ω (the operation called operator exchange -op.ex.). In the approach of Moszynski et al. [33] one utilizes the normalized CC wave function as Ψ 0

Ψ 0 = e T Φ e T |e T 1/2 , (7)
and the following Ansatz for the first-order wave function Ψ (1) , being a linear response to the perturbation X oscillating with a frequency ω

Ψ (1) = (Ω 0 (ω) + Ω(ω))e T Φ e T |e T 1/2 . ( 8
)
In Eq. ( 8) Ω 0 (ω) is a number (which does not enter the final equations and does not need to be considered any further) and Ω(ω) is an excitation operator defined by Eq. (4). Additionally, an excitation operator S, appearing for the first time in expectation-value CC theory for the first-order properties [32], will beis used in the XCC polarization propagator formula. This operator is uniquely defined by

It has been proven [32,42] that the S operator is connected (it can be constructed from a finite number of multiple commutators containing T and T † ), what makes it a useful quantity for constructing of size-extensive models. The S operator correct through the third W order is given by the equations [32]

S 1 = T 1 + P 1 ([T † 1 , T 2 ] + [T † 2 , T 3 ]), S 2 = T 2 + 1 2 P 2 ([[T † 2 , T 2 ], T 2 ]), S 3 = T 3 . (10)
Starting from Eqs. ( 6)-( 9) the following formula has been derived in Ref. [33] for the CC polarization propagator

X; Y ω = e -S e T † Y e -T † e S | P(e S † Ω(ω) e -S † ) + g.c.c. (11)
In Eq. ( 11) P stands for the superoperator P 1 + P 2 + P 3 + . . ., defined such that P n "projects" the n-tuple excitation part of an operator Z, i.e.,

P n (Z) = (n!) -2 i1...in a1...an e i1...in an...a1 Z e an...a1 i1... in , (12)
and e an...a1 i1... in = a a1 • • • a an a in • • • a i1 is the n-body spinorbital excitation operator (a ak and a ik are a virtual spinorbital creation operator and an occupied spinorbital annihilation operator, respectively). The form of Eq. ( 11) makes it straightforward to apply the expansion from Eq. ( 5). The resulting formula contains a finite, although extremely long series of multiple commutators of the S † and T operators. Most of these terms give a very small contribution to the total value of the polarization propagator and can be safely neglected. The resulting truncated models are also size-extensive since they are constructed as the sums of multiple commutators of connected quantities. The usual criteria for the importance of such terms are based on the analysis of the orders of T 's and W 's [33]. These criteria have been used to develop several truncation schemes for the CCSD case [34], from which the CCSD(3) model has been selected as the optimal one for many-electron cases. The CCSD(3) model includes a minimum number of terms, which are necessary to account for all contributions up to the O(W 3 ) order arising from Eq. ( 11), if T , Ω, and Ξ are truncated after single and double excitations. It is important to emphasize that this model does not include all contributions of the O(W 3 ) order entering the exact XCC polarization propagator. Similarly as for the TD-CCSD case, missing terms arise from triple excitation operators T 3 , Ω 3 (ω), and Ξ 3 (-ω).

All components of the XCCSDT polarization propagators giving a contribution of the O(W 3 ) or lower order are shown in Table 2. As in Table 1, the leading MP order and the presence of a given term in the XCCSD (equivalent to the CCSD(3) model from Ref. [34]) formula are also indicated in this table.

Some contributions presented in Tables 1 and2 look quite different. However, it can be shown by an order-by-order MP analysis that they contain the same or very similar contributions through the third MP order.

CC polarization propagator correct through the third MP order

It is clear from Tables 1 and2 that triple excitation operators are necessary, if a propagator correct through the third MP order is required. Obviously, the presence of triple excitations increases the computational cost of the CC propagator from the sixth to at least the seventh power of the molecular size N , what in turn reduces substantially sizes of molecules for which such calculations can be performed. In this section it will be shown how to introduce the corrections to XCCSD and TD-CCSD polarization propagators, which account for missing O(W 3 ) contributions, but do not contain the T 3 , Ω 3 (ω), or Ξ 3 (-ω) excitation operators.

Møller-Plesset expansion of excitation operators

An analysis of the contributions from Tables 1 and2 is based on the Møller-Plesset expansion of the excitation operators. Such expansions can be most easily obtained from the following iterative formulae [32,33,43]

T n = R n,0 (e -T W e T ), Ω n (ω) = R n,ω (e -T Xe T + [e -T W e T , Ω(ω)]). ( 13
)
The action of the resolvent superoperator R n,ω on ansome operator ZA is given by the equation

R n,ω ( @ @ ZA) = (n!) -2 i1...in a1...an (ǫ i1 +. . .+ǫ in -ǫ a1 -. . .-ǫ an -ω) -1 e i1.
.. in an...a1 @ @ ZA e an...a1 i1... in ,

where ǫ p is a spinorbitalan energy forof a spinorbital p. As shown in Ref. [33], the R n,ω superoperator fulfills the identity (for any operator A)

[F, R n,ω (Z Z XA)] + ω R n,ω (Z Z XA) + P n (Z Z XA) = 0. ( 15
)
By iterating Eq. ( 13) the following expressions can be obtained for consecutive orders of T (k) and Ω (k) (ω), where k denotes the power of the fluctuation operator W T

(1)

2 = R 2,0 (W ), T (2) n = R n,0 ([W, T (1)
2 ]), n = 1, 2, 3,

T (3) n = R n,0 ([W, T (2)
1 + T (2) 2 + T (2) 3 ] + (1 -δ 1n ) 1 2 [[W, T (1)
2 ], T

2 ]), n = 1, 2, 3,

Ω (0) 1 (ω) = R 1,ω (X), Ω (1) n (ω) = R n,ω ([X, T (1)
2 ] + [W, Ω (0) 1 (ω)]), n = 1, 2, Ω (2) n (ω) = R n,ω ([X, T (2) n + T (2) n+1 ] + [W, Ω (1)
1 (ω) + Ω (1) 2 (ω)] + [[W, T (1)
2 ], Ω (0)

1 (ω)]), n = 1, 2, Ω (2)
3 (ω) = R 3,ω ([X, T (2)
3 ] + [W, Ω (1)
2 (ω)] + [[W, T (1)
2 ], Ω (0)

1 (ω)] + 1 2 [[X, T (1)
2 ], T

2 ]), Ω

1 (ω) = R 1,ω ([X, T (3)
1 + T (3) 2 ] + [W, Ω (3)
1 (ω) + Ω (2) 2 (ω) + Ω (2) 3 (ω)] + + [[W, T (2)
2 ], Ω (0)

1 (ω)] + [[W, T (2)
1 ], Ω (0)

1 (ω)] + [[W, T (1)
2 ], Ω

1 (ω)]). ( 16) The formulae for Ω (k) (ω) (k ≤ 2) and for T (k) have been already presented elsewhere and are only repeated for the sake of completeness. Note that in Eq. ( 41) of Ref. [33] the last term is missing in the definition of Ω

F
(2)

3 (ω). The expressions for Ξ (k) n (-ω) are analogous to Ω (k) n (ω) and will be not listed here. The Λ n operator can be expanded in terms of the n-body excitation operators e a1...an i1 ... in , denoted here shortly as τ K(n) ,

Λ n = K(n) Λ K(n) τ K(n) . (17)
The iterative formula for the expansion coefficient Λ K(n) has the following form

Λ K(n) = (∆ǫ K(n) ) -1 1 + Λ|[W , τ K(n) ] , (18)
where ∆ǫ K(n) is the energy denominator from Eq. ( 14) for ω = 0. By iterating Eq. ( 18) one can easily verify that the Λ and T operators have the same MP expansion up to the second MP order and that differences between them appear first on the MP3 level

Λ (1) 2 = T (1) 2 , Λ (2) n = T (2) n , n = 1, 2, 3, Λ (3)
K(1) = (∆ǫ K(1) ) -1 [[W, T (2) 1 ], τ K(1) ] + T (1) 2 |[[W, T (1) 2
], τ K(1) ] + + T (2) 1 |[W, τ K(1) ] + T (2) 2 |[W, τ K(1) ] . (19)
It is important to observe that no contribution from the triple excitation operator T

(2) 3

enters the expression for Λ

1 . Since for the XCC propagator also the S operator is needed, its MP expansion is also of interest here. It can be deduced from Eq. ( 10) that the S n and T n operators are identical through the O(W 2 ) order. FThe first differences between S and T occur in the third MP order, when terms resulting from the commutator expansion containing the T † deexcitation operator start to contribute to S, S

1 = T (3) 1 + P 1 ([(T (2) 1 ) † , T (1) 2 ] + [(T (1) 2 ) † , T (2) 3 ]), S (3) 2 = T (3) 2 + 1 2 P 2 ([[(T (1) 2 ) † , T (1) 2 ], T (1) 2 ]). ( (3)
) 20
From an analysis of Eqs. ( 16)-( 20) one can identify fragments of single and double excitation operators dependent on the triple excitation operators T

3 and Ω

(2)

3 (ω). The explicit expressions for these operator fragments are listed below (the small letter t in superscript denotes that only triple-related parts of the operator are taken into account) As already mentioned, Λ

T (3t) n = R n,0 ([W, T (2) 3 ]), n = 1, 2, Ω (2t) 2 (ω) = R 2,ω ([X, T (2) 3 ]), Ω (3t) 1 (ω) = R 1,ω ([X, T (3t) 1 + T (3t) 2 ] + [W, Ω (2t) 2 (ω) + Ω (2) 3 (ω)]), S (3t) 1 = T (3t) 1 + P 1 ([(T (1) 2 ) † , T (2) 3 ]). (21)
= 0. By substituting operators from Eq. ( 21) into the terms from Tables 1 and2 it is now easy to select the contributions which depend on triple excitation operators and are therefore not included into the XCCSD and TD-CCSD propagators.

The elimination of triple excitation operators from the MP3 contributions of the XCC polarization propagator

The MP3 contributions of the XCC polarization propagator, missing within the XCCSD treatment are collected in Table 3. In Ref. [32] it has been demonstrated that long triple denominators from the MP-expanded expression for the first-order CC properties orders can be eliminated on the operator level. Similar manipulations with double excitation operators have been utilized for the XCC propagator in order to separate true and apparent correlation contributions up to the second MP order [33]. The elimination of long energy denominators can be most easily explained on the algebraic level (see e.g. Ref. [44]). If a sum of two similar expressions is computed, which differ only by various placements of the energy denominators, and if the denominators for the first and second contribution have the form a -1 (a + b) -1 and b -1 (a + b) -1 , respectively (a and b denote sums of orbital energies), then the sum of these two denominators equals simply to a -1 b -1 , i.e. the long (a + b) -1 denominator dissapears. Such techniques have been used to derive various orbital expressions, like polarizabilities for the second order of MP theory [45]. An important advantage of the approach presented in Refs. [32,33] consists in performing such eliminations on the operator level, i.e. without the necessity of keeping track of a multitude of algebraic expressions, which would appear if the Wick theorem or the contraction theorem is applied to the contributions collected in Table 3. In the remainder of this section this approach will be extended and applied to eliminate the T

and Ω

3 (ω) operators from the MP3 contributions to the XCCSDT and TD-CCSDT propagators. In order to facilitate the search for appropriate pairs of terms from Table 3 a diagrammatic analysis can be utilized. A short description of this technique, as well as diagrammatic representations of most important derivation steps from the next Section, are presented in the Appendix. 3 (ω) operator can be rewritten in a short form as

Ω (2) 3 (ω) = R 3,ω (Z), (22)
where the operator Z is given by the equation

Z = [X, T (2)
3 ] + [W, Ω (1)
2 (ω)] + [[W, T (1)
2 ], Ω (0)

1 (ω)] + 1 2 [[X, T (1)
2 ], T

]. ( (1) 2
) 23
The operator Ω

(2)

3 (ω) appears in the terms Id and VII. The term Id can be rewritten as follows

Id = Y | R 1,ω ([W, Ω (2) 3 (ω)]) = R 1,ω * (Y )|[W, Ω (2) 3 (ω)] = = P 2 (W ) R 1,ω * (Y )|Ω (2) 3 (ω) = P 2 (W ) R 1,ω * (Y )| R 3,ω (Z) . (24)

Molecular Physics w3˙rev

Molecular Physics

Here the "Hermicity" of the superresolvent operator has been used (A and B denote some general operators)

Z Z XA| R n,ω (@ @ Y B) = R n,ω * (Z Z XA)|@ @ Y B , (25)
as well as the facts that (i) the commutator [W, Ω

3 (ω)] can be replaced by the product W Ω The second term containing the Ω

3 (ω) operator can be rewritten as

VII = T (1) 2 |[Y, Ω (2)
3 (ω)] = T

(1)

2 P 1 (Y )|Ω (2) 3 (ω) = R 2,0 (W ) P 1 (Y )| R 3,ω (Z) , (26)
if one makes use of the Hermicity of the Y operator and notices that in this case Y plays a role of a single deexcitation operator. In the last step the definitions of the T

(1) 2

and Ω

(2) 3 (ω) operators (Eqs. ( 16) and ( 22), respectively) have been used. Note that the operators on the bra side of Eq. ( 26) commute, since they are both the excitation operators.

In order to continue the following identity should be employed

[F, R A 1n,ω (Z Z XA) R A 2k,0 (@ @ Y B)]+ω R A 1n,ω (Z Z XA) R A 2k,0 (@ @ Y B)+ P A 1n (Z Z XA) R A 2k,0 (@ @ Y B)+ R A 1n,ω (Z Z XA) P A 2k (@ @ Y B) = 0. (27)
Eq. ( 27) has been derived by performing an addition of two equations: the first one obtained by the multiplication of Eq. ( 15) by R A 2k,0 (@ @ Y B) from the right and the second one retrieved from Eq. ( 15) by setting Z Z XA = @ @ Y B, n = A 2k, and ω = 0, and by the multiplication of the resulting expression by R n,ω (Z Z XA) from the left. Applying Eq. ( 27) to the sum of the terms Id and VII gives

Id + VII = R 1,ω * (Y ) P 2 (W ) + P 1 (Y ) R 2,0 (W )| R 3,ω (Z) = = -R 1,ω * (Y ) R 2,0 (W )|[F, R 3,ω (Z)] -ω R 1,ω * (Y ) R 2,0 (W )| R 3,ω (Z) . (28)
In view of Eq. ( 15) the commutator in Eq. ( 28) can be replaced by -ω R 3,ω (Z) -P 3 (Z). As a result the terms multiplied by ω cancel each other and finally the following expression is obtained

Id + VII = R 1,ω * (Y ) R 2,0 (W )|Z = Ξ (0) 1 (ω * )T (1) 2 |Z , (29)
where in the last step the definitions of the Ξ (0) 1 (ω * ) and T

(1) 2

operators have been used.

In this way the first triple excitation operator has been eliminated. After applying Eq. ( 23) the formula (29) separates into four fragments

Id + VII = Ξ (0) 1 (ω * )T (1) 2 |[X, T (2) 3 ] + Ξ (0) 1 (ω * )T (1) 2 |[W, Ω (1) 2 (ω)] + + Ξ (0) 1 (ω * )T (1) 2 |[[W, T (1) 2 ], Ω (0) 1 (ω)] + 1 2 Ξ (0) 1 (ω * )T (1) 2 |[[X, T (1)
2 ], T

2 ] . (30)

The last three terms of Eq. ( 30) already have the desired form, i.e. they do not contain any triple excitation operators. The diagrammatic analysis of possible con- tractions allows to rewrite the first part of the r.h.s. of this equationEq. (30) as a sum of two fragments

Ξ (0) 1 (ω * )T (1) 2 |[X, T (2)
3 ] = = P 1 (XΞ to replace the MP operators by their converged CCSD counterparts, what results in the following formula,

(0) 1 (ω * ))T (1) 2 |T (2) 3 + Ξ (0) 1 (ω * ) P 2 ([X, T (1)
∆(XCC) = Ξ 1 (ω * )T 2 |[W, Ω 2 (ω)] + [[W, T 2 ], Ω 1 (ω)] + 1 2 [[X, T 2 ], T 2 ] + + Σ 1 |Y Ω 1 (ω) + XΞ 1 (ω * )|Σ 1 + + [W, T 2 ]|Ω 2 (ω)Ξ 1 (-ω) + [W, T 2 ]|Ω 1 (ω)Ξ 2 (-ω) + g.c.c., (52)
where the Σ 1 operator is defined as

[F, Σ 1 ] + P 1 (T † 2 [W, T 2 ]) = 0. ( 53
)
Eq. ( 52) is the final result of this section. It includes all O(W 3 ) contributions of the XCC polarization propagator, but does not contain triple excitation operators and can be calculated forat the O(N 6) cost (a more detailed analysis reveals that the most expensive parts of this formula scale as O(o 3 v 3 ), where o and v are dimensions of occupied and virtual orbital spaces, respectively). If the correction ∆(XCC) is added to the XCCSD polarization propagator, a low-cost propagator correct through the third MP order is obtained.

The elimination of triple excitation operators from the MP3 contributions of the TD-CC polarization propagator

A similar correction formula can be proposed for the TD-CCSD polarization propagator. The MP3 terms missing in the TD-CCSD formula are collected in Table 4.

Since they are identical to some terms from Table 3, the same roman labels are preserved. It is important to mention that the component 3 from Table 1 does not contribute to Table 4, because the operator Λ (3t) 1

is equal to zero, while a similar component 3 from Table 2 gives rise to the terms II and III of Table 3.

An identification of the correction terms for the TD-CCSD case is, however, not as straightforward, as it appears from the first look at Table 4, since for each of listed contributions its op.ex. counterpart instead of the g.c.c. counterpart should be added. Since in the derivations presented in Sections 3.2.1-3.2.4 the g.c.c. of some terms from Table 3 have been used in several places, these cases should be reexamined once more in order to check if the new appropriate pairing terms can be found.

Terms Ia, Ib, Ic, Id, VI, and VII for the TD-CCSD case

Let us first analyse the terms from the first line of Table 4, which have appeared in Sections 3.2.1 and 3.2.3. In Section 3.2.1 the g.c.c. of the term IVa (Eq. ( 36)) is utilized as the counterpart of one of the terms arising during the derivation. A closer look at Table 4 reveals that this term is identical to the op.ex. of the term VI. Since all other terms from Section 3.2.1 are explicitly listed in Table 4, the final result of this section is valid also for the TD-CC case. The remaining terms from Table 4 are the terms Ib and Ic, while in Section 3.2.3 the following terms: the g.c.c. of the term Ib, the g.c.c. of the term Ic, IVb, and V have been considered. Alternatively, one could make the derivations from Section 3.2.3 for the g.c.c. of the abovementioned terms, i.e. for the terms Ib, Ic, the g.c.c. of IVb, and the g.c.c. of V (they are included in Eq. (51) anyway under the abbreviation g.c.c.). After some operator manipulations it can be shown that the g.c.c. of IVb is identical to the op.ex. of Ic and the g.c.c. of V is equal to the op.ex. of Ib, so the g.c.c. of the final result of Section 3.2.3 can be used also for the TD-CCSD triples correction. term can be rewritten as [W, Λ 2 ]|Ω 1 (ω)Ξ 2 (-ω) , i.e. it represents exactly the last term of Eq. ( 52).

Summary

The analysis of the most important contributions for the TD-CC and XCC polarization propagators has been performed by means of MP theory. The TD-CCSD and XCCSD polarization propagators, which in their original formulations are rigorously correct through the second order of MP theory only, have been augmented with noniterative corrections ∆(XCC) and ∆(TD), respectively, accounting for the missing triple contributions on the MP3 level. In the final formulae triple excitation operators have been eliminated by utilizing techniques which enable the removal of long energy denominators on the operator level. Since missing triple contributions of the leading MP3 order are responsible for a less satisfactory performance of the TD-CCSD and XCCSD polarization propagators as compared to the CCSD correlation energy [which is rigorously correct through the O(W 3 ) order], it can be expected that the proposed correction terms will improve the agreement between the CCSD polarization propagators and the full configuration interaction (FCI) benchmarks. The scaling of the computational time for Eqs. ( 52) and ( 55) is at least one order of magnitude lower than for the case of an explicit use of triple excitation operators (e.g. the computational time for the TD-CC3 polarization propagator scales as iterative O(o 3 v 4 ) [19], while the proposed corrections can be obtained with the noniterative O(o 3 v 3 ) cost). The low cost of the ∆(XCC) and ∆(TD) corrections will allow to perform calculations for CC polarization propagators accurate through the O(W 3 ) order for much larger molecules and orbital basis sets than it is possible with e.g. the TD-CC3 method. The comparison of the XCC and TD-CC approaches on the operator level reveals that that some terms, which appear in the XCC model only when triple excitations are included, are already incorporated in the original TD-CCSD polarization propagator, either with help of the double commutator term Λ 2 |[[H, Ω 1 (ω)], Ξ 2 (-ω)] , or through the Λ 1 operator. However, numerical tests show that errors of both methods with respect to the FCI benchmarks are of the same order of magnitude, although on average they are somewhat smaller for the TD-CCSD case [34]. The implementation of Eqs. ( 52) and ( 55) is planned in the near future.

Appendix: On the use of diagrams for the derivation of the ∆(XCC) correction

The derivations presented in Sections 3.2.1-3.2.4 can be greatly simplified if the diagrammatic approach is utilized as an auxiliary tool. For this purpose it is enough to employ the arrowless (nonoriented) Hugenholtz diagrams [46]. In the following it will be demonstrated on one selected example from Section 3.2.1 how one can utilize diagrams in order to find the counterparts of some term under study. The following drawing conventions will be used for operator vertices: a small filled circle denotes the W operator, an empty (white) circle -the T operator, a cross -the X operator, a crossedsmall circle with the letter Y inside -the Y operator, a black squarethe Ω operator and finally a large filled circle -the Ξ operator. Thin vVertical lines denote energy denominators from Eq. ( 14), with a corresponding frequency listed above or below the denominator line. Additionally, auxiliary operators, like the W T n operators or the Z operator defined in Eq. ( 23) are represented by an ellipse with the corresponding symbol inside. A left-to-right drawing convention is used, i.e. for a formula ABC . . . the order of vertices A, B, C, etc. corresponds to their order in the formula. As usual, if operators commute the order of their vertices in a diagram can be changed and e.g. these vertices can be drawn one over another, if convenient. Let us reiterate the search for the complementary terms, which has been the main issue in Sections 3.2.1 to 3.2.4, by utilizing the diagrammatic methods.

Let us now take as an example the first term of the r.h.s. of Eq. ( 31). This term is shown on Figure 1 ) . This term represents exactly the term Ia from Table 3. In this way the counterpart of the first term of the r.h.s. of Eq. ( 31) has been identified. Other cases from Sections 3.2.1-3.2.4 can be proceeded in a similar way.

The first term examined in Section 3.2.1 is the term Id (see Eq. ( 24)), which is shown in Figure 1(a). Its rightmost energy denominator can be eliminated, if this diagram is considered together with a diagram presented in Figure 1(b), having the reverse pattern of the short (i.e. corresponding to single or double excitations) denominator lines. The frequency in the latter diagram (i.e. zero) has been obtained from the condition that the total frequency should be equal to the frequency in the long energy denominator. In the diagram shown in Figure 1 24)- (29). In order to clarify the correspondence between the diagrammatic and algebraic approaches, let us write explicitly the algebraic expression obtained from the diagram presented in Figure 1(a). For this purpose this diagram has been decorated with arrows and line labels (indices labelling the occupied and virtual spinorbitals are denoted by letters i, j, k and a, b, c, respectively). The algebraic form of the diagram from Figure 1(a) (with amplitudes Y a i etc. implicitly defined) is given below

Id = ijkabc Y a i (ǫ i -ǫ a -ω) -1 W bc jk (ǫ i + ǫ j + ǫ k -ǫ a -ǫ b -ǫ c -ω) -1 Z ijk abc . (57)
It can be easily verified that the expression for the neighbour diagram (see Figure 1(b)) differs only by the absence of the term (ǫ iǫ aω) -1 and the presence of (ǫ j + ǫ kǫ bǫ c ) -1 . The rest of the procedure for the elimination of the triple energy denominator has been described at the beginning of Section 3.2.

In the following the diagrams containing the T

(2) 3

operator will be examined in detail. For clarity reasons the full diagrammatic representation of this operator (with the reduced resolvent R 3,0 and the [W, T

2 ] commutator, cf. Figure 2(a)) will be not used here. Instead, either a simple T vertex for T

. Components of the TD-CCSDT propagator giving a contribution in a low W order. For the components containing Y and Ω(ω) and for the term 10 the respective op.ex. terms should be added (the op.ex. term is obtained by a simultaneous exchange of the operators X and Y and replacing ω by -ω). The lowest MP order to which a given term contributes is presented in the second column. The last column contains an information if these terms are fully included on the CCSD level (the abbreviation Y denotes that the whole term is included, N -that the whole term is missing, P -that some parts of this term are missing).

Nr MP order Term in TD-CCSD?

1 0 Y Ω 1 (ω) P 2 1 [[W, Ω 1 (ω)], Ξ 1 (-ω)] Y 3 2 Λ 1 |[Y , Ω 1 (ω)] Y 4 2 Λ 2 |[Y , Ω 2 (ω)] P 5 2 Λ 2 |[Y , Ω 1 (ω)] Y 6 2 Λ 2 |[[W , Ω 1 (ω)], Ξ 1 (-ω)] Y 7 3 Λ 1 |[Y , Ω 2 (ω)] Y 8 3 Λ 2 |[Y , Ω 3 (ω)] N 9 3 Λ 1 |[[W , Ω 1 (ω)], Ξ 1 (-ω)] Y 10 3 Λ 2 |[[W , Ω 1 (ω)], Ξ 2 (-ω)] Y Table 2.
Components of the XCCSDT propagator giving a contribution in a low W order. For each term the corresponding g.c.c. counterpart should be added. For explanations of the second and fourth columns see Table 1.

Nr MP order Term in XCCSD? (ω) from Eq. ( 21) (the order of terms from Eq. ( 21) is preserved). Similarly, an insertion of Ω

1 0 Y |Ω 1 (ω) P 2 1 [Y, S 2 ]|Ω 1 (ω) P 3 2 [Y, S 1 ]|Ω 1 (ω) P 4 2 [Y, S 2 ]|Ω 2 (ω) P 5 2 [[T † 2 , Y ], S 2 ]|Ω 1 (ω) Y 6 3 Y |[S † 1 , Ω 2 (ω)] Y 7 3 Y |[S † 2 , Ω 3 (ω)] N 8 3 [[T † 1 , Y ], S 2 ]|Ω 1 (ω) Y 9 3 [[T † 2 , Y ], S 3 ]|Ω 1 (ω) N 10
2 (ω) (Eq. ( 16)) leads to the terms IVa and IVb.

3 (ω)]) , where Eq. ( 22) has been utilized for Ω

(2)

3 (ω); b) its guessed counterpart T (1) 2 |[Y, Ω (2)
3 (ω)] , which can be identified as the term VII; c) the T |Y Ω (0) 1 (ω) , where the diagrammatic representation of the operator (T (3t) 1

) † , cf.

3

|Ω (0) 1 (ω) b P 2 ([Y, T (1)

  e T † e T Φ e T |e T .

3. 2 . 1 .

 21 Terms Ia, Id, IVa, VII Let us start from two terms containing the perturbed triple excitation operator Ω

  )W generates a pentuple excitation, which obviously cannot be deexcited by a single R 1,ω * (Y ) operator), and (ii) the W operator should deexcite two out of three excitations produced by the Ω

  * ) replaced by R 1,ω * (Y ). From Eq. (34) it follows that the appropriate "pairing" counterpart for this term has the reverse pattern of the denominator lines, like depicted on Figure2(a). It is easy to see that the diagram from Figure 2(a) includes the T (3t) 1 operator, presented alone on Figure 2(b). This in turn allows to identify the whole diagram from Figure 2(a) as Y | R 1,ω (XT (3t) 1

  (b) the Hermitian adjoint of the doubly-excited cluster operator T (1) 2 can be easily seen (see Figure 1(c)), while the remaining part of the diagram contains the Y operator contracted to the Ω (2) 3 (ω) operator, so that finally this diagram can be identified as the term VII. The sum of both terms gives the term depicted in Figure 1(d), i.e. it represents the diagram where both short energy denominator lines are present, but the troublesome long line vanishes. The same result on the operator level has been obtained in a sequence of Eqs. (

Figure 2 (

 2 b), or a simplified version of diagram shown in Figure 2(b), as given in Figure 2(c), will be utilized.

F

Figure 1 .

 1 Figure 1. Diagrams replacing the derivation steps from Eq. (24) to Eq. (29): a) the term Id , i.e. Y | b R 1,ω ([W, Ω

Figure 2 .Figure 4 .) 2 (

 242 Figure 2. Diagrammatic representations of the T (2) 3operator.

Figure 5 .Figure 6 .

 56 Figure 5. The diagrams explaining Eqs. (31) and (36): a) the second part of the r.h.s. of Eq. (31) with the expanded form of Ξ (0) 1 (ω * ) (see Eq. (16)); b) its counterpart b R 2,-ω * ([X, T (1) 2 ])|[Y, T(2)3 ] , which can be identified as g.c.c. of the term IVa.

Figure 4 (Figure 7 .Figure 8 .) 2 (

 4782 Figure 7. The operators T (3t) 2 (a) and Ω (2t) 2 (ω) (b), see Eq. (21).

Figure 9 .

 9 Figure 9. The diagrams explaining Eqs. (48)-(50): a) the term VIII, T

  i.e. g.c.c. of the term VI. Page 25 of 25 URL: http://mc.manuscriptcentral.com/tandf/tmph Molecular Physics

Table 3 .

 3 Terms arising from triple excitation operators of the XCCSDT polarization propagator contributing at the MP3 order. For each term its g.c.c. counterpart should be added. Terms Ia, Ib, Ic, and Id are obtained by inserting the definition of Ω

	(3t)
	1

URL: http://mc.manuscriptcentral.com/tandf/tmph Molecular Physics

where we use the fact that the commutator [X, T

3 ] can be replaced by the operator product XT

(2) 3 in this case, and that only hole-hole and particle-particle parts of the X operator contribute to Eq. (31).

Let us now examine the term Ia. A sequence of transformations listed below changes this term to the form suitable for our purposes

If the result of derivations presented in Eq. (32) and the first term of the r.h.s. of Eq. ( 31) are summed up, the following expression results

(2) 3

+ R 1,0 (XΞ

2 |[W, T

.

In the derivation of Eq. ( 33) the identity [33] P k (@ @ Y B) R n,ω (Z Z XA) + R k,-ω (@ @ Y B) P n (Z Z XA) + [F, R k,-ω (@ @ Y B) R n,ω (Z Z XA)] = 0, (34) has been employed for ω = 0, k = 1, and n = 2.. Eq. ( 34) ishas been obtained by using the same procedure as for the case of Eq. (27). Eq. ( 33) does not contain any triple excitation operators, but can be still somewhat simplified. If the operator T

(1) 2

is moved to the ket side and the Hermicity of the R n,0 superoperator is used, then the r.h.s. of Eq. ( 33) can be rewritten as

1 , where the auxiliary operator Σ

(3) 1 is defined as

Note that the Σ

1 operator is identical with the second term of the Λ

1 operator (see Eq. ( 19)). This part of the derivation is represented diagrammatically in the Appendix.

Finally, it can be found from the examination of Table 3 that the counterpart (in a sense of Eq. ( 34)) for the second fragment from the r.h.s. of Eq. ( 31) is the g.c.c. of the term IVa, which can be rewritten as

The sum of these two terms gives R 2,-ω * ([X, T

.

Terms II and III

The techniques presented in Section 3.2.1 can be applied to the terms II and III giving

and

3 |T

(1)

As before, by applying a sequence of Eqs. ( 34) and ( 15) the sum of these two terms can be transformed to the following expression

Finally, if the T

(1) 2

operator in Eq. ( 39) is moved to the bra side and the definition of the operator Σ

(3) 1 is again used (see Eq. ( 35)), then the final form of Eq. ( 39) can be written as

and of the g.c.c. of the term Ib

results in the following expression

3 ]|Ω

WhenIf the operator W is moved to the ket side and the diagrammatic analysis of possible contractions is performed, the following two terms are obtained

The first of the terms on the r.h.s. of Eq. ( 44) can be paired with the term IVb,

while the second term of Eq. ( 44) makes a pair with the g.c.c. of the term Ic,

respectively. Summarizing, the final result of this Section is listed below

2 ]|Ω

Terms VI and VIII

The remaining components in Table 3 are the terms VI and VIII. The g.c.c. of the term VI can be rewritten as

and similarly the term VIII can be transformed to

in a few steps explained in detail in the previous sections. It is easy to see that the summation of these two terms leads to the result

2 ]|Ω

2 ]) ,

what completes the derivation of the correction terms for the XCCSD polarization propagator.

The final correction formula for the XCCSD polarization propagator

When the partial results from Sections 3.2.1-3.2.4 are collected together and the definitions of the Ω (1) 2 (ω) and Ξ

(1) 2 (-ω) operators from Eq. ( 16) are utilized, then the final correction formula for the triple contributions of the third MP order can be written as

2 ], T

]|Ω

(1)

2 ]|Ω

where g.c.c. denotes, as usual, that for each term its g.c.c. counterpart should be added in order to obtain the complete formula: note that in the derivations presented in Sections 3.2.1-3.2.4 either the term k, or the g.c.c. of the term k (but never both terms simultaneously) has been employed. Eq. ( 51) includes the lowest nonvanishing MP orders of the T , Ω, and Ξ operators. It can be therefore proposed

It should be stressed that in this case both terms have been utilized together with their op.ex. counterparts, so the g.c.c. of the r.h.s. of Eq. ( 47) should appear only once in the final formula.

The final correction formula for the TD-CCSD polarization propagator

After collecting the results from Sections 3.2.1 and 3.2.3 the following triples correction for the TD-CC polarization propagator can be proposed

2 ], T

+ Ω

(1)

2 ] + op.ex.,

where op.ex. denotes that for each term its op.ex. counterpart should be added in order to restore the complete formula. After replacing the MP operators by their CCSD counterparts the following expression is obtained for the triple-excitation correction for the TD-CCSD case,

Eq. ( 55) is the main result of this section. Similarly as for the XCCSD case, if the correction ∆(TD) is added to the TD-CCSD propagator, the expression accurate through the O(W 3 ) order is obtained at the O(N 6 ) computational cost.

A comparison of the XCCSD and TD-CCSD correction formulae

The formulae (51) and (54) (or alternatively (52) and ( 55)) are similar, but not identical, e.g. it is evident that the second equation contains less terms than the first one. In particular, one term with the Σ

1 operator and one term containing the product of the first-order single and double excitation operators are absent in the second expression. In order to explain this fact let us first analyse triples contributions to the S (3) 1 operator (see Eq. ( 21)). One of the terms from Eq. ( 21) originates from the T (3) 1 operator, while the second one -from the T † T part of the S operator. The sum of these two terms projected on some arbitrary single excitation operator τ 1 gives

where once more the sequence of Eqs. ( 34) and ( 15) have been employed. The result of these derivations is identical with the second term of Eq. ( 19) for the Λ 1 . However, if the XCCSD polarization propagator is calculated (the CCSD(3) model from Ref. [34]), an approximate form of S 1 is utilized, which is constructed from the T 1 and T 2 operators of CCSD theory. The components absent for the CCSD case (see Eq. ( 56)) give rise to the terms II and III from Table 3 and are responsible for one contribution containing Σ

(3) 1 in Eq. (51). The second "missing" term in Eq. ( 55) can be identified as the term 10 from Table 1. After some manipulations it can be demonstrated that this

The next identity which can be clarified by the diagrammatic approach is the result presented in Eq. (31). Since in the l.h.s. of this equation the one-electron operator X should be connected to the T (2) 3 vertex (because of the presence of the commutator with T

(2)

3 ), the remaining line of the X vertex connects either to the Ξ (0) 1 (ω * ), or to the T (1) 2 vertex, giving two terms presented in Figure 3. The first term of the r.h.s. of Eq. ( 31) is shown in Figure 4(a) with the T 33)).

We proceed similarly for the second term from the r.h.s. of Eq. ( 31). Its diagram is presented in Figure 5(a) with the Ξ (0) 1 (ω * ) operator replaced by R 1,ω * (Y ), while the complementary diagram is depicted in Figure 5(b). The latter diagram can be identified as g.c.c. of the term IVa. For the sake of brevity the final diagrams obtained from the summation of the diagram pairs will be not shown in the rest of this Appendix, since the way of the elimination of long energy denominators is strictly analogous to this presented in Figure 4.

In Section 3.2.2 a counterpart for the term II is looked for. The term II is shown in Figure 6 and Ω (2t) 2 (ω). These operators are presented in Figure 7. For the sake of brevity the diagrammatic representation of Eqs. ( 41)- (43), showing the elimination of the double energy denominator of T (3t) 2 is skipped here and only the result of these derivations, already splitted into two possible diagrams, is given in Figures 8(a


Joonsuk Huh

Michael Neff

Guntram Rauhut

Robert Berger

wavefunctions

Keywords: Franck-Condon, VSCF/VCI, photodetachment-photoelectron spectra, CCSD

teaching and research institutions in France or abroad, or from public or private research centers.

Bowman et al. 14 used vibrational self-consistent field (VSCF) and vibrational configuration interaction (VCI) methods together with the Watson Hamiltonian to compute FC factors (FCFs). They employed, however, one set of primitive harmonic oscillator basis functions for both states, which renders the computation of their corresponding overlap integrals trivial but typically requires much longer VCI expansions to describe accurately the intensity profile near the origin of the vibronic band (0'-0 transition region).

Luis et al. 13,22 used either vibrational perturbation theory or VSCF/VCI functions, apparently only for the final electronic state, and determined the Franck-Condon factors by solving subsequently a homogeneous linear system of equations. They employed, however, only a simple diagonal kinetic energy operator without vibrational angular momentum terms and without the Watson correction term. These authors have reported applications even to systems with nine atoms, although with a significantly reduced set of mode-mode couplings.

Rodriguez-Garcia et al. 15 also used this simplified kinetic energy operator and compared two types of approaches, namely (i) the use of harmonic oscillator basis functions adapted to each electronic state involved, which thereby require to compute overlap integrals between distorted, displaced and rotated harmonic oscillators, and (ii) VCI calculations using a common set of basis functions for both electronic states. The vibronic spectra of the systems studied, however, were largely dominated by the 0'-0 transition which thus seems to be a particularly favourable situation for a common basis set.

The photodetachment process of disulfanide can (depending on the final electronic state) be accompanied by a relatively large change in the equilibrium structure due to the shortening of the S-S bond length (∆r S-S ∼ 0.1 Å), as a relatively long vibrational progression is visible in the corresponding photoelectron spectrum of the anion, thereby posing a problem which appears somewhat better suited for computational methods that employ individual vibrational basis functions for the vibrational motion in the initial and final electronic states.

We thus present herein our integrated vibronic structure approach for the calculation of the FC profile of the photodetachment-photoelectron spectra. This approach, which we will outline in the next section, shares some features of the one by Mok et al. 12 , in that we are using separate sets of basis functions in the two electronic states and employing the Watson rovibrational Hamiltonian including the pseudopotential-like Watson correction term and the vibrational angular momentum terms for the solution of the vibrational problem. As in most previous work, all additional angular momentum coupling terms like electron spinrotational, electron spin-electron orbit and nuclear rotation-electron orbit couplings were dropped. However, these are necessary for high-accuracy benchmark calculations of vibrational wavenumbers, see for example the work of Werner and others on the water cation 23 .

We deviate from the ansatz of Mok et al. by using a VSCF/VCI ansatz for the description of the numerous vibrational wavefunctions and by employing for each electronic state a separate set of distributed Gaussians as primitive basis functions for the VSCF calculations. The latter functions allow, due to their locality, simple prescreening strategies in the calculation of the FC integrals of n-dimensional anharmonic oscillators to be exploited. Potential energy surfaces are obtained by making use of a fully automated surface construction code 24,25 (see below). For the present application, the potential energy surfaces were computed at the (U)CCSD(T)-F12a/vtz-f12 level.

II. THEORY

As outlined above our approach for obtaining vibrational wavefunctions is based on the Watson Hamiltonian 26 for polyatomic non-linear molecules.

Ĥ = 1 2 αβ πα µ αβ πβ - 1 8 α µ αα - 1 2 i ∂ 2 ∂q 2 r + V (q 1 , . . . , q 3N -6 ) (1)
Like the potential, the 2nd term, the so-called Watson correction term, is represented in a many-mode expansion and is added as a pseudopotential-like contribution to the potential V (q 1 , . . . , q 3N -6 ) expressed in normal coordinates q i . For details see Ref. 24. We used a distributed Gaussian basis χ µ

χ µ (q) = 4 2A µ π exp{-A µ (q -q µ ) 2 } (2)
for representing the one-mode wavefunctions φ i (q i ) (modals) within the VSCF and VCI approaches. The parameters A µ were determined as described in the seminal paper of Hamilton and Light 27 . As normal modes are distinguishable, the VSCF many-mode wavefunction Φ can be expressed as a simple Hartree product of the modals, i.e. where n denotes the occupation number vector (ONV) with elements n i and we assume the expansion coefficients of the modals to be real numbers. The corresponding VCI solution is a simple linear combination of such wavefunctions with different ONVs. Based on a separability approach, the 3N -6 dimensional eigenvalue problem is reduced to 3N -6 one-dimensional problems within the VSCF method. Once a polynomial representation of the potential has been chosen, the resulting one-dimensional effective polynomial, p(i) r , can be expressed as 28

Φ (n) (q) = i φ (n i ) i (q i ) with φ (n i ) i (q i ) = µ C (n i ) µi χ µ (q i ) (3)
p(i) r = p (i) r + j s X js p (ij) rs + k t X kt 1 2 p (ijk) rst + . . . (4)
where r denotes the order of the polynomial. The effective polynomial exclusively depends on the fitting coefficients p and the one-dimensional potential integrals:

Q r µν = χ µ |q r | χ ν and X ir = µν C µi C νi Q r µν (5)
Note that the integrals Q r µν are mode-independent, provided that the same Gaussian basis set is employed for all modes. The VSCF/VCI approach outlined here was successfully used in many applications but usually was limited to the calculation of fundamental modes or low lying overtones and combination bands [29][30][31] .

The basic building block of the vibronic structure calculation with the VSCF and VCI wavefunctions is the overlap integral between rotated, distorted and displaced Gaussian functions. The normal coordinates of the two different electronic states are connected by the Duschinsky relation 10

q ′ = Sq + d ( 6
)
where S is the Duschinsky rotation matrix and d is the displacement vector. With the quantities defined by Doktorov et al. 32 , i.e.

R = (I + J t J) -1 J t , P = J(I + J t J) -1 J t , J = Ω ′ SΩ -1 and δ = Ω ′ d (7) in which Ω = diag((2A µ1 ), . . . , (2A µ(3N -6) )) 1/2 and Ω ′ = diag((2A ′ µ ′ 1 ), . . . , (2A ′ µ ′ (3N -6) )) 1/2
, where prime ( ′ ) is used to specify the parameters belonging to the final electronic state, the multi-dimensional overlap integral (characterised by the vectors µ ′ and µ with elements µ ′ i and µ i , respectively) of the Gaussian basis functions is given as The norm of the overlap integral has a trivial upper bound in the limiting case of identical Gaussian exponent parameters between the initial and final electronic state Gaussian basis set functions and in the absence of Duschinsky rotation, which is

i χ ′ µ ′ i (q ′ i )χ µ i (q i ) dq = 2 (3N -6)/2 |det(R)| 1/2 exp - 1 2 δ t (I -P)δ . (8
i χ ′ µ ′ i (q ′ i )χ µ i (q i ) dq ≤ exp - 1 4 ||δ|| 2 (9)
and only depends on the norm of the delta vector. This simple bound condition is exploited in our overlap integral prescreening and integrals (with norm estimate) below 10 -8 were ignored.

The VSCF wavefunction overlap integral is expressed in terms of the Gaussian overlap integrals

Φ ′ (n ′ ) (q ′ )Φ (n) (q) dq = µ,µ ′ i C (n ′ i ) µ ′ i C (n i ) µ i i χ ′ µ ′ i (q ′ i )χ µ i (q i ) dq. ( 10
)
The VCI wavefunction of the m-th vibrational state is expanded in terms of the VSCF wavefunction with real expansion coefficients

B ( n) m Ψ m (q) = n B (n) m Φ (n) (q). (11)
With this, the corresponding overlap integrals of the VCI wavefunctions reads as

Ψ ′ m ′ (q ′ )Ψ m (q) dq = n,n ′ B ′ ′ (n ′ ) m ′ B (n) m µ,µ ′ i C ′ (n ′ i ) µ ′ i C (n i ) µ i i χ ′ µ ′ i (q ′ i )χ µ i (q i ) dq. ( 12
)
In the evaluation of the VCI overlap integral the terms of B

′ (n ′ ) m ′ B (n)
m below 10 -8 were dropped in the summation.

III. COMPUTATIONAL DETAILS A. Electronic Structure Calculations

Electronic structure calculations were performed at the (U)CCSD(T)-F12a/vtz-f12 level.

Explicitly correlated coupled-cluster calculations are much better suited for the generation of accurate potential energy surfaces than conventional (U)CCSD(T) calculations, because less demanding orbital basis sets can be used. This reduces the I/O bottleneck in the corresponding Hartree-Fock calculations, which may be severe once many ab initio calculations need to be performed. The RMP2-F12 calculations, which are the first step in (U)CCSD(T)-F12a calculations, were performed using the MP2-F12/3C(FIX) method 33,34 .

In this method the numerous two-electron integrals are computed using robust density fitting (DF) approximations 35,36 (DF was not used in the Hartree-Fock and UCCSD(T)-F12 calculations). The aug-cc-pVTZ/MP2FIT basis set of Weigend et al. 37 was used as auxiliary basis. For evaluating the Fock matrix (which is needed in the AO and RI basis sets) the vtz/OPTRI basis sets 38 was employed. The complementary auxiliary basis set (CABS) approach was employed, i.e., the union of the AO and RI basis sets was used to approximate the resolution of the identity. The perturbative CABS singles correction as described in Refs. 39,40 was applied in all F12 calculations. This significantly reduces the Hartree-Fock basis set error. Careful tests showed that the choice of the density fitting basis has only a negligible effect on the vibrational wavenumbers. The choice of the RI basis is more critical; however, we found that employing larger RI basis sets does not reduce the statistical errors.

Tight CABS thresholds (10 -9 ) were used in order to minimise the number of functions that are deleted due to near linear dependencies. The number of functions was kept constant for all displacements. This guarantees smooth potential energy surfaces. For an application of explicitly correlated coupled-cluster theory, CCSD(T)-F12, on the calculation of vibrational transitions see Ref. 31.

B. Calculation of Potential Energy Surfaces

Potential energy surfaces were represented by a multi-mode expansion in terms of normal coordinates rather than a Taylor series expansion. These surfaces were generated in a fully automated fashion based on an iterative interpolation algorithm as described in detail elsewhere 24,25 . In a first step all surfaces were represented by 20 grid points in each direction but a subsequent transformation to polynomials up to 8th order was used for the calculation of the fundamental modes 28 .

C. VSCF/VCI Calculations

Once the potential had been determined -either in terms of grid points or polynomialsthe vibrational self-consistent field approximation was used for calculating the fundamental modes of HS distributed Gaussians. In the VSCF calculations, vibrational angular momentum terms were included perturbatively using the converged VSCF wavefunction and a constant µtensor as discussed by Carbonniere and Barone 30,41 . Vibrational correlation effects were included subsequently by vibrational configuration interaction calculations including single, double and triple excitations, leading to a total of 216 VCI states 28 . In the VCI calculations, vibrational angular momentum terms were included variationally.

D. Vibronic Structure Calculations

Franck-Condon factors (FCFs) for the ( X

1 A ′ )HS - 2 → ( X 2 A ′′ )HS 2 , ( X 1 A ′ )HS - 2 → ( Ã 2 A ′ )HS 2 , ( X 1 A ′ )DS - 2 → ( X 2 A ′′ )DS 2 and ( X 1 A ′ )DS - 2 → ( Ã 2 A ′
)DS 2 photodetachmentphotoelectron spectra were computed using the VSCF/VCI vibrational wavefunctions obtained as described above and, for comparison, using approximate three-dimensional harmonic oscillator wavefunctions. As these wavefunctions were obtained using normal modes specific for each electronic state, Duschinsky mode mixing effects 10 had to be taken into account. To this end, a development version of the MOLPRO program package 42 used for the electronic and vibrational structure calculations reported herein was interfaced to a development version of the hotFCHT software 9,[START_REF] Berger | [END_REF]44 which was employed for subsequent vibronic structure studies. VSCF based vibronic spectra reported below were computed by using the modals obtained in a VSCF calculation for the vibrational ground state of each electronic state. Excited vibrational states were approximated within an uncoupled Hartree framework by products involving excited modals obtained in the vibrational ground state SCF procedure and the energies of the various states (including the ground state) were approximated by sums of the corresponding modal energies. None of these vibrational state energies were subsequently corrected perturbatively for angular momentum contributions. Thus, shifts between VSCF and VCI vibronic band positions are due to vibrational correlation effects and the vibrational angular momentum terms included in the latter approach but neglected in the former. As the changes in equilibrium structures and normal modes are comparatively large in the X 1 A ′ → X 2 A ′′ transitions, the use of specific normal modes appears beneficial as it reduces the length of the VCI expansions (number of VCI basis functions) required to accurately describe the intensity pattern in the 0'-0 transition region. The price to pay is the increased cost in the calculation of individual overlap integrals between basis func- This situation is significantly different from a previous study along similar lines reported by Rodriguez-Garcia et al. 15 , where the vibronic spectra were largely dominated by the 0'-0 transition band and thus structural changes were also very small. Besides prescreening on the level of multidimensional overlap integrals between Gaussian basis functions (Eq. ( 10))

and on the level of the VCI expansion coefficients, one can also envisage an additional prescreening stage for the overlap integrals between VSCF wavefunctions (in Eq. ( 12)) that will also impact on the cross-over point for the common basis set and individual basis set approaches. Low-cost a priori estimates of the respective numerical effort would be desirable, but are beyond the scope of the present work.

IV. RESULTS AND DISCUSSION

Bond lengths, harmonic wavenumber and the fundamentals of HS I andII. With respect to the S-S bond lengths, our UCCSD(T)-F12a/vtz-f12 results deviate slightly from those of Peterson et al. 3 . However they do agree nicely with those obtained at the complete basis set limit but without core correlation and and high-order correlation effects, which are not included in our calculations. Except for the H-S-stretching modes our VCI results are in excellent agreement with the computed anharmonic wavenumbers of Peterson et al. 3 This holds true in particular for the deuterated species. In general, we consider the results of these authors to be more accurate than our's, because they did not only account for high-order electron correlation effects (in terms of CCSDTQ calculations), but also incorporated core-correlation effects and scalar-relativistic corrections. As our intention is not to reproduce the spectroscopic properties of these molecules most accurately but to study routes for determining the vibrational finestructure of their photodetachment-photoelectron spectra, we did not correct our potential energy surfaces for these effects. However, a comparison of the results of Peterson et al. error bar in the experiments. We believe that the accuracy of our calculations is sufficient in order to simulate the vibrational structure of the photodetachment-photoelectron spectra.

The computed Franck-Condon (FC) profiles of the computed ( X 1 A ′ ) → ( X 2 A ′′ ) photodetachment-photoelectron of HS - 2 and DS - 2 are displayed in Figs. 1 and2, respectively. For each of these transitions we show the stick representation for the FC profile (corresponding to the individual FCFs of the various transitions) as well as the FC profile as obtained after convolution with Lorentzian lineshape function with full width at half maximum (FWHM) of 200 cm -1 . The full band listing up to 4000 cm -1 above the 0'-0 transition wavenumber is reported in the appendix.

For both species the computed main spectral feature of this transition is a pronounced progression 3 n 0 in the S-S stretching mode ν 3 that extends (visibly) up to about n ′ = 5. This progression has also been identified by Entfellner et al. [START_REF] Entfellner | [END_REF] in the experimental low resolution photodetachment-photoelectron spectrum, with the progression extending as well up to about n ′ = 5. The FCFs of the 3 1 0 and 3 2 0 transition are the largest of the spectrum, with the former being less than 5 % (less than 13%) larger than the latter in HS 2 (in DS 2 ). The FCFs of the 0 0 0 transition and the 3 3 0 transition are of similar magnitude and both are predicted to be half the size of the 3 1 0 within the VSCF and VCI framework. The harmonic approximation instead predicts about a 4:7 ratio ((0 0 0 , 3 3 0 ) : (3 1 0 , 3 2 0 )) and also produces a longer progression in the stretching mode (up to about n ′ = 6 or n ′ = 7). A direct comparison of the intensity pattern computed for this progression with experiment 1 is hampered by the somewhat low experimental resolution and overlapping signals due to additional photoprocesses, that have been attributed to a concomitant photodetachment-photoelectron spectrum of S - 2 , which is formed by photodissociation of the parent anionic compound.

The second prominent feature, which remained, however, unresolved in the experimental spectrum 1 , is the progression 2 1 0 3 n 0 that builds upon the excitation of the bending mode ν 2 . The FCF for the 2 1 0 transition is in the deuterated species more pronounced than in the non-deuterated parent compound due to the larger projection of the change in equilibrium structure, that is caused by the electron detachment, onto the DS 2 bending mode.

Excitation in ν 1 , the H-S or D-S stretching, mode are connected with very small FCFs (at most one per mil of the FCF of 0 0 0 ). This is primarily due to the almost negligible change of the corresponding bond length upon electron detachment, which is predicted by theory.

The computed FC profile of the photoelectron spectrum for the photodetachment process that leads to the energetically lowest excited doublet state of HS 2 and DS 2 (( X 1 A ′ ) → ( X 2 A ′′ ), see Figs. 3 and4) reveals significantly less features. The main progression 2 n 0 is due to the bending mode ν 2 , but is very short and ends essentially at n ′ = 2. The most intense band is the 0 0 0 transition, with the other bands falling off rapidly in intensity. This is in good agreement with the experimental findings of Entfellner et al.. [START_REF] Entfellner | [END_REF] Excitations of the S-S stretching mode ν 3 are strongly suppressed due to the almost equal S-S bond length in the initial and final electronic state. The FCF of the 3 1 0 transition is only about 1/60 (1/15) of the FCF of 0 0 0 in HS 2 (DS 2 ). These bands were not resolvable in the experimental spectrum.

Even weaker is the computed 1 n 0 progressions in the H-S and D-S stretching mode ν 1 . The FCF obtained for 1 [START_REF] Entfellner | [END_REF] 0 is only about one per mil of 0 0 0 . Entfellner et al. [START_REF] Entfellner | [END_REF] assigned, however, a signal in the corresponding wavenumber range to this transition. On the basis of our calculations, we can neither support nor fully rule out this assignment. While the wavenumber range matches, our computed intensity appears too low to justify the assignment. Due to the approximate nature of our calculation, however, where we neither included core-valence correlation effects nor relativistic corrections, the change in bond length upon electron detachment might be underestimated which will result in a too low bandstrength for this transition.

Overall, VSCF and VCI give quite similar FC profiles for the four photoelectronphotodetachment spectra. The VCI wavefunctions are for the most prominent bands usually dominated by a single reference, which also impacts on the resulting intensity pattern. A slightly more pronounced multi-configurational character is observed in particular when higher excitations of the bending mode are involved. Some strongly perturbed transitions are present in the photodetachment spectra of the deuterated species, where a close resonance situation is predicted for 1 1 0 3 2 0 /2 1 0 3 4 0 (and various progressions built on these) in the electron detachment spectrum that involves the electronic ground state of the radical. A similar situation occurs for 1

1 0 2 4 0 /1 1 0 2 3 0 3 1 0 in the corresponding ( X 1 A ′ )DS - 2 → ( Ã 2 A ′ )DS 2 electronic transition.
Differences between harmonic and VSCF/VCI profiles are most prominent in the 3 n 0 progressions (excitations of the S-S stretching mode) in the ( X 1 A ′ ) → ( X 2 A ′′ ) electron detachment process. The anharmonic progression is shorter than the harmonic one, which is not unexpected, as the S-S bond contracts upon removal of the electron. The details of the profile depend on the shape of the potential function. We have to emphasise, however, that the underlying electronic structure calculations are difficult to converge for such large excursions from the equilibrium structure so that in particular the computed vibronic bands in the larger wavenumber region are to be considered as crude approximations. The corresponding data are thus intended primarily for facilitating the reproduction of our results.

Given the current experimental resolution 1 for the photodetachment-photoelectron spectra of HS - 2 and DS - 2 , even the harmonic approximation might be considered a good starting point for the analysis.

V. SUMMARY, CONCLUSIONS AND OUTLOOK

We have reported herein the computed Franck-Condon profile for the photodetachmentphotoelectron spectra of different isotopomers of disulfanide upon which the disulfanyl in its doublet ground and energetically lowest doublet excited state is formed. Three-dimensional anharmonic potential energy hypersurfaces have been computed in a fully automated fashion. Electronic structure calculations were performed at the explicitly correlated coupled clusters singles and doubles level with perturbative triples correction using a correlationconsistent polarized valence triple zeta atomic basis set 38,39 . The approach employed for the vibronic structure calculations combines VSCF/VCI wavefunctions in basis sets adapted to each electronic state, namely displaced Gaussian functions along the normal modes of each electronic state. The non-orthogonal basis functions require the calculation of multidimensional overlap integrals, but have the potential of supporting shorter VCI expansions and facilitating straight forward prescreening techniques.

The ( X 1 A ′ ) → ( X 2 A ′′ ) is dominated by excitations in the S-S stretching mode, whereas various weak additional features are due to an excitation of the bending mode. The computed FC profile is in reasonable agreement with the low-resolution spectrum obtained experimentally. As a result of smaller structural changes upon removal of the electron, the ( X 1 A ′ ) → ( Ã 2 A ′ ) transition is relatively featureless and contains a short progression in the bending mode as its main feature, which is in line with the experimental spectrum. Excitations in the S-S stretching mode are significantly less pronounced and remained unresolved in experiment. The single signal assigned in the experiment to an excitation of the H-S and D-S stretching mode, respectively, is computed to be extremely weak. Nevertheless,

̺ FCW (ν) /cm (ν -ν0 )/cm -1
Harmonic VSCF VCI ν 1 ν 2 ν 3 a ν -ν0 FCF ν 1 ν 2 ν 3 b ν -ν0 FCF ν 1 ν 2 ν 3 (%) c
ν -ν0 FCF 0 0 0 0.00 1.33E-01 0 0 0 0.00 1.31E-01 0 0 0 (99.9) 0.00

Harmonic

VSCF VCI

ν 1 ν 2 ν 3 a ν -ν0 FCF ν 1 ν 2 ν 3 b ν -ν0 FCF ν 1 ν 2 ν 3 (%) c
ν -ν0 FCF 0 0 0 0.00 6.25E-01 0 0 0 0.00 6.17E-01 0 0 0 (99.9) 0.00

Harmonic

VSCF VCI

ν 1 ν 2 ν 3 a ν -ν0 FCF ν 1 ν 2 ν 3 b ν -ν0 FCF ν 1 ν 2 ν 3 (%) c ν -ν0 FCF 0 0 0 0.
Harmonic VSCF VCI ν 1 ν 2 ν 3 a ν -ν0 FCF ν 1 ν 2 ν 3 b ν -ν0 FCF ν 1 ν 2 ν 3 (%) c ν -ν0 FCF 0 0 0 0.00 5

  integrals. The authors employed the Watson rovibrational Hamiltonian in their calculations, which have, up to now, primarily been performed for various non-linear triatomic molecules.

6

 6

2 .

 2 Each modal was represented by a linear combination of 30 7

  the basis functions are non-orthogonal and the corresponding overlap integrals are not separable into simple products of one-dimensional integrals. Prescreening techniques for multidimensional overlap integrals of distributed Gaussians still allow rapid calculation of these integrals such that the reduced expansion lengths may outweight the higher cost.

  at the VCI level are shown in Tables

  to our's shows that these corrections are usually quite small. In most cases deviations with respect to the experimental values 1 are significantly larger which we attribute to the large 9

10

 10

11

 11

  not be fully ruled out by virtue of the inherent approximation in the calculations.While our initial application was limited to a three-dimensional vibronic problem, ongoing work will focus on the extension to medium sized systems.VI. ACKNOWLEDGMENTSRB acknowledges the Volkswagen Foundation for financial support and the Center for Scientific Computing (CSC) Frankfurt for computer time.

Figure 1 : 2 → ( X 2 A 18 Page 19 of 25 URL- 1 Figure 2 : 2 → ( X 2 A 19 Page 20 of 25 URL- 1 Figure 3 : 2 → ( Ã 2 A 20 Page 21 of 25 URL- 1 Figure 4 : 2 → ( Ã 2 A 21 Page 22 of 25 URL

 1221825122219251322202514222125 Figure 1: Franck-Condon profiles of the ( X 1 A ′ )HS - 2 → ( X 2 A ′′ )HS 2 photodetachmentphotoelectron spectrum computed in the harmonic approximation, the vibrational self-consistent field (VSCF) approximation and the vibrational configuration interaction (VCI) approximation. The lower graph shows a stick representation of the averaged Franck-Condon weighted density of states ̺ FCW (ν) whereas the upper graph shows the corresponding profiles obtained after convolution with a Lorentzian lineshape function with full width at half maximum (FWHM) of 200 cm -1 .

  (a) Peak assignment based on harmonic oscillator wavefunctions (b) Peak assignment based on VSCF wavefunctions (c) Peak assignment for the VCI wavefunctions based on the (maximum) percentage contribution (numbers in parentheses) of the corresponding VSCF wavefunction.

  ). (a) Peak assignment based on harmonic oscillator wavefunctions (b) Peak assignment based on VSCF wavefunctions (c) Peak assignment for the VCI wavefunctions based on the (maximum) percentage contribution (numbers in parentheses) of the corresponding VSCF wave function.

  (a) Peak assignment based on harmonic oscillator wavefunctions (b) Peak assignment based on VSCF wavefunctions (c) Peak assignment for the VCI wavefunctions based on the (maximum) percentage contribution (numbers in parentheses) of the corresponding VSCF wavefunction.

Table I :

 I Computed and observed fundamental vibrational wavenumbers (in cm -1 ) of HS 2 and HS - 2 .

		Molecular Physics
		X2 A ′′	Ã2 A ′	X1 A ′
		This work PMF a Exp. b This work PMF a	Exp. b This work
	r(HS)	1.3499 1.3482 1.3523	1.3435 1.3417		1.3472
	R(SS)	1.9646 1.9608 1.9603	2.0807 2.0752		2.0883
	θ(HSS)	101.58 101.52 101.74	93.28 93.25		101.68
	ω 1 ω 2 ω 3 F o 2602.4 2598.3 918.7 919.0 605.0 605.0 r	2677.3 2674.1 771.9 770.5 510.5 513.0		2598.1 831.5 489.7
	ν 1 ν 2 ν 3	2482.4 2474.5 2463 898.2 900.8 904±8 P 598.4 598.4 595±4 e e r	2569.7 2563.2 2550±200 752.6 752.8 750±200 504.6 506.9 504±4	2469.6 809.7 478.4
		R
		e
		v i e
			w
			O n l
				y

a) Taken from Peterson, Mitrushchenkov and Francisco, see Ref. 3 b) Experimental data taken from Refs. 45-48 16 Page 17 of 25 URL: http://mc.manuscriptcentral.com/tandf/tmph

Table II :

 II Computed and observed fundamental vibrational wavenumbers (in cm -1 ) of DS 2 and

				Molecular Physics		Page 18 of 25
		1.0×10 -3
	DS -2 .	8.0×10 -4	X2 A ′′		Ã2 A ′	VCI VSCF X1 A ′ Harmonic
	̺ FCW,L (ν)/cm	This work PMF a Exp. b This work PMF a 1869.6 1866.8 1922.4 1920.1 667.2 667.3 559.3 558.4 4.0×10 ω 1 ω 2 -4 6.0×10 -4	Exp. b This work 601.4 1866.5
	ω 3 ν 1 ν 2 ν 3 F 1807.5 1803.2 603.4 603.6 656.6 657.9 696±20 597.1 597.1 591±10 0.0 r o 2.0×10 -4	508.6 511.0 1866.1 1863.0 1830±160 549.0 549.0 540±160 503.0 505.5 502±15	484.6 1800.9 589.6 478.9
		a) Taken from Peterson, Mitrushchenkov and Francisco, see Ref. 3 P 0 1000 2000 3000 4000 (ν -ν0 )/cm -1	5000
			b) Experimental data taken from Refs. 45-48 e
			e r			VCI
		-1 2.5×10		R			VSCF Harmonic
		2.0×10 -1		e
		5.0×10 -2 1.0×10 -1 1.5×10 -1		v w e i
		0.0	0	1000	2000	O l 3000 4000 n	5000
							y
				17

URL: http://mc.manuscriptcentral.com/tandf/tmph

Table IV :

 IV Computed Franck-Condon factors (for a wavenumber range extending up to ν -ν0 < 4000 cm -1 ) and corresponding band assignments of the photodetachment process from HS - 2 ( X1 A ′ ) to the lowest excited doublet state of HS 2 ( Ã2 A ′ ).

	F o r
		P
		e
			e r
				R
				e
	0 0 1 0 1 0 0 0 2 0 1 1 0 0 3	604.98 918.73 1209.96 1523.71 1814.94	2.28E-01 9.12E-03 2.22E-01 1.67E-02 1.59E-01	0 0 1 0 1 0 0 0 2 0 1 1 0 0 3	599.14 916.96 1193.05 v 2.78E-01 8.96E-03 2.70E-01 i e 1516.10 2.01E-02 1781.72 1.58E-01 w 0 0 1 (99.8) 0 1 0 (98.9) 0 0 2 (99.8) 0 1 1 (98.4) 0 0 3 (99.7)	598.38 898.29 1191.51 1492.16 1779.39	1.31E-01 2.78E-01 1.27E-02 2.68E-01 2.53E-02 1.57E-01
	0 2 0	1837.46	9.38E-04	0 2 0	1852.80	8.68E-04	0 2 0 (94.7)	1785.99	4.83E-04
	0 1 2 0 0 4 0 2 1 1 0 0 0 1 3 0 3 0 0 0 5 0 2 2 1 0 1	2128.69 2419.92 2442.44 2602.40 2733.67 2756.19 3024.90 3047.42 3207.38	1.71E-02 9.37E-02 1.73E-03 2.54E-06 1.29E-02 7.30E-05 4.77E-02 1.80E-03 3.99E-06	0 1 2 0 0 4 0 2 1 1 0 0 0 1 3 0 3 0 0 0 5 0 2 2 1 0 1	2110.01 2365.14 2451.95 2472.92 2698.68 2802.89 2943.29 3045.85 3072.06	2.08E-02 6.26E-02 1.98E-03 1.12E-06 1.30E-02 5.96E-05 1.78E-02 2.07E-03 2.10E-06	0 1 2 (98.0) 0 0 4 (99.6) 0 2 1 (93.1) 1 0 0 (96.9) 0 1 3 (97.2) 0 3 0 (87.5) 0 0 5 (99.7) O n 2080.82 2362.08 2375.21 2482.42 2664.25 2666.94 2941.73 l 0 2 2 (91.7) 2959.32 1 0 1 (96.3) 3082.35 y	2.26E-02 6.17E-02 1.19E-03 5.53E-06 1.21E-02 1.72E-06 1.74E-02 1.12E-03 5.90E-06
	0 1 4	3338.65	7.95E-03	0 1 4	3282.10	5.56E-03	0 1 4 (97.0)	3242.78	4.10E-03
	0 3 1	3361.17	1.39E-04	1 1 0	3389.88	5.36E-06	0 3 1 (84.8)	3251.56	2.31E-05
	1 1 0	3521.13	6.82E-06	0 3 1	3402.03	1.44E-04	1 1 0 (90.4)	3365.17	1.20E-06
	0 0 6	3629.88	2.18E-02	0 0 6	3516.18	3.76E-03	0 2 3 (90.4)	3538.26	6.49E-04
	0 2 3	3652.40	1.37E-03	0 2 3	3634.52	1.32E-03	0 4 0 (80.7)	3579.96	1.02E-06
	0 4 0	3674.92	6.00E-06	1 0 2	3665.97	1.78E-06	1 0 2 (95.8)	3676.99	1.97E-06
	1 0 2	3812.36	3.57E-06	0 4 0	3763.19	3.78E-06	0 1 5 (96.6)	3826.60	9.63E-04
	0 1 5	3943.63	4.22E-03	0 1 5	3860.25	1.72E-03	0 3 2 (81.4)	3831.30	1.48E-05
	0 3 2	3966.15	1.48E-04	1 1 1	3989.02	1.13E-05	1 1 1 (89.6)	3960.72	1.84E-06
				0 3 2	3995.94	1.59E-04
					22

Table V :

 V Computed Franck-Condon factors (for a wavenumber range extending up to ν -ν0 < 4000 cm -1 ) and corresponding band assignments of the photodetachment process from DS - 2 ( X1 A ′ ) to ground state DS 2 ( X2 A ′′

	0 0 1 0 1 0 0 0 2 0 1 1 0 0 3 F 1021.04 510.52 771.91 1282.43 1531.56 o r	8.15E-03 2.94E-01 4.65E-04 2.83E-03 1.09E-05	0 0 1 0 1 0 0 0 2 0 1 1 0 0 3	504.86 776.33 1005.23 1281.18 1501.07	9.35E-03 2.99E-01 1.92E-04 3.42E-03 1.39E-06	0 0 1 (99.8) 0 1 0 (98.8) 0 0 2 (99.8) 0 1 1 (98.7) 0 2 0 (94.8)	504.53 752.59 1004.52 1254.83 1495.01	6.17E-01 1.02E-02 3.06E-01 6.49E-05 2.01E-03 5.84E-02
	0 2 0	1543.82	5.80E-02	0 2 0	1585.30	5.88E-02	0 0 3 (99.7)	1499.93	1.04E-06
	0 1 2 0 0 4 0 2 1 0 1 3 0 3 0 0 0 5 0 2 2 1 0 0	1792.95 2042.08 2054.34 P 1.79E-04 4.42E-07 3.63E-04 2303.47 3.43E-06 2315.73 6.08E-03 e 2552.60 1.19E-08 2564.86 2.81E-05 2677.34 5.73E-04 e 0 1 2 0 0 4 0 2 1 0 1 3 0 3 0 0 0 5 1 0 0 0 2 2 r	1781.56 1992.37 2090.16 2277.40 2418.17 2479.37 2557.85 2590.54	6.19E-05 2.64E-09 4.55E-04 3.29E-07 7.47E-03 3.31E-11 1.86E-04 6.99E-06	0 1 2 (98.5) 0 0 4 (99.6) 0 2 1 (94.4) 0 3 0 (88.2) 0 1 3 (98.2) 0 0 5 (99.6) 0 2 2 (94.0) 1 0 0 (98.1)	1752.29 1990.79 1994.97 2232.08 2244.93 2478.07 2489.93 2569.62	9.15E-07 8.97E-08 2.03E-04 4.44E-03 1.65E-08 1.95E-10 2.46E-07 5.36E-04
	0 1 4	2813.99	1.47E-07	0 1 4	2768.69	1.77E-11	0 3 1 (87.4)	2730.01	2.83E-05
	0 3 1 0 0 6 0 2 3 0 4 0 1 0 1 0 1 5 0 3 2 1 1 0 0 0 7 0 2 4 0 4 1 1 0 2 0 1 6	2826.25 3063.12 3075.38 3087.64 3187.86 3324.51 3336.77 3449.25 3573.64 3585.90 3598.16 3698.38 3835.03	1.86E-05 4.23E-10 3.92E-07 3.44E-04 2.71E-06 3.43E-09 2.24E-06 2.00E-03 1.20E-11 1.92E-08 1.61E-07 2.28E-07 1.25E-10	0 3 1 0 0 6 1 0 1 R 2923.02 2963.91 3062.71 0 2 3 3086.38 0 1 5 3255.69 e 0 4 0 3268.03 1 1 0 3334.18 0 3 2 3423.40 0 0 7 3453.01 1 0 2 3563.08 v 3.29E-05 1.03E-12 3.89E-07 2.33E-08 1.96E-11 8.78E-04 1.28E-03 4.21E-07 6.92E-13 3.02E-09 i e 0 2 4 3577.67 1.28E-10 0 1 6 3740.24 1.26E-12 0 4 1 3772.89 1.78E-06 w 0 1 4 (97.8) 0 2 3 (93.4) 0 4 0 (83.3) 1 0 1 (98.0) 0 1 5 (97.3) 0 3 2 (86.3) 1 1 0 (93.0) 0 2 4 (93.1) 0 4 1 (82.6) 1 0 2 (97.9) 0 3 3 (85.8) 0 5 0 (74.9) 1 1 1 (92.5)	2733.31 2979.84 3037.67 3073.72 3220.20 3222.95 3303.62 3465.95 3534.95 3573.32 3710.80 3801.82 3805.88	1.41E-07 6.37E-09 8.43E-05 4.66E-06 9.36E-12 2.20E-07 7.19E-04 4.96E-09 6.53E-06 6.92E-08 1.39E-09 1.06E-06 1.22E-06
	0 3 3 0 5 0 1 1 1	3847.29 3859.55 3959.77	1.69E-08 8.73E-06 1.68E-05	1 1 1 0 3 3 0 0 8	3839.03 3919.24 3960.83	1.26E-05 5.87E-10 2.04E-13	0 2 5 (93.0) O n 3953.55 l	8.96E-10
							y
					23

Table VI :

 VI Computed Franck-Condon factors (for a wavenumber range extending up to ν -ν0 < 4000 cm -1 ) and corresponding band assignments of the photodetachment process from DS - 2 ( X1 A ′ ) to the lowest excited doublet state of DS 2 ( Ã2 A ′ ).

	Page 25 of 25			Molecular Physics
	00 603.44 667.15 1206.88 1270.59 1334.30 F 0 0 1 0 1 0 0 0 2 0 1 1 0 2 0 o r	1.30E-01 2.06E-01 2.88E-02 1.86E-01 4.85E-02 4.97E-03	0 0 0 0 0 1 0 1 0 0 0 2 0 1 1 0 2 0	0.00 598.38 666.53 1191.91 1264.91 1340.97	1.28E-01 2.48E-01 2.91E-02 2.20E-01 5.87E-02 4.81E-03	0 0 0 (99.9) 0 0 1 (99.7) 0 1 0 (99.2) 0 0 2 (99.6) 0 1 1 (98.0) 0 2 0 (96.1)	0.00 597.08 656.57 1189.24 1249.11 1308.27	1.28E-01 2.48E-01 3.91E-02 2.21E-01 6.59E-02 6.96E-03
	0 0 3	1810.32	1.24E-01	0 0 3	1780.63	1.18E-01	0 0 3 (99.5)	1776.42	1.21E-01
	1 0 0 0 1 2 0 2 1 0 3 0 0 0 4 1 0 1 0 1 3 1 1 0	1869.60 1874.03 1937.74 P 3.56E-07 4.61E-02 8.68E-03 2001.45 7.17E-04 2413.76 6.80E-02 e 2473.04 1.32E-07 2477.47 3.22E-02 2536.75 3.95E-06 e 1 0 0 0 1 2 0 2 1 0 3 0 0 0 4 1 0 1 0 1 3 1 1 0 r	1802.77 1858.45 1939.35 2022.01 2364.55 2401.14 2447.16 2469.30	6.05E-08 5.48E-02 1.00E-02 6.20E-04 4.27E-02 4.71E-09 3.15E-02 3.01E-06	1 0 0 (97.4) 0 1 2 (96.8) 0 2 1 (92.9) 0 3 0 (90.8) 0 0 4 (99.2) 1 0 1 (92.6) 0 1 3 (93.0) 1 1 0 (91.7)	1807.43 1837.10 1896.05 1956.35 2359.02 2404.36 2420.77 2455.84	1.71E-04 5.08E-02 1.04E-02 8.95E-04 4.50E-02 7.71E-04 2.32E-02 6.05E-05
	0 2 2	2541.18	8.51E-03	0 2 2	2532.88	9.77E-03	0 2 2 (89.8)	2479.75	7.06E-03
	0 3 1 0 4 0 0 0 5 1 0 2 0 1 4	2604.89 2668.60 3017.20 3076.48 3080.91	1.29E-03 9.33E-05 3.24E-02 5.40E-10 1.85E-02	0 3 1 0 4 0 0 0 5 R 2620.39 2708.44 2943.71 1 0 2 2994.68 0 1 4 3031.08 e	1.36E-03 6.67E-05 1.11E-02 2.13E-07 1.24E-02	0 3 1 (85.5) 0 4 0 (85.6) 0 0 5 (98.9) 1 0 2 (59.8) 0 1 4 (61.4)	2539.76 2618.72 2938.92 2994.23 3003.10	1.21E-03 7.83E-05 1.22E-02 2.47E-03 4.74E-03
	1 1 1 0 2 3 1 2 0 0 3 2 0 4 1	3140.19 3144.62 3203.90 3208.33 3272.04	5.42E-06 6.13E-03 1.57E-06 1.30E-03 1.73E-04	1 1 1 0 2 3 1 2 0 0 3 2 0 4 1	3067.68 3121.60 3143.74 3213.92 3306.82 v 4.96E-06 5.90E-03 1.25E-06 1.41E-03 1.56E-04 i e	1 1 1 (77.0) 0 2 3 (78.7) 1 2 0 (80.9) 0 3 2 (80.3) 0 4 1 (78.0)	3046.89 3060.10 3101.68 3119.49 3201.42	2.97E-04 2.65E-03 1.32E-05 7.27E-04 9.65E-05
	0 5 0 0 0 6 1 0 3	3335.75 3620.64 3679.92	1.12E-05 1.38E-02 4.77E-08	0 5 0 0 0 6 2 0 0	3399.15 3518.15 3548.87	5.88E-06 2.15E-03 1.99E-07 w 0 5 0 (78.8) 2 0 0 (95.3) 0 1 5 (53.5)	3276.91 3554.86 3579.13	8.12E-06 4.65E-06 9.46E-04
	0 1 5	3684.35	9.13E-03	1 0 3	3583.40	4.04E-07	1 0 3 (54.0)	3589.49	8.57E-04
	2 0 0 1 1 2 0 2 4 1 2 1 0 3 3 1 3 0	3739.20 3743.63 3748.06 3807.34 3811.77 3871.05	5.16E-08 4.23E-06 3.61E-03 2.41E-06 9.64E-04 3.87E-07	0 1 5 1 1 2 0 2 4 1 2 1 0 3 3 1 3 0	3610.24 3661.21 3705.52 3742.11 3802.64 3824.78	3.58E-03 3.63E-06 2.47E-03 2.29E-06 9.04E-04 2.90E-07	0 2 4 (46.9) 1 1 2 (48.0) 1 2 1 (53.7) 0 3 3 (57.5) 1 3 0 (70.1) 0 4 2 (72.9) O n 3629.74 3641.86 3687.48 3697.04 3746.91 3779.90 l	4.16E-04 4.49E-04 7.35E-05 2.03E-04 4.08E-07 5.49E-05
	0 4 2 0 5 1	3875.48 3939.19	1.79E-04 2.13E-05	0 4 2 0 5 1	3900.36 3997.53	1.73E-04 1.49E-05	0 5 1 (73.2)	3865.65 y	1.01E-05
					24

URL: http://mc.manuscriptcentral.com/tandf/tmph

URL: http://mc.manuscriptcentral.com/tandf/tmph Molecular Physics

the author for peer review, but cannot be