19 December 1997
Chemical Physics Letters 281 Ž1997. 161–167
On the quasidiabatic character of average natural orbitals V.M. Garcıa ´ a, M. Reguero a, R. Caballol a, J.P. Malrieu
b
a b
Departament de Quımica Fısica i Inorganica, UniÕersitat RoÕira i Virgili, Plac¸a Imperial Tarraco 1, 43005 Tarragona, Spain ´ ´ ` IRSAMC, Laboratoire de Physique Quantique, UniÕersite´ Paul Sabatier, 118, route de Narbonne, 31062 Toulouse Cedex, France Received 7 August 1997; in final form 17 September 1997
Abstract The nearly diabatic character of average natural orbitals obtained from linear combinations of density matrices associated with different states is rationalized. These orbitals may be used to build low-dimensional model spaces and effective Hamiltonians which give a nearly diabatic description of the set of interacting states. q 1997 Elsevier Science B.V.
1. Introduction The solutions of the Schrodinger equation provide ¨ adiabatic wavefunctions and potential energy surfaces. Nearly diabatic descriptions of electronic states can be useful for a qualitative understanding of wavefunctions and how they change when their molecular conformation is modified, as well as being suitable for dynamical treatments. Numerous procedures have been proposed to obtain such nearly diabatic descriptions, most of which involve some external criteria w1–8x. We show first in this paper that the eigenfunctions of mean density matrices have a nearly diabatic character. These types of eigenfunctions can as well be obtained by diagonalizing appropriate linear combinations of the density matrices associated with different states lying close in energy, and can be used to define nearly diabatic configurations that span model spaces which are suitable for building low-dimensional effective Hamiltonians. Although extensive configuration interactions are necessary to obtain quantitative potential energy surfaces, it is usually possible to find appropriate one-
electron orbitals in such a way that the excited eigenstate of interest Cm is dominated by very few configurations, in general one or two. The most rational way to obtain this representation is to diagonalize the one-electron density matrix Rm associated with Cm
Rm s ²Cm < R
Rs Ý p
Ý aqp a q ,
Ž 1.
q
where p and q belong to a set of orthonormal monoelectronic functions. Around avoided crossings, the nature of the eigenstate Cm changes rapidly with the nuclear conformation, and a configuration F 1 predominates in a certain domain D 1 of nuclear coordinates, while another configuration F 2 predominates in an adjacent domain D 2 of conformations. In these cases some natural orbitals of R1 will change when going from D 1 to D 2 . If the interaction between F 1 and F 2 is weak, the change will be rapid on the border between D 1 and D 2 and the adiabatic natural orbitals may undergo rapid changes. It has been suggested elsewhere w9x that for energy differences between two Žor several. states to be
0009-2614r97r$17.00 q 1997 Elsevier Science B.V. All rights reserved. PII S 0 0 0 9 - 2 6 1 4 Ž 9 7 . 0 1 0 7 3 - 7
V.M. Garcıa ´ et al.r Chemical Physics Letters 281 (1997) 161–167
162
consistently calculated, it is worth using a common set of molecular orbitals ŽMOs.. This leads to the use of average density matrices
Rs
Ž R1 q R 2 . 2
cisely, the valence part of the ground state may be written
C 0 s lsg sg y msu su , ,
Ž 2.
and the corresponding average natural orbitals are obtained. The starting matrices R1 and R2 may be obtained from any set of MOs. The average MOs are expected to give a better simultaneous description of C 1 and C 2 . Section 2 shows that in many cases the mean natural orbitals ŽMNO. obtained in this way are nearly diabatic. This property, that has been pointed out in interesting papers by Ruedenberg and Atchity w10,11x as a numerical result, is rationalized here. Once nearly diabatic MOs have been obtained, they may be exploited to produce nearly diabatic leading configurations that define a simple model space on which an effective Hamiltonian may be built. The diagonal elements represent the effective energies of nearly diabatic states, and the off-diagonal terms represent the couplings between them. Section 3 discusses the domain of applicability of this technique and shows that in some problems the use of differences of density matrices may be helpful to remove near degeneracies of occupation numbers and obtain well-localized nearly diabatic MOs.
2. Nearly diabatic eigenvectors of mean density matrices: a few examples This section will illustrate and rationalize the nearly diabatic character of mean natural orbitals through simple examples from the spectroscopy of diatomic systems. 2.1. The double well of the E,F 1Sgq potential curÕe of H2 The E,F 1 Sq g state is a text-book example of a double well potential curve. It is interpreted as being due to the crossing between a state generated by a Rydberg single excitation from the ground state and the valence excited state which is dominated by a doubly excited determinant of essentially ionic character in the valence bond ŽVB. sense. More pre-
l)m ,
Ž 3.
where sg and su are the valence bonding and antibonding MOs, respectively, built of 1s type atomic orbitals. The valence excited state belonging to the same symmetry would be
Cv) s msg sg q lsu su ,
Ž 4.
and the Rydberg configuration
C R) s
Ž sg sgR q sgR sg . '2
,
Ž 5.
where sgR is the diffuse Rydberg MO of Sq g symmetry. The minimum associated with C R) is found at shorter distances than the minimum associated with Cv) . The diagonal elements of the density matrices for these states, expressed in the basis of the three MOs, sg , su and sgR , can be given schematically as
sg
sgR
su
R0 :
2 l2
0
2 m2
Rv : RR :
2 m2 1
0 1
2 l2 0
Ž 6. .
The interaction between C R) and Cv) is expected to mix them, giving
C 1 s aCv) q bC R) ,
Ž 7.
C 2 s ybCv) q aC R) .
Ž 8.
If the basis set, and consequently the Cl are larger, minor components appear in C 1 and C 2 . The associated density matrices will be approximately:
R1 f
sg
sgR
su
2 m2a 2 q b 2
2 mabr'2
0 0 2 a 2l2
2 mabr'2 0
2
b 0
Ž 9.
V.M. Garcıa ´ et al.r Chemical Physics Letters 281 (1997) 161–167
163
and
R2 f
sg
sgR
su
2 m2b 2 q a 2
y2 mab r'2
0 0 . 2 b 2l2
y2 mab r'2
2
a 0
0
Ž 10 .
If a mean density matrix from R1 and R2 is calculated, the off diagonal terms ² sg < R 1 < sgR : and ² sg < R 2 < sgR : have nearly opposite values, 2 mabr '2 and y2 mabr '2 respectively, so
¦s < R < s ;f 0 , g
Ž 11 .
gR
and the following occupation numbers are expected
R:
sg
sgR
1 q 2 m2
1
2
2
su l2
.
Ž 12 .
If m2 is not negligible Žas is the case for the lengthened H–H bond., the difference between the diagonal matrix elements associated to the valence sg and the Rydberg sgR orbitals is large, 2 m2 , and both orbitals will not mix in the natural orbitals. However, when m2 is small a mixture of both orbitals is expected. In this case, the physical identity of the mean natural orbitals is perfectly preserved when the average of three density matrices is considered
Rs
Ž R 0 q R1 q R 2 . 3
,
Ž 13 .
H eff s
X X Ý Ž PsCi . :Ei² Ž PsCi .
,
Ž 16 .
is0,2
where Ei is the exact eigenenergy, and Ci the associated FCI eigenvector
Ž 17 .
and
Ž 2 R 0 q R1 q R 2 .
,
4
Ž 14 .
since the interference effect cancels the off-diagonal terms of the mean density matrix and the diagonal terms from Eq. Ž13. are different
R:
in a small basis set w7s4p2dr5s4p2dx. The mean density matrix of the last combination Žequation 14. was diagonalized. Whatever the interatomic distance, the first sg natural MO has a valence character and the second is a Rydberg orbital. The three configurations sg sg , su su and Ž sg sgR q sgR sg . define a model space of projector Ps on which the exact effective Hamiltonian can be built according to the definition of des Cloizeaux w12x
H
or
Rs
Fig. 1. Potential energy curves for the adiabatic Ž – – – . and . excited 1 Sq diabatic Ž g states of H 2 .
sg
sgR
su
3
1
2
3
3
3
.
Ž 15 .
Hence the MNOs of valence and Rydberg orbitals will not mix. The three lowest 1 Sq g states of H 2 were calculated by a Complete Interaction Configuration ŽFCI.
PsCiX 4 s Sy1 r2 PsCi 4 .
Ž 18 .
A more conventional two-state modelization may be obtained by using the C R) and Cv) configurations defined in Eqs. Ž4. and Ž5., respectively. The diabatic and adiabatic potential energy curves of the two lowest excited 1 Sq g states are shown in ˚ as it was in the Fig. 1. The crossing is close to 1.7 A, exact treatment of Kolos and Wolniewicz w13x. 2.2. Neutralr ionic crossing in the ground state of polar diatomic molecules: LiH and LiF The ground state of polar molecules such as LiH and LiF is dominated by an ionic VB component
V.M. Garcıa ´ et al.r Chemical Physics Letters 281 (1997) 161–167
164
LiqHy, LiqFy around their equilibrium geometries. They dissociate into neutral atoms, so the asymptotic content of the singlet ground state is neutral in Valence Bond ŽVB. terms. A curve crossing occurs between the ionic and neutral diabatic states. For LiF the ionicrneutral curve crossing takes place at a very large interatomic distance Ž10 to 12 bohr., the crossing is very weakly avoided and the change in the function is sudden. For LiH the diabatic curve crossing takes place at a much shorter interatomic distance and the crossing is more strongly avoided. In these LiX ŽX s H, F. problems, two atomic orbitals play the main roles, namely the 2s orbital of Li and one of the s X orbitals of X, the 1s orbital of H or the 2p z orbital of F. From these atomic orbitals two main configurations can be defined, one of which has an ionic VB character, F i , and the other has a neutral VB character, Fn
Fn s
'2
,
Ž 19 .
F i s
Ž 20 .
The diagonal terms of the density matrices associated to the neutral open shell configuration, Fn , and the ionic closed shell configuration, F i , are Rn : Ri :
sX 1 2
2s Ž Li . 1 0
Ž 21 . .
Whatever the interatomic distance, the two lowest 1 q S states are dominated by orthogonal combinations of Fn and F i
C 1 s lFn q mF i ,
Ž 22 . Ž 23 .
C 2 s ymFn q lF i , for which the density matrices are respectively
sX
2s
R1 f
l2 2 mlr'2 2s
R2 f
m
2
y2 mlr'2
2 mlr'2 2
l q2m
2
Ž 24 .
,
sX y2 mlr'2
m 2 q 2 l2
.
Ž 25 .
So, when the mean density matrix R s Ž R1 q R2 .r2
is considered, the occupation numbers are close to 1r2 for 2s and 3r2 for s X . In contrast to this large difference between the diagonal terms of R , it should be pointed out that the cancellation that takes place in the off-diagonal matrix elements of R ²2s < R 1 < s X : s 2 lmr'2 ,
Ž 26 .
²2s < R 2 < s X : s y2 lmr'2 ,
Ž 27 .
²2s < R < s X : s 0 ,
Ž 28 .
will still be approximately valid when working with non minimal basis sets and large CI expansions, i.e. when
C 1 s lFn q mF i q remainder , C 2 s ymFn q lF i q remainder ,
Ž 29 . Ž 30 .
where ‘‘remainder’’ represents the components of the external space but the projections of the two wave functions in the model space spanned by F i and Fn still remain nearly orthogonal, so that
l f lX , mf mX and X X ²2s < R < s X : f lm y lm f0 .
Ž 31 .
As a consequence, the MNOs have still a local character. A FCI calculation was performed on the two lowest 1 Sq states of LiH. The basis set for lithium was a DZP w14x to which one d polarization function Ž zd s 0.2. was added. The basis set for hydrogen was a w5s2pr3s2px w15x. The mean density matrix gave MNOs with the expected character. From these natural orbitals, whose occupation numbers were close to 3r2 and 1r2, the two main configurations Fn and F i were defined and the des Cloizeaux effective Hamiltonian was built as before on the model space given by these two last configurations. The diabatic and adiabatic potential energy curves are shown in Fig. 2. The curve crossing is close to 6 bohr, as it was in the treatment of Boutalib et al. w16x where an explicit diabatization was performed. The value of the off-diagonal matrix element of the effective Hamiltonian against the interatomic distance increases when the interatomic distance decreases, as expected for an effective charge-transfer ²1s H < F <2s Li : integral Ž F is the Fock operator., which behaves as the overlap between the two non-orthogonal atomic orbitals.
V.M. Garcıa ´ et al.r Chemical Physics Letters 281 (1997) 161–167
165
clearly in the zoom. At this distance, the neutralrionic coupling is very weak, 0.047 eV. The diabatic coupling versus the interatomic distance shows a very similar behaviour to the one reported by Bauschlisher et al. w14x The charge-transfer interaction presents a regular exponential decrease at large interatomic distances and its maximum is close to 3.0 bohr.
3. Using other combinations of density matrices
Fig. 2. Potential energy curves for the adiabatic Ž – – – . and . lowest 1 Sq states of LiH. The energy is diabatic Ž shifted by y8.0 au.
For the LiF calculation, the DZP basis set used by Bauschlisher et al. w14x in their FCI calculation was used. The two lowest 1 Sq g states were calculated with the iterative DDCI method, considering two active orbitals for the starting CAS. The adiabatic curve for the ground state, already published w9x, was in perfect agreement with the FCI results. The des Cloizeaux effective Hamiltonian was built on a model space given by the Fn and F i configurations. Fig. 3 shows the adiabatic and diabatic curves in the range of 5 to 14 bohr. The crossing between the diabatic curves at 11.8 bohr can be seen more
It should not be forgotten that the eigenvectors of mean density matrices are not always nearly diabatic. Let us consider for instance the interaction between the two configurations F 1 and F 2 that are singly excited from the closed shell ground state F 0 and which are the leading configurations of two states C 1 and C 2 . Defining 1
q q < : '2 Ž a i a p " a i a p F 0 . ,
q q < : '2 Ž a j a q " a j a q F 0 . ,
1
Ž 32 . Ž 33 .
where p and q are occupied and i and j virtual MOs, and
C 1 s lF 1 q mF 2 ,
Ž 34 .
C 2 s ymF 1 q lF 2 ,
Ž 35 .
the diagonal elements of the density matrices relative to F 1 and F 2 are R1 : R2 :
p 1 2
q 2 1
i 1 0
j 0 1
Ž 36 . .
So those of C 1 and C 2 will be p
Fig. 3. Potential energy curves for the adiabatic Ž – – – . and . lowest 1 Sq states of LiF. The energy is diabatic Ž shifted by y106.84 au.
2
RC1 :
l
RC 2 :
2 m2
q 2m
l2
i 2
j
2
m2
m2
l2
l
Ž 37 . .
Then the diagonal matrix elements of the mean density matrix R s R1 q R2 are degenerate for both the occupied Ž R p p s R q q s 2 m 2 q l2 . and the virtual MOs Ž R i i s R j j s l2 q m2 ., with reference to F 0 . Since F 1 and F 2 differ by two orbitals, the off-di-
V.M. Garcıa ´ et al.r Chemical Physics Letters 281 (1997) 161–167
166
agonal elements of the mean density matrix are zero. Hence the mean natural orbitals remain degenerate and have no reason to be localized. As a physical example, let us consider two identical molecules A and B in weak interaction, such as the H 2 PPP H 2 system. Let us suppose that the bond lengths are different, i.e. r y E for molecule A and r q E for molecule B. The excitation energy will be smaller for the lengthened molecule than for the shortened one. A natural diabatic description would use orbitals pA and i A localized on A, and q B and j B localized on B. In the lowest eigenstate C 1
C 1 s lFpA i A q mFq B j B ,
Ž 38 .
The excitation is essentially localized on B for E ) 0 Ž m 4 l., while it moves to A when E - 0 Ž l 4 m .. However, in such a problem the mean natural orbitals are not localized. Considering now another combination of density matrices, namely the difference Rs R 1 y R 2 , it may be easily seen that its diagonal elements are p R:
2
m yl
q 2
2
l ym
i 2
2
l ym
j 2
2
m y l2
.
Ž 39 . If the diagonalization is performed by restricting by blocks of the occupied MOs p and q on one side and of the virtual MOs i and j on the other side, the corresponding eigenvectors will not rotate p and q
nor i and j and will therefore localize on the separate subsystems, except when l s m , i.e. for E s 0. This model problem ŽH 2 . 2 was calculated using a DZP basis set w17x, with the following collinear geometry
and several values of the distortion E between E s ˚ The states of interest are built "0.01 and 0.1 A. from local 1 Sq states. Even for the lowest value of E u the natural MOs of R are perfectly localized and they lead to the diabatic curves plotted in Fig. 4 and to an almost constant diabatic coupling which has an amplitude of 0.0047 au. The adiabatic curves are also shown in the figure. The off-diagonal element follows a perfect quadratic dependence on E . This E-dependence of the coupling can be easily rationalized by taking into account that the interaction between the two localized singlet states is a dipoletransition interaction ™ ™ pA i A q B j B
m
²Fp i < H
s
Ž ryE . Ž rqE . R3
m
R3 r2yE 2 s
R3
.
Ž 40 .
Fig. 4. Potential energy curves for the adiabatic Ž – – – . and . lowest excited 1 Sq states of H 2 –H 2 . The diabatic Ž energy is shifted by y1.8 au.
The strategy presented on this model problem can be generalized to systems involving inactive core and inactive virtual MOs. After a preliminary definition of the three sets of MOs Žcore, active, virtual. by diagonalizing the mean density matrix R , the MOs can be rotated according to the following steps: Ži. C 1 and C 2 and the corresponding density matrices R1 and R2 are calculated; Žii. the block of the R s R1 y R2 matrix corresponding to active orbitals is separeted into two subblocks of occupied and virtual MOs with respect to F 0 ; Žiii. both blocks are separately diagonalized. The resulting MOs are well localized as in the model problem.
V.M. Garcıa ´ et al.r Chemical Physics Letters 281 (1997) 161–167
4. Conclusions In this Letter, Atchity and Ruedenberg’s w10x observation that the eigenvectors of mean density matrices have a quasidiabatic character has been rationalized. Several examples ŽH 2 1 Sq g excited states, ionicrneutral 1 Sq curve crossing in LiH and LiF. confirm this property and explain its occurrence. A counter example has been found for the interaction between two states that are combinations of two singly excited configurations. In such cases diabatic MOs may be obtained by carrying out a difference of density matrices. It may be recommended in general: Ži. To build and diagonalize the mean density matrix relative to C 1 and C 2 . The mean natural orbitals can be classified into a few groups: those with occupation numbers close to 2 are inactive occupied MOs that do not participate in the electronic change between C 1 and C 2 . The inactive virtual orbitals are those with small, nearly zero, occupation numbers. The active holes have occupation numbers larger than 1 and significantly smaller than 2. The active particles have non negligible occupation numbers lower than 1. In many cases these mean natural orbitals already have a diabatic character and may be used to construct the small-dimensional model space and effective Hamiltonian. Žii. If the active natural MOs do not have a diabatic character, other linear combinations of density matrices can be built. The difference of the density matrices relative to C 1 and C 2 can be calculated, and the block of the active holes on one side, and the active particles on the other side diagonalized. The corresponding eigenvectors have been shown to become nearly diabatic in the difficult Žand frequent. case of interaction between two singly excited states which differ by two orbitals. The advantage of these diabatization procedures is that they do not call for extra information, such as projections of valence bond guess functions, orbital following from the asymptotes etc. The method only manipulates matrices of the adiabatic states, which are intrinsic. A previous work showed the computational efficiency of considering mean density matrices to im-
167
prove the orbitals in correlated calculations of eigenenergies. It is interesting to see that the corresponding eigenvectors have some qualitative properties, that can make a physical interpretation of the processes occurring on the potential energy surfaces possible.
Acknowledgements This work has been supported by DGICYT of the Ministerio de Educacion ´ y Ciencia of Spain Žproject No. PB95-0847-CO2-02., by CIRIT of the Generalitat de Catalunya Žgrant SGR95-426. and by the European Commission through the TMR network contract ERB FMRX-CT96-0079 ŽQuantum Chemistry of Excited States.. VMG thanks the CIRIT of the Generalitat de Catalunya for a predoctoral fellowship. The Laboratoire de Physique Quantique is Unite´ Mixte de Recherche 5626 du CNRS.
References w1x w2x w3x w4x w5x w6x w7x w8x w9x w10x w11x w12x w13x w14x w15x w16x w17x
F.T. Smith, Phys. Rev. 179 Ž1969. 111. V. Sidis, H. Lefebvre Brion, J. Phys. B 4 Ž1971. 1040. J.C. Tully, J. Chem. Phys. 50 Ž1973. 5122. C. Kubach, V. Sidis, Phys. Rev. A 14 Ž1976. 152. R. Cimiraglia, M. Persico, Mol. Phys. 70 Ž1979. 1774. C. Petrogolo, R.J. Buenker, S.D. Peyerimhoff, J. Chem. Phys. 78 Ž1983. 7284. F. Spiegelmenn, J.P. Malrieu, J. Phys. B: At Mol. Phys. 17 Ž1984. 1259. F.X. Gadea, M. Pelissier, J. Chem. Phys. 93 Ž1990. 545. ´ V.M. Garcıa, ´ O. Castell, R. Caballol, J.P. Malrieu, Chem. Phys. Letters 238 Ž1995. 229. G.J. Atchity, K. Ruedenberg, J. Chem. Phys. 99 Ž1993. 3790. K. Ruedenberg, G.J. Atchity, J. Chem. Phys. 99 Ž1993. 3799. J. des Cloizeaux, Nucl. Phys. 20 Ž1960. 321. W. Kolos, L. Wolniewicz, J. Chem. Phys. 49 Ž1968. 404. C.W. Bauschlicher Jr., S.R. Langhoff, J. Chem. Phys. 89 Ž1988. 4246. G. Jeung, J.P. Daudey, J.P. Malrieu, Chem. Phys. Letters 98 Ž1983. 433. A. Boutalib, F.X. Gadea, J. Chem. Phys. 97 Ž1992. 1144. C.W. Bauschlicher, P.R. Taylor, J. Chem. Phys. 85 Ž1986. 6510.