- Email: [email protected]
An efficient relativistic multireference configuration interaction method
CHEMICAL PHYSICS LETTERS
Volume 215, number 1,2,3
26 November 1993
An efficient relativistic multireference configuration interaction method M. Vijayakumar, S. Roszak Departmentof Received 1
The tion
and IL
and Biochemistry, 1993; in
form 20
State University,
AZ 85287-1604,
1993
of the confiiation interaction and energy procedures. The for the of energies properties. The correction accounting unlinked quadrupole
Introduction Relativistic effects are known to play an important role on chemical and spectroscopic properties of molecules containing very heavy atoms [ 11. The relativistic effective core potentials (RECP) employed in conjunction with the multiconfiguration self-consistent field/multireference configuration interaction (CASSCF/MRCI ) method have provided very good description of molecules containing heavy atoms [ 21. Pitzer and co-workers [ 3,4] have provided reliable RECPs derived from the relativistic Dirac-Fock calculations. The spin-orbit operator itself can be expressed in terms of RECP as shown both by Hafner and Schwartz [ 5 ] and Ermler et al. [ 6 1. The radial part of the spin-orbit was molded in a Gaussian analytical form convenient for use in molecular structure calculations. Pitzer and Winter [ 7 ] have developed the codes to evaluate the spin-orbit integrals in a Gaussian basis set. The spin-orbit integrals thus obtained can be introduced either perturbatively or variationally. One of the authors [ 81 (KB) has developed an approach that takes advantage of recent advances in the nonrelativistic MCSCF and direct configuration methods and incorporates ’ Permanent address: Institute of Physical and Theoretical Chemistry, Technical University, Wyb. Wyspiadskiego 27, 50370 Wrociaw, Poland. 0009-2614/93/S 06.00 0 1993 Elsevier Science Publishers B.V.
has been by employing tests for and Pb2 the usefulness developed method includes both is called relativistic multireference
confgurathe proposed and the
all the relativistic effects, including the spin-orbit term. In the RCI calculations, the spin-orbit integrals are transformed in the natural orbital basis obtained from the CASSCF/CI methods and then added to the one-electron Hamiltonian matrix at the RCI step. Consequently, the spin-orbit coupling term is included variationally rather than perturbatively. The RCI code developed up to now is restricted to calculations that include up to 25000 configurations. The configurations included in the CI process contribute in a different way to the total energy of the system. We propose a selection procedure based on the perturbation energy expression to calculate an energy contribution of any particular configuration when added to the reference configuration set. Based on a chosen threshold, the configurations are selected in the new relativistic multireference configuration interaction (RMRCI) calculation. Contributions of the unselected configurations are calculated, assuming that there are no direct couplings between contributions. The contributions of unselected configurations depend linearly on the assumed threshold and an extrapolation procedure has been developed to account for them. In addition the Davidson correction for unlinked quadrupole clusters is also included in the RMRCI. It is shown that the results of the RMRCI method are quite comAll rights reserved.
87
Volume 2 15, number 1,2,3
CHEMICAL PHYSICS LETTERS
parable to fully variational RCI expansions of much larger sizes.
26 November 1993
figurations. The efficiency of the program can be improved by applying the configuration selection procedure described in section 3.
2. General CASSCF/MRSDCI/RCI method 3. Method of configuration selection The calculations follow three main steps. Molecular orbitals for the MRSDCI calculations are generated by the complete active space MCSCF method (abbreviated as CASSCF). In this method a set of the most important electrons for chemical bonding (active electrons) are distributed in all possible ways among a chosen set of orbitals referred to as the internal or active set of orbitals. The CASSCF method provides a starting set of orbitals for the inclusion of higher-order correlation effects. Higher-order electron correlation effects not included in the CASSCF, are taken into account using the configuration interaction method. The CI calculations are multireference singles + doubles CI (MRSDCI). The MRSDCI calculations include a subset of reference configurations in the SOCI, determined by the most important configurations in the CASSCF. Subsequently, single and double excitations from each of the reference configurations are allowed in the MRSDCI. The RCI calculations are carried out following the CASSCF/MRSDCI methods. The CASSCF/ MRSDCI calculations described above include electron correlation corrections as accurately as possible, but neglect the spin-orbit coupling term entirely, which is far from negligible for molecules containing heavy atoms. In the RCI method the spin-orbit integrals are transformed into the natural orbital basis obtained from the CASSCF/MRSDCI methods and then added to the one-electron Hamiltonian matrix at the RCI step. The spin-orbit coupling term is thus included variationally. The RCI calculations include all configurations that mix in the presence of the spin-orbit coupling term as reference conflgurations. The inclusion of spin-orbit coupling term in the Hamiltonian changes the normal point group symmetry of a molecule into the spin double group. Thus, all configurations that have the same symmetry in a spin double group could mix in the RCI. The RCI code developed up to now does not have the capabilities of the direct CI procedures and is restricted to calculations that include up to 25000 con88
The goal of the configuration selection procedure is to select configurations which contribute significantly to the energy when included variationally in the CI space. One convenient way of doing this is to select a set of most important configurations (main reference space) and to check the energy lowering when a tested configuration is added to the reference set. The energy lowering generated by a tested configuration #’ by interacting with configuration em which is a member of the main reference set can be calculated by the perturbation expression [ 9 ]
The energy lowering implied by eq. ( 1) was first used by Whitten and Hackmeyer [ 10,111 as a contiguration selection technique in non-relativistic CI computations. Different selection techniques based on eq. ( 1) have been developed and used successfully by Bender and Davidson [ 12 1, Shavitt [ 131 and most extensively by Buenker and Peyerimhoff [ 14,151. In this work, configurations are tested against every configuration which is a member of the main reference set using eq. ( 1). The configuration is being selected when AE exceeds the assumed threshold for at least one main configuration.
4. Accounting of unselected con&u&ions Expression ( 1) may also be used for the estimation of contribution of unselected configurations,
(2) where the contribution of the unselected contiguration is calculated against the zeroth-order CI wavefunction obtained within the reference space. Eq. (2) does not account for couplings between configurations which are outside the reference space. The ex-
Volume 215, number 1,2,3
CHEMICAL. PHYSICS LETTERS
act energy (threshold equal to zero) can be expressed as E(O)=E(T)+I(T)K(T).
(3)
For a particular selection threshold T, the function K(T) has to be scaled to account for internal couplings and possible relaxation processes in the CI performed with no threshold. The proportionality factor 1, in general, may be a function of the threshold. Bunker and Peyerimhoff [ 16 ] proposed that A is independent of the threshold for non-relativistic CI computations. However our test calculations indicate a linear dependence on the threshold T namely, A(T)=A+BT.
(4)
The performed RMRCI calculations are of the second-order type. The difference between the RMRCI and the full configuration interaction calculations may be estimated by the unlinked quadrupole cluster correction formula proposed by Davidson et al. [ 17,18 1, extended later for a multireference case by Peyerimhoff and Buenker [ 191. The correction expression is given as
where E,, is a reference secular equation and EC1 the second-order energy. Eq. (5) has been shown to be a good correction for the lack of size extensivity of the energy with single and double excitations [ 20,2 11. Since the full relativistic CI calculations are not yet available for comparison, the correction should be considered cautiously as an estimated measure of the missing contribution to the total energy.
5. Computational details To test the present method, we have performed calculations for some molecules with heavy elements such as Pbz and Biz. These two dimers were considered before using the RCI method [ 22,23 1. The natural orbitals for various points in the potential energy curves have been generated from the first-order configuration interaction (FOCI) method following the complete-active-space multicontiguration self-
26 November 1993
consistent field (CAS-MCSCF) [ 8 ] calculations. All CASSCF/FOCI calculations are performed using the modified version [ 8 ] of the ALCHEMY II #l code. In all calculations, we have employed the relativistic effective core potentials (RECPs) [ 241 for both the Pb and Bi atoms. Thus the valence space contains 6s26p2 and 6s26p3 shells for the Pb and Bi atoms, respectively. The ( 3s3p 1d ) valence Gaussian basis set has been used. The spin-orbit integrals were computed using Pitzer’s modified ARGOS codes [ 7 1. In the RMRCI calculations, the low-lying electronics states which yield an n state of the desired symmetry have been included as reference contigurations. Along with the reference configurations, the singly- and doubly-excited configurations from these reference configurations have also been included in the RMRCI calculations. For the Pb2 dimer, the 0: state calculations included the reference configurations from the 3Z; ground state and the excited states such as ‘Ez, 311g,‘Zz (II) and ‘II, (II) states. For the 0: state calculations of Bi2, the reference configurations which we included are the ‘El ground state and the 0: states arising from the excited states such as 31’Ig,‘E:, ‘II, (II), ‘IS; and ‘El (II). The 1U state calculations of Biz are performed by using the 3Z,*( l,), 3AU( 1,) and ‘Zy ( 1,) states as the reference configurations,
6. The effect of configuration selection procedure on energies and properties As seen from tables 1 and 2 the dramatic reduction of selected configurations introduces very small absolute error to the total energy. The threshold of 1 phartree gives errors lying on the border of precision of calculations, despite the reduction of number of configurations to one third of the original number of CSFs. The largest threshold ( 100 uhartree) reduces the number of configurations almost ten times, but the relative error to the total energy is less than 0.05% in all cases studied here. The potential energy surfaces (figs. 1 and 2) as well as geometries (see table 3)) are also well reproduced, even for the highest threshold employed. The impact of c1 The major authors of ALCHEMY II are B. Lengsfield, B. Liu and M. Yoshimine.
89
Volume 2 15, number 1,2,3
CHEMICAL PHYSICS LETTERS
26 November 1993
Table 1 Number of selected configurations, total energy and absolute error for different selection thresholds in hartree atomic units for the 0: ground and excited states of Bis (2.70 A) Selection threshold (Wree)
Number of selected configurations
0: energy
AEX 10-s
energy
AEX lo-2
0 1.0 20.0 40.0 60.0 80.0 100.0
18113 7393 3843 3059 2579 2291 2058
-10.71161 -10.71154 - 10.71031 - 10.70958 - 10.70871 - 10.70829 - 10.70748
0.0073 0.1299 0.2029 0.2899 0.3319 0.4129
-
0.0085 0.1015 0.1710 0.2440 0.2892 0.3561
0: (II)
K(T)XlO_2
10.63308 10.63300 10.63206 10.63137 10.63063 10.63019 10.62951
-0.0036 -0.1528 -0.2291 -0.3126 -0.3603 -0.4416
Table 2 Number of selected configurations, total energy and absolute error (in atomic units) for different selection thresholds for 1, states of Bi2(2.70 A) Selection threshold ( phartree)
Number of selected configurations
energy
AEX lo-2
energy
AEX lo-2
energy
AEX lo-2
0 1.0 20.0 40.0 60.0 80.0 100.0
13780 4377 2405 1965 1733 1581 1455
-
0.0006 0.1088 0.1783 0.2506 0.2835 0.3600
-
0.0006 0.1066 0.1688 0.2410 0.2741 0.3474
-10.61180 -10.61172 -10.61068 - 10.60999 - 10.60926 - 10.60892 - 10.60816
0.0008 0.115 0.1813 0.2540 0.2882 0.3641
1.
l”(II)
10.66207 10.66201 10.66098 10.66028 10.65956 10.65923 10.65847
K(T)XlO_2
l”(III)
10.64028 10.64022 10.63921 10.63859 10.63787 10.63754 10.63681
-0.0024 -0.1121 -0.1786 -0.2365 -0.2635 -0.3243
-0.800,
-0.*20\ 2.5
2.75
3 Distance
3.25
3.5
&)
Fig. 1, Potential energy curves for the ground state (0: ) of Bis for different selection thresholds, (-) T=O phartree, (---) T=20 phartree, (---) T= 100 bhartree.
90
2.7
2.9
2.9
3 Distance(
3.1
3.2
3.3
‘A)
Fii. 2. Potential energy curves for the ground state of Pbs for different selection thresholds, (-) T=O phartree, (---) T=20 phartree, ( ---) T= 100 phartree.
CHEMICAL PHYSICS LETTERS
Volume 215, number 1,2,3
Table 3 Ground state bond distances (A) of Bia and Pb2 calculated with and without spin-orbit coupling for different selection thresholds T (in phartree ) Molecule
Spin-orbit
T=O
T=20
T=lOO
Bi2 Bi2 Pbt Pbs
Yes no Yes no
2.1564 2.7614 2.9443 2.8855
2.7683 2.7734 2.9425 2.8902
2.7586 2.7746 2.9346 2.9115
Table 4 The energy lowering (in eV) due to spin-orbit coupling for the ground state of Biz and Pbs for different thresholds T (in @hattree) Molecule
TsO
T=20
T= 100
Bi2 Pbs
0.88 1.45
0.89 1.46
0.91 1.47
the spin-orbit coupling on the energy change is also not significantly changed, even for a selection threshold of 100 uhartree (table 4). The reproduction of the wavefunction is critical for a proper description of several electronic properties. The main reference configuration sets which have been chosen (73 reference configurations for Bi2 and 56 reference configurations for Pbz) were reasonably extended, leading to the total weight (sum of the squares of the coefficients of reference configurations) as 0.98 1 for Bil and 0.960 for Pbz. The wavefunction structure is almost unchanged for different assumed thresholds. For instance, the wavefunction of the O,+ state of Biz is composed of ‘E: (74.2%) and ‘II, (13.8W). When T is equal to 20 uhartree (and 100 uhartree), the wavefunction of the 0: state is made of ‘C: (74.1%) and 311p (13.91%) (‘ZZ - 74.2% and ‘II, - 14.0%). The wavefunctions of the 1, state of Bi2, with T=O, 20 and 100phartree, are composed of 3G (80.8%) and 3A,, (14.4%), ‘E; (80.9W) and ‘A,, (14.4%) and 3Ey (81.1%) and 3A, ( 14.3%), respectively.
7. Energy extrapolation and the Davidson correction
results are known and J(T) function is calculated (fig. 3). The function is linear for a wide range of threshold values. Error bars for nonlinearity are obvious, since configurations belong to discribed sets rather than being contiuously distributed according to the energy thresholds. The function K(T) reproduces closely the energy difference E( 0) -E(T) (tables 1 and 2). Using eqs. (3) and (4) the values for at least three different thresholds are needed to determine the A and B coefficients of L(T) function. The value of the correction A( T)K(T) is added to E(T) for calculations with the smallest energy threshold studied. The averaged value of 1 ( T= 20) using energy thresholds in the range between 20 and 100 uhartree (fig. 3 ) resulted as 0.8 1 for the 0: state and 0.94 for the 1, state. These values are close to the exact I of 0.85 (0: ) and 0.99 (1,). Our linear form for 1( T)is more general than the constant value of 1 (independent of T) proposed by Buenker and Peyerimhoff [ 16 1. None of the cases studied in this work could be reasonably described by a parameter which is independent of T. The full configuration interaction energy was estimated by a standard expression (5). The results for a number of excited states of Biq are presented in table 5. Due to an extended reference space (Cc: for a reference space is larger than 0.98 for all states) the correction is very small.
0.60
I 0
The total MRSDCI energy may be calculated from eq. (3), assuming that the proportionality factor A( T) is known. In the Bi2 case, exact (T= 0)RMRCI
26 November 1993
20
1 40
60
I 60
100
T WI
Fig. 3. The A(T) scaling factor values for different electronic states of Bis as a function of the selection threshold.
91
Volume 215, number 1,2,3
CHEMICAL PHYSICS LETTERS
Table 5 Full configuration interaction energy estimation Q-corrected energy for different states of Bir Electronic state
AEX 10”
Q-corrected energy
0:
- 0.408 -0.475 -0.158 -0.166 -0.198
-10.71191 - 10.63355 - 10.66222 - 10.64045 - 10.61200
0: (II) 1” l”(II) l”(III)
8. Conclusions The efficiency of the relativistic configuration interaction (RCI) method has been improved by employing the individualized configuration selection procedure. The error introduced by the configuration selection procedure is well controlled by the calculation of the contribution of unselected configurations to the total energy. The performed tests of energy as well as properties indicate the usefulness of the proposed approach. The molecular wavefunction is also well reproduced for the electronic states considered here. The calculated contribution of unselected configurations is used for the estimation of full second-order CI energy. The linear dependence of the scaling factor A on the threshold Tallows the estimation of full SOCI energies for at least three different thresholds. The averaged A ( TMIN)value leads to a more precise value of the full SOCI energy. The contribution of unlinked quadruple clusters to the total energy (Q correction) is also calculated using the Davidson correction.
Acknowledgement This research was supported in part by US Deof Energy under grant DEFGpartment 0286ER13558.
92
26 November 1993
References [ 1] IL Balasubramanian and KS. Pitzer, in: Ab initio methods in quantum chemistry, ed. K-P. Lawley (Wiley, New York, 1987). [2] K. Balasubramanian, J. Phys. Chem. 93 (1989) 6585. [3] Y.S. Lee, W.C. Ermler and K.S. Pitzer, J. Chem. Phys. 67 (1977) 5861. [4] P.A. Christiansen, Y.S. Lee and KS. Pitzer, J. Chem. Phys. 71 (1979) 445. [ 51 P. Hafner and W.H.E. Schwartz, Chem. Phys. Letters 65 (1979) 537. [6] WC. Ermler, Y.S. Lee, P.A. Christiansen and K.S. Pitzer, Chem. Phys. Letters 8 1 ( 198 1) 70. [7] R.M. Pitzer and N.W. Winter, J. Phys. Chem. 92 (1988) 3061. [8] K. Balasubramanian, J. Chem. Phys. 89 (1988) 5731. [9] J.C. Slater, Quantum theory of atomic structure, Vol. 1 (McGraw Hill, New York, 1960). [ lo] J.L. Whitten and M. Hackmeyer, J. Chem. Phys. 51 (1969) 5584. [ 111 M. Hackmeyer and J.L. Whitten, J. Chem. Phys. 54 (1971) 3739. [ 121 C.F. Benderand E.R. Davidson, Phys. Rev. 183 (1929) 23. [ 13 ] I. Shavitt, in: Modem theoretical chemistry, VoL 3. Methods of electronic structure theory, ed. H.F. Schaefer (Plenum Press, New York, 1977) p. 189. [ 141 R.J. Buenker and S.D. Peyerimhoff, Theoret. Chim. Acta 35 (1974) 33. [ 151 R.J. Buenker and S.D. Peyerimhoff, in: New horizons of quantum chemistry, eds. P.O. Lawdin and B. Pullman (Reidel, Dordrecht, 1983) p. 183. [ 161 R.J. Buenker and S.D. Peyerimhoff, Theoret. Chim. Acta 39 (1975) 217. [ 171 E.R. Davidson, The world of quantum chemistry, eds. R. Daudel and B. Pullman (Reidel, Dordrecht, 1973) p. 17. [ 181 S.R. Langhoff and E.R. Davidson, Intern. J. Quantum. Chem. 8 (1974) 61. [ 191 SD. Peyerimhoff and R. Buenker, Chem. Phys. 42 (1979) 167. [20] J. Paldus, P.E. Warmer, F. Vissar and A. Van der Avoird, Proceedings of the 5th Seminar on Computational Methods in Quantum Chemistry, Groningen, The Netherlands (1983) p. 13. [21] D.B. Knowles, J.R. Alvarez-Collado, G. Hirsch and R.J. Buenker, 3. Chem. Phys. 92 (1990) 585. [22] K. Balasubramanian, Chem. Rev. 90 (1990) 93. [23] K. Balasubramanian and D.-W. Liao, J. Chem. Phys. 95 (1991) 3064. [ 241 R.B. Ross, J.M. Powers, T. Atashroo, W.C. Ermler, L.A. LaJohn and P.A. Christiansen, J. Chem. Phys. 93 (1990) 6654.
Volume 215, number 1,2,3
26 November 1993
An efficient relativistic multireference configuration interaction method M. Vijayakumar, S. Roszak Departmentof Received 1
The tion
and IL
and Biochemistry, 1993; in
form 20
State University,
AZ 85287-1604,
1993
of the confiiation interaction and energy procedures. The for the of energies properties. The correction accounting unlinked quadrupole
Introduction Relativistic effects are known to play an important role on chemical and spectroscopic properties of molecules containing very heavy atoms [ 11. The relativistic effective core potentials (RECP) employed in conjunction with the multiconfiguration self-consistent field/multireference configuration interaction (CASSCF/MRCI ) method have provided very good description of molecules containing heavy atoms [ 21. Pitzer and co-workers [ 3,4] have provided reliable RECPs derived from the relativistic Dirac-Fock calculations. The spin-orbit operator itself can be expressed in terms of RECP as shown both by Hafner and Schwartz [ 5 ] and Ermler et al. [ 6 1. The radial part of the spin-orbit was molded in a Gaussian analytical form convenient for use in molecular structure calculations. Pitzer and Winter [ 7 ] have developed the codes to evaluate the spin-orbit integrals in a Gaussian basis set. The spin-orbit integrals thus obtained can be introduced either perturbatively or variationally. One of the authors [ 81 (KB) has developed an approach that takes advantage of recent advances in the nonrelativistic MCSCF and direct configuration methods and incorporates ’ Permanent address: Institute of Physical and Theoretical Chemistry, Technical University, Wyb. Wyspiadskiego 27, 50370 Wrociaw, Poland. 0009-2614/93/S 06.00 0 1993 Elsevier Science Publishers B.V.
has been by employing tests for and Pb2 the usefulness developed method includes both is called relativistic multireference
confgurathe proposed and the
all the relativistic effects, including the spin-orbit term. In the RCI calculations, the spin-orbit integrals are transformed in the natural orbital basis obtained from the CASSCF/CI methods and then added to the one-electron Hamiltonian matrix at the RCI step. Consequently, the spin-orbit coupling term is included variationally rather than perturbatively. The RCI code developed up to now is restricted to calculations that include up to 25000 configurations. The configurations included in the CI process contribute in a different way to the total energy of the system. We propose a selection procedure based on the perturbation energy expression to calculate an energy contribution of any particular configuration when added to the reference configuration set. Based on a chosen threshold, the configurations are selected in the new relativistic multireference configuration interaction (RMRCI) calculation. Contributions of the unselected configurations are calculated, assuming that there are no direct couplings between contributions. The contributions of unselected configurations depend linearly on the assumed threshold and an extrapolation procedure has been developed to account for them. In addition the Davidson correction for unlinked quadrupole clusters is also included in the RMRCI. It is shown that the results of the RMRCI method are quite comAll rights reserved.
87
Volume 2 15, number 1,2,3
CHEMICAL PHYSICS LETTERS
parable to fully variational RCI expansions of much larger sizes.
26 November 1993
figurations. The efficiency of the program can be improved by applying the configuration selection procedure described in section 3.
2. General CASSCF/MRSDCI/RCI method 3. Method of configuration selection The calculations follow three main steps. Molecular orbitals for the MRSDCI calculations are generated by the complete active space MCSCF method (abbreviated as CASSCF). In this method a set of the most important electrons for chemical bonding (active electrons) are distributed in all possible ways among a chosen set of orbitals referred to as the internal or active set of orbitals. The CASSCF method provides a starting set of orbitals for the inclusion of higher-order correlation effects. Higher-order electron correlation effects not included in the CASSCF, are taken into account using the configuration interaction method. The CI calculations are multireference singles + doubles CI (MRSDCI). The MRSDCI calculations include a subset of reference configurations in the SOCI, determined by the most important configurations in the CASSCF. Subsequently, single and double excitations from each of the reference configurations are allowed in the MRSDCI. The RCI calculations are carried out following the CASSCF/MRSDCI methods. The CASSCF/ MRSDCI calculations described above include electron correlation corrections as accurately as possible, but neglect the spin-orbit coupling term entirely, which is far from negligible for molecules containing heavy atoms. In the RCI method the spin-orbit integrals are transformed into the natural orbital basis obtained from the CASSCF/MRSDCI methods and then added to the one-electron Hamiltonian matrix at the RCI step. The spin-orbit coupling term is thus included variationally. The RCI calculations include all configurations that mix in the presence of the spin-orbit coupling term as reference conflgurations. The inclusion of spin-orbit coupling term in the Hamiltonian changes the normal point group symmetry of a molecule into the spin double group. Thus, all configurations that have the same symmetry in a spin double group could mix in the RCI. The RCI code developed up to now does not have the capabilities of the direct CI procedures and is restricted to calculations that include up to 25000 con88
The goal of the configuration selection procedure is to select configurations which contribute significantly to the energy when included variationally in the CI space. One convenient way of doing this is to select a set of most important configurations (main reference space) and to check the energy lowering when a tested configuration is added to the reference set. The energy lowering generated by a tested configuration #’ by interacting with configuration em which is a member of the main reference set can be calculated by the perturbation expression [ 9 ]
The energy lowering implied by eq. ( 1) was first used by Whitten and Hackmeyer [ 10,111 as a contiguration selection technique in non-relativistic CI computations. Different selection techniques based on eq. ( 1) have been developed and used successfully by Bender and Davidson [ 12 1, Shavitt [ 131 and most extensively by Buenker and Peyerimhoff [ 14,151. In this work, configurations are tested against every configuration which is a member of the main reference set using eq. ( 1). The configuration is being selected when AE exceeds the assumed threshold for at least one main configuration.
4. Accounting of unselected con&u&ions Expression ( 1) may also be used for the estimation of contribution of unselected configurations,
(2) where the contribution of the unselected contiguration is calculated against the zeroth-order CI wavefunction obtained within the reference space. Eq. (2) does not account for couplings between configurations which are outside the reference space. The ex-
Volume 215, number 1,2,3
CHEMICAL. PHYSICS LETTERS
act energy (threshold equal to zero) can be expressed as E(O)=E(T)+I(T)K(T).
(3)
For a particular selection threshold T, the function K(T) has to be scaled to account for internal couplings and possible relaxation processes in the CI performed with no threshold. The proportionality factor 1, in general, may be a function of the threshold. Bunker and Peyerimhoff [ 16 ] proposed that A is independent of the threshold for non-relativistic CI computations. However our test calculations indicate a linear dependence on the threshold T namely, A(T)=A+BT.
(4)
The performed RMRCI calculations are of the second-order type. The difference between the RMRCI and the full configuration interaction calculations may be estimated by the unlinked quadrupole cluster correction formula proposed by Davidson et al. [ 17,18 1, extended later for a multireference case by Peyerimhoff and Buenker [ 191. The correction expression is given as
where E,, is a reference secular equation and EC1 the second-order energy. Eq. (5) has been shown to be a good correction for the lack of size extensivity of the energy with single and double excitations [ 20,2 11. Since the full relativistic CI calculations are not yet available for comparison, the correction should be considered cautiously as an estimated measure of the missing contribution to the total energy.
5. Computational details To test the present method, we have performed calculations for some molecules with heavy elements such as Pbz and Biz. These two dimers were considered before using the RCI method [ 22,23 1. The natural orbitals for various points in the potential energy curves have been generated from the first-order configuration interaction (FOCI) method following the complete-active-space multicontiguration self-
26 November 1993
consistent field (CAS-MCSCF) [ 8 ] calculations. All CASSCF/FOCI calculations are performed using the modified version [ 8 ] of the ALCHEMY II #l code. In all calculations, we have employed the relativistic effective core potentials (RECPs) [ 241 for both the Pb and Bi atoms. Thus the valence space contains 6s26p2 and 6s26p3 shells for the Pb and Bi atoms, respectively. The ( 3s3p 1d ) valence Gaussian basis set has been used. The spin-orbit integrals were computed using Pitzer’s modified ARGOS codes [ 7 1. In the RMRCI calculations, the low-lying electronics states which yield an n state of the desired symmetry have been included as reference contigurations. Along with the reference configurations, the singly- and doubly-excited configurations from these reference configurations have also been included in the RMRCI calculations. For the Pb2 dimer, the 0: state calculations included the reference configurations from the 3Z; ground state and the excited states such as ‘Ez, 311g,‘Zz (II) and ‘II, (II) states. For the 0: state calculations of Bi2, the reference configurations which we included are the ‘El ground state and the 0: states arising from the excited states such as 31’Ig,‘E:, ‘II, (II), ‘IS; and ‘El (II). The 1U state calculations of Biz are performed by using the 3Z,*( l,), 3AU( 1,) and ‘Zy ( 1,) states as the reference configurations,
6. The effect of configuration selection procedure on energies and properties As seen from tables 1 and 2 the dramatic reduction of selected configurations introduces very small absolute error to the total energy. The threshold of 1 phartree gives errors lying on the border of precision of calculations, despite the reduction of number of configurations to one third of the original number of CSFs. The largest threshold ( 100 uhartree) reduces the number of configurations almost ten times, but the relative error to the total energy is less than 0.05% in all cases studied here. The potential energy surfaces (figs. 1 and 2) as well as geometries (see table 3)) are also well reproduced, even for the highest threshold employed. The impact of c1 The major authors of ALCHEMY II are B. Lengsfield, B. Liu and M. Yoshimine.
89
Volume 2 15, number 1,2,3
CHEMICAL PHYSICS LETTERS
26 November 1993
Table 1 Number of selected configurations, total energy and absolute error for different selection thresholds in hartree atomic units for the 0: ground and excited states of Bis (2.70 A) Selection threshold (Wree)
Number of selected configurations
0: energy
AEX 10-s
energy
AEX lo-2
0 1.0 20.0 40.0 60.0 80.0 100.0
18113 7393 3843 3059 2579 2291 2058
-10.71161 -10.71154 - 10.71031 - 10.70958 - 10.70871 - 10.70829 - 10.70748
0.0073 0.1299 0.2029 0.2899 0.3319 0.4129
-
0.0085 0.1015 0.1710 0.2440 0.2892 0.3561
0: (II)
K(T)XlO_2
10.63308 10.63300 10.63206 10.63137 10.63063 10.63019 10.62951
-0.0036 -0.1528 -0.2291 -0.3126 -0.3603 -0.4416
Table 2 Number of selected configurations, total energy and absolute error (in atomic units) for different selection thresholds for 1, states of Bi2(2.70 A) Selection threshold ( phartree)
Number of selected configurations
energy
AEX lo-2
energy
AEX lo-2
energy
AEX lo-2
0 1.0 20.0 40.0 60.0 80.0 100.0
13780 4377 2405 1965 1733 1581 1455
-
0.0006 0.1088 0.1783 0.2506 0.2835 0.3600
-
0.0006 0.1066 0.1688 0.2410 0.2741 0.3474
-10.61180 -10.61172 -10.61068 - 10.60999 - 10.60926 - 10.60892 - 10.60816
0.0008 0.115 0.1813 0.2540 0.2882 0.3641
1.
l”(II)
10.66207 10.66201 10.66098 10.66028 10.65956 10.65923 10.65847
K(T)XlO_2
l”(III)
10.64028 10.64022 10.63921 10.63859 10.63787 10.63754 10.63681
-0.0024 -0.1121 -0.1786 -0.2365 -0.2635 -0.3243
-0.800,
-0.*20\ 2.5
2.75
3 Distance
3.25
3.5
&)
Fig. 1, Potential energy curves for the ground state (0: ) of Bis for different selection thresholds, (-) T=O phartree, (---) T=20 phartree, (---) T= 100 bhartree.
90
2.7
2.9
2.9
3 Distance(
3.1
3.2
3.3
‘A)
Fii. 2. Potential energy curves for the ground state of Pbs for different selection thresholds, (-) T=O phartree, (---) T=20 phartree, ( ---) T= 100 phartree.
CHEMICAL PHYSICS LETTERS
Volume 215, number 1,2,3
Table 3 Ground state bond distances (A) of Bia and Pb2 calculated with and without spin-orbit coupling for different selection thresholds T (in phartree ) Molecule
Spin-orbit
T=O
T=20
T=lOO
Bi2 Bi2 Pbt Pbs
Yes no Yes no
2.1564 2.7614 2.9443 2.8855
2.7683 2.7734 2.9425 2.8902
2.7586 2.7746 2.9346 2.9115
Table 4 The energy lowering (in eV) due to spin-orbit coupling for the ground state of Biz and Pbs for different thresholds T (in @hattree) Molecule
TsO
T=20
T= 100
Bi2 Pbs
0.88 1.45
0.89 1.46
0.91 1.47
the spin-orbit coupling on the energy change is also not significantly changed, even for a selection threshold of 100 uhartree (table 4). The reproduction of the wavefunction is critical for a proper description of several electronic properties. The main reference configuration sets which have been chosen (73 reference configurations for Bi2 and 56 reference configurations for Pbz) were reasonably extended, leading to the total weight (sum of the squares of the coefficients of reference configurations) as 0.98 1 for Bil and 0.960 for Pbz. The wavefunction structure is almost unchanged for different assumed thresholds. For instance, the wavefunction of the O,+ state of Biz is composed of ‘E: (74.2%) and ‘II, (13.8W). When T is equal to 20 uhartree (and 100 uhartree), the wavefunction of the 0: state is made of ‘C: (74.1%) and 311p (13.91%) (‘ZZ - 74.2% and ‘II, - 14.0%). The wavefunctions of the 1, state of Bi2, with T=O, 20 and 100phartree, are composed of 3G (80.8%) and 3A,, (14.4%), ‘E; (80.9W) and ‘A,, (14.4%) and 3Ey (81.1%) and 3A, ( 14.3%), respectively.
7. Energy extrapolation and the Davidson correction
results are known and J(T) function is calculated (fig. 3). The function is linear for a wide range of threshold values. Error bars for nonlinearity are obvious, since configurations belong to discribed sets rather than being contiuously distributed according to the energy thresholds. The function K(T) reproduces closely the energy difference E( 0) -E(T) (tables 1 and 2). Using eqs. (3) and (4) the values for at least three different thresholds are needed to determine the A and B coefficients of L(T) function. The value of the correction A( T)K(T) is added to E(T) for calculations with the smallest energy threshold studied. The averaged value of 1 ( T= 20) using energy thresholds in the range between 20 and 100 uhartree (fig. 3 ) resulted as 0.8 1 for the 0: state and 0.94 for the 1, state. These values are close to the exact I of 0.85 (0: ) and 0.99 (1,). Our linear form for 1( T)is more general than the constant value of 1 (independent of T) proposed by Buenker and Peyerimhoff [ 16 1. None of the cases studied in this work could be reasonably described by a parameter which is independent of T. The full configuration interaction energy was estimated by a standard expression (5). The results for a number of excited states of Biq are presented in table 5. Due to an extended reference space (Cc: for a reference space is larger than 0.98 for all states) the correction is very small.
0.60
I 0
The total MRSDCI energy may be calculated from eq. (3), assuming that the proportionality factor A( T) is known. In the Bi2 case, exact (T= 0)RMRCI
26 November 1993
20
1 40
60
I 60
100
T WI
Fig. 3. The A(T) scaling factor values for different electronic states of Bis as a function of the selection threshold.
91
Volume 215, number 1,2,3
CHEMICAL PHYSICS LETTERS
Table 5 Full configuration interaction energy estimation Q-corrected energy for different states of Bir Electronic state
AEX 10”
Q-corrected energy
0:
- 0.408 -0.475 -0.158 -0.166 -0.198
-10.71191 - 10.63355 - 10.66222 - 10.64045 - 10.61200
0: (II) 1” l”(II) l”(III)
8. Conclusions The efficiency of the relativistic configuration interaction (RCI) method has been improved by employing the individualized configuration selection procedure. The error introduced by the configuration selection procedure is well controlled by the calculation of the contribution of unselected configurations to the total energy. The performed tests of energy as well as properties indicate the usefulness of the proposed approach. The molecular wavefunction is also well reproduced for the electronic states considered here. The calculated contribution of unselected configurations is used for the estimation of full second-order CI energy. The linear dependence of the scaling factor A on the threshold Tallows the estimation of full SOCI energies for at least three different thresholds. The averaged A ( TMIN)value leads to a more precise value of the full SOCI energy. The contribution of unlinked quadruple clusters to the total energy (Q correction) is also calculated using the Davidson correction.
Acknowledgement This research was supported in part by US Deof Energy under grant DEFGpartment 0286ER13558.
92
26 November 1993
References [ 1] IL Balasubramanian and KS. Pitzer, in: Ab initio methods in quantum chemistry, ed. K-P. Lawley (Wiley, New York, 1987). [2] K. Balasubramanian, J. Phys. Chem. 93 (1989) 6585. [3] Y.S. Lee, W.C. Ermler and K.S. Pitzer, J. Chem. Phys. 67 (1977) 5861. [4] P.A. Christiansen, Y.S. Lee and KS. Pitzer, J. Chem. Phys. 71 (1979) 445. [ 51 P. Hafner and W.H.E. Schwartz, Chem. Phys. Letters 65 (1979) 537. [6] WC. Ermler, Y.S. Lee, P.A. Christiansen and K.S. Pitzer, Chem. Phys. Letters 8 1 ( 198 1) 70. [7] R.M. Pitzer and N.W. Winter, J. Phys. Chem. 92 (1988) 3061. [8] K. Balasubramanian, J. Chem. Phys. 89 (1988) 5731. [9] J.C. Slater, Quantum theory of atomic structure, Vol. 1 (McGraw Hill, New York, 1960). [ lo] J.L. Whitten and M. Hackmeyer, J. Chem. Phys. 51 (1969) 5584. [ 111 M. Hackmeyer and J.L. Whitten, J. Chem. Phys. 54 (1971) 3739. [ 121 C.F. Benderand E.R. Davidson, Phys. Rev. 183 (1929) 23. [ 13 ] I. Shavitt, in: Modem theoretical chemistry, VoL 3. Methods of electronic structure theory, ed. H.F. Schaefer (Plenum Press, New York, 1977) p. 189. [ 141 R.J. Buenker and S.D. Peyerimhoff, Theoret. Chim. Acta 35 (1974) 33. [ 151 R.J. Buenker and S.D. Peyerimhoff, in: New horizons of quantum chemistry, eds. P.O. Lawdin and B. Pullman (Reidel, Dordrecht, 1983) p. 183. [ 161 R.J. Buenker and S.D. Peyerimhoff, Theoret. Chim. Acta 39 (1975) 217. [ 171 E.R. Davidson, The world of quantum chemistry, eds. R. Daudel and B. Pullman (Reidel, Dordrecht, 1973) p. 17. [ 181 S.R. Langhoff and E.R. Davidson, Intern. J. Quantum. Chem. 8 (1974) 61. [ 191 SD. Peyerimhoff and R. Buenker, Chem. Phys. 42 (1979) 167. [20] J. Paldus, P.E. Warmer, F. Vissar and A. Van der Avoird, Proceedings of the 5th Seminar on Computational Methods in Quantum Chemistry, Groningen, The Netherlands (1983) p. 13. [21] D.B. Knowles, J.R. Alvarez-Collado, G. Hirsch and R.J. Buenker, 3. Chem. Phys. 92 (1990) 585. [22] K. Balasubramanian, Chem. Rev. 90 (1990) 93. [23] K. Balasubramanian and D.-W. Liao, J. Chem. Phys. 95 (1991) 3064. [ 241 R.B. Ross, J.M. Powers, T. Atashroo, W.C. Ermler, L.A. LaJohn and P.A. Christiansen, J. Chem. Phys. 93 (1990) 6654.