10 September 1999
Chemical Physics Letters 310 Ž1999. 530–536 www.elsevier.nlrlocatercplett
A full-configuration benchmark for the N2 molecule Elda Rossi a , Gian Luigi Bendazzoli b, Stefano Evangelisti
b,c,)
, Daniel Maynau
c
a
CINECA, Via Magnanelli 6 r 3, I-40033 Casalecchio di Reno (BO), Italy Dipartimento di Chimica Fisica e Inorganica, UniÕersita` di Bologna, Viale Risorgimento 4, I-40136 Bologna, Italy Laboratoire de Physique Quantique, UMR 5626, UniÕersite´ Paul Sabatier, 118 Route de Narbonne, 31062 Toulouse Cedex, France b
c
Received 10 June 1999; in final form 12 July 1999
Abstract A full-configuration interaction ŽFCI. calculation has been performed for the nitrogen molecule using an ANO w4s3p1dx basis set. The FCI space for such a system contains about 9.68 = 10 9 symmetry-adapted Slater determinants. The FCI results are compared with several approximate methods, both of single- and multi-reference type, in order to test their accuracy. q 1999 Elsevier Science B.V. All rights reserved.
1. Introduction Full-configuration interaction ŽFCI. calculations have an important role in quantum chemistry, since they provide benchmark results useful for assessing the accuracy of approximate methods. Several FCI algorithms exist w1–5x. An FCI algorithm has been proposed and implemented in Bologna w6,7x, and used for the study of ab initio and semi-empirical Hamiltonians. The algorithm has been implemented on massively parallel computers w8–11x, where huge calculations can be performed. This permits a comparison of the FCI result with approximate calculations, like perturbation theory ŽPT. or configuration interaction ŽCI. calculations. In particular, FCI results permit us to investigate the size consistency of approximate methods, i.e., whether the energy of two non-interacting systems is indeed the sum of the energies of the two individual systems. ) Corresponding author: Fax: q33-5-61556065; e-mail:
[email protected]
In this Letter, we present the results of a FCI calculation on the ground state of the nitrogen dimer N2 . The calculation was done using a valence atomic natural orbital ŽANO. basis set. The size of the FCI space is 9 678 561 408 symmetry-adapted Slater determinants. Due to the large dimension of the FCI vector, only one internuclear distance, close to the equilibrium geometry of the molecule, was investigated. The size-consistency property of FCI was used in order to obtain the FCI energy of the dissociated system, so that the ‘dissociation energy’ of the molecule Ži.e., the energy difference between the dissociated system and the molecule at the approximated equilibrium distance. could be obtained. The FCI results are compared with several approximate methods, both of single- ŽSR. and multireference ŽMR. type, in order to test their accuracy. In particular, the total energies at equilibrium can be compared in the case of SR methods. Since SR methods cannot follow the dissociation of the molecule, only MR methods can be used to compute the dissociation energy of the molecule.
0009-2614r99r$ - see front matter q 1999 Elsevier Science B.V. All rights reserved. PII: S 0 0 0 9 - 2 6 1 4 Ž 9 9 . 0 0 7 9 1 - 5
E. Rossi et al.r Chemical Physics Letters 310 (1999) 530–536
2. Computational details We used the parallel version of our FCI algorithm that is adapted to CRAY T3E computers. In its present implementation, the method requires only two vectors of the dimension of the FCI space: one for the coefficients of the wavefunction
531
is of about 80 Gbytes, with a global IrO traffic per iteration of more than 500 Gbytes. The wall-clock time is about 26 500 s per iteration, mainly because of the huge IrO load. In fact, the CPU time needed to produce Y from X is about 13 500 s, while 13 000 s are spent in IrO operations. In order to reduce the overall IrO time, we performed the first iterations neglecting those components of X and Y that are below a given threshold Žin our calculation, we chose 10y8 .. It is only when the iterative process is almost at convergence that a small number of zero-threshold FCI iterations is performed. In this way, the real FCI process is reduced to a very small number of iterations, while the preliminary iterations needed to compute the approximate vector require limited IrO resources. For the present calculation on nitrogen dimer, the ANO basis set optimized by Widmark et al. w14x was used, adopting a w4s3p1dx contraction. The actual group used in all calculations was the D 2h subgroup of the full group of the molecule D`h . We optimized the molecular orbitals ŽMO. through a valence complete-active-space self-consistent field ŽCAS-SCF. calculation, with 10 active electrons in eight MOs. The two lowest MOs, of sg and su symmetries, respectively, corresponding to the two 1s orbitals of the nitrogen atoms, were kept frozen in the FCI calculation at the CAS-SCF level. This leaves a total of 34 MOs and 10 electrons for the FCI calculation. The orbital symmetries of the FCI orbitals are 7sg , 7su , 4 pu , 4 pg , 1 dg , 1 du in D`h Ž8a g , 4b 3u , 4b 2u , 1b 1g , 8b1u , 4b 2g , 4b 3g , 1a u in D 2h ..
3. I r O considerations As mentioned before, the full-CI code is heavily IrO bounded. In fact, two coefficient vectors must be stored on disks. In order to decrease storage requirements, symmetry considerations are taken into account Žthe vectors are blocked diagonal Ž8 blocks for the actual case. and symmetric.. In the present case, with 10 billion determinants, about 10 billion real numbers must be stored on disk, resulting in about 80 Gbytes of disk occupation. Moreover, during each iteration, the whole
532
E. Rossi et al.r Chemical Physics Letters 310 (1999) 530–536
WRITE operations, with 40 Gbytes transferred each, for a total data traffic of about 500 Gbytes per iteration. This is obviously a real problem, particularly if we consider that the CINECA T3E is equipped with a scratch file system of 130 Gbytes of DD-308 disks with a sustained bandwidth of 8–12 Mbyters. A normal IrO transfer of 500 Gbytes at a speed of 10 Mbyters would require 50 000 s only for IrO activity for every single iteration. In order to solve this problem, we adopted a heavy IrO optimization, focused to ‘parallelize’ as much as possible the IrO activity on a number of independent disks connected to the system through high bandwidth channels able to host many of them. We can distinguish all the optimization interventions between ‘normal’ optimization and ‘parallel’ optimization. As far as the first part is concerned, we take advantage of the fact that we need to handle a small number of very large records Žabout 10 million words in the present case.. Using large transfer sizes alleviates the longer system–call processing time. When system intervention is avoided, data are transferred directly between user buffer in memory and disks, bypassing system buffers. This happens only if the IrO request is ‘well formed’, i.e., Ži. the size of the request is a multiple of the disk sector Ž512 words in this case., and Žii. the user area in memory is located on a 8-word boundary Ža specific compiler directive exists for this.. To meet these requirements, our memory buffers were rounded to the nearest multiple of the disk sector, with a very limited increase Žonly 0.02%. of the memory occupancy. Direct access type IrO Žasynchronous queued IrO. was adopted to let multiple processors access the same file. In this case, it is important to coordinate all the processors performing IrO on the same file, to let the IrO requests be issued sequentially. In the other case, when the processors issue the request at their own speed, the host will interpret this as a non-sequential extension of a file, and the sustained bandwidth can be reduced to one-half. Adopting these measures, we were able to reach almost the peak bandwidth of the disks. In order to go over this limit, we decided to ‘parallelize’ the IrO activity, taking advantage of the structure of the /tmp file system. It is composed of a set of DD-308 disks arranged on five indepen-
dent ‘fibre channel loops’, each loop hosting from 3 to 10 disks. The five loops are connected, through a FCN-1 interface, to the system nodes with a GigaRing IrO channel. Different IrO requests issued on different disks at the same time can proceed concurrently so that, in principle, throughput can scale accordingly. DD-308 disks are characterized by a 10 Mbyters rate average, fibre channel is a serial channel that has a 100 Mbyters peak bandwidth, the GigaRing channel is a counter-rotating, dual-ring channel that can deliver 1 Gbyters of data. The calculation on the nitrogen dimer could not use a so large a number of disks because of production rules and limitations. Four files were used instead and a speed of 36 Mbyters average Žout of a peak of 40. was reached. This result, even if partial, allowed us to reduce considerably the elapsed time for each single iteration, at least for the last ones, those requiring the complete treating of the two coefficient vectors. In particular, instead of 50 000 s, only 13 000 s were spent for IrO activity, resulting in a total wall clock time of about 7.5 h for every single final iteration.
4. Results and discussion The main purpose in performing FCI calculations is to compare the FCI result with the analogous result obtained with approximate methods. A first comparison can be done on the total energy; a good approximate method should reproduce the FCI energy as closely as possible. However, for chemical applications, the relevant quantity is usually the energy difference between different geometries. Therefore, we are interested not only in the absolute error of the approximate methods with respect to FCI, but also in the error of the dissociation energy of the molecule Žassuming a fixed value for the equilibrium distance, in the present case R s 2.10 bohr.. Therefore, we performed calculations at R s 2.10 bohr, near the equilibrium bond length, and at ‘infinite’ Ži.e., very large. internuclear distance. For the dissociated molecule, we chose the value of R s 40 bohr, a distance at which the energy differs from the asymptotic limit by less than 10y8 hartree. In principle, CAS-SCF is size consistent, unless symmetry breaking occurs. We verified that we do
E. Rossi et al.r Chemical Physics Letters 310 (1999) 530–536 Table 1 Total energies for the different methods we report: the energy at 2.10 bohr Ž E2.10 ., the energy at infinite distance Ž E` .. E2.10 and E` are in hartree Method
E2.10
E`
SCF CAS-SCF MRCI FCI
y108.97987753 y109.13214769 y109.31840560 y109.325905r5
– y108.80149555 y108.99624144 y109.00517879
not have such a breaking for N2 with our choice of the active space, and that the CAS-SCF energy at 40 bohr is indeed twice the corresponding energy of an isolated nitrogen atom. In a given basis set, FCI is size consistent. However, since we performed a valence FCI calculation, one must be sure that the frozen orbitals of the dimer and the two monomers span the same orbital space. Again, this is true if no symmetry breaking occurs at the CAS-SCF level. Since the CAS-SCF solution is size consistent, the Žvalence. FCI energy at dissociation is twice the FCI energy of a single atom, even in presence of uncorrelated Žfrozen. orbitals. Our benchmark can be used for a comparison with approximate calculations. Therefore we report in Table 1 the FCI total energies at R s 2.10 bohr and R s `. In order to make the comparison with our results easier, we also report in Table 1 the SCF, CAS-SCF and MR-CI results for the two distances. The approximated calculations we performed belong to two different classes: SR and MR calculations. SR calculations Žat SCF and post-SCF levels. are suitable for describing the molecule close to equilibrium, but in general fail to describe dissociation. In Table 2 we report, for R s 2.10 bohr, the energy difference from FCI obtained using SR methods. These are SCF, Moller–Plesset PT at second-, third- and ¨ fourth-order ŽMP2, MP3, and MP4, respectively., single- and double-CI and CC ŽCISD and CCSD.. For the CCSD case, the effect of non-iterative inclusion of triple excitations has also been considered. In the case of SCF, the error is obviously very large. The error goes down from 346 mhartree for SCF to 40 and 16 mhartree for CISD and CCSD, respectively. By adding the perturbative contribution of triple excitations in CC, the error is significantly reduced, going down to less than 0.5 mhartree for
533
the CCSDwTx method. The MP2 result is extremely good, with an error similar to CCSD. MP3 has a larger error and so is for truncated MP4. Only for the full MP4 is the error substantially reduced. In Table 3, the energy differences with respect to FCI for MR methods are shown, both at R s 2.10 bohr Ž D2.10 . and R s ` Ž D` ., for PT- and CI-like methods. The perturbative methods considered are second-order PT over a CAS-SCF wavefunction ŽCAS-PT2, as implemented in MOLCAS w15x, and the analogous approach RS, as implemented in MOLPRO w16x.. The CAS-SCF results are about 200 mhartree higher than the corresponding FCI values, with an error of about 10 kcal in the dissociation energy. Second-order PT reduces the absolute errors D2.10 and D` by an order of magnitude Ždown to 7–30 mhartree, depending on the different choices of H0 .. However, the dissociation error DD is reduced by a much smaller factor: about 2r3, and we have both positive and negative errors. ŽHowever, the use of larger basis sets would certainly give a comparison less favorable to CAS-SCF.. The differences between the MOLCAS and MOLPRO results are due to a different definition of the unperturbed Hamiltonian Ž‘diagonal’ versus ‘not diagonal’., and to a different treatment of single excitations. In the RS case, we report the results obtained with four possi-
Table 2 SR methods: energy difference from FCI Žin mhartree. at Rs 2.10 bohr Ž D2.10 . Method
D2.10
SCF
346.028
MP-2 ŽD. MP-3 ŽD. MP-4 ŽDQ. MP-4 ŽSDQ. MP-4 ŽSDTQ.
16.846 25.187 19.609 14.190 -3.741
CISD CISDqDavidson CISDqPople
40.356 8.175 12.930
CCSD CCSDŽT. CCSDwTx CCSD-T
16.320 1.340 0.480 1.746
E. Rossi et al.r Chemical Physics Letters 310 (1999) 530–536
534
Table 3 MR methods: energy difference from FCI Žin mhartree. at Rs 2.10 bohr Ž D2.10 ., Rs` Ž D` ., and difference between these two values, DD s D` y D2.10 Method
D2.10
D`
CAS-SCF
193.757
203.683
9.926
CAS-PT2 RS2 g0 RS2 g1 RS2 g2 RS2 g3 RS2 g4
14.653 20.468 20.837 20.610 20.481 19.827
6.669 13.875 27.454 28.643 28.643 27.534
y7.984 y6.593 6.617 8.033 8.162 7.707
7.052 6.688 7.045 7.050 6.804
4.851 2.116 3.253 3.253 2.185
y2.201 y4.572 y3.793 y3.798 y4.620
C-MRCI C-MRCIqDavidson C-MRCIqPople C-MRCI ACPF
9.690 y0.639 1.144 1.047
11.604 y2.583 y0.261 y0.392
1.914 y1.944 y1.404 y1.439
MRCI MRCI ACPF MRCI AQCC ŽSC. 2-MRCI
7.500 y1.707 0.520 2.712
8.937 y3.746 y0.379 3.124
1.438 y2.039 y0.899 0.412
RS3 g0 RS3 g1 RS3 g2 RS3 g3 RS3 g4
DD
ble definitions of H0 . The default ‘g0’ partition, the three possibilities ‘g1’, ‘g2’, ‘g3’, w17x, as well as a fourth ‘g4’ partition, which gives size-consistent results in the case of dissociation into high-spin openshell systems w16x. Not surprisingly, the third-order PT ŽRS3. gives definitely better results, with an absolute error of 7 mhartree or below, and a dissociation error of about 4 mhartree. As a general remark, we can say that perturbation methods are very sensitive to the choice of H0 , and the differential effects in the non-dynamical correlation are difficult to reproduce. If we compare CI and PT2 methods, at both SR and MR levels, it appears that the MR result is much better than the SR one at CI level Ža reduction of about a factor 5 in the error at equilibrium., while the PT results are of the same quality at SR and MR levels. This fact, together with the high variability of the MR-PT result depending of the choice of H0 , may give rise to some doubts about the splitting of the Hamiltonian into unperturbed and perturbative parts.
We come now to CI-like methods. We have considered both contracted ŽC-MRCI. w18x and uncontracted ŽMRCI. approaches. Standard uncontracted MRCI is 7.5 mhartree higher than FCI at equilibrium, and about 9 mhartree at dissociation, with a DD of about 1.5 mhartree. The effect of contracting the CI wavefunction is relatively small, giving a slightly larger DD , about 2 mhartree. We compare the MRCI results with three modified methods that try to correct the size-consistency error of MRCI: the averaged coupled pair functional ŽACPF. w19x; the averaged quadratic coupled cluster ŽAQCC. w20–22x; and the size-consistent self-consistent CI ŽŽSC. 2-CI. w23–28x. In the present implementation of the ŽSC. 2 method, the dressing is computed with respect to a single Slater determinant. In the case of N2 at dissociation, many different determinants are dominant in the MRCI wavefunction. In order to solve this problem, we computed the ŽSC. 2 energy at infinity for a high spin multiplicity, i.e., we chose the value S z s 3 Žwe must note, however, that this procedure is unable to give the whole dissociation curve for this molecule.. The size-consistency corrected results with all the three methods are much closer to the FCI result at equilibrium. However, the amount of the correction can be different at equilibrium and dissociation so the error in the dissociation energy can be relatively large. In some cases, it can even become larger than the corresponding MRCI error Žthis is the case for ACPF.. Of the three methods we considered, ŽSC. 2 gives the smallest error for the dissociation energy. In order to test whether the good agreement between the ŽSC. 2 dissociation energy and the FCI result is indeed related to a correct description of unlinked terms, we computed the size-consistency error at dissociation. This is defined as the difference
Table 4 Size-consistency error: we report the difference 2 EŽN. – E`ŽN2 . Žin mhartree. Method
2 EŽN. – E `ŽN2 .
MRCI MRCI ACPF MRCI AQCC ŽSC. 2 -MRCI FCI
5.257 y2.604 y1.935 0.788 0.000
E. Rossi et al.r Chemical Physics Letters 310 (1999) 530–536 Table 5 Population of the FCI vector Class Threshold Population 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
0.1 Eq00 0.1 Ey01 0.1 Ey02 0.1 Ey03 0.1 Ey04 0.1 Ey05 0.1 Ey06 0.1 Ey07 0.1 Ey08 0.1 Ey09 0.1 Ey10 0.1 Ey11 0.1 Ey12 0.1 Ey13 0.0
Total
3 146 3 064 60 644 893 583 9 029 290 70 875 347 342 964 883 1 158 615 041 2 279 641 912 2 752 930 763 2 010 453 996 828 103 800 195 031 342 29 957 594
5 51 1.17156007 2.55456485 8.46870300 12.85891813 22.91372024 22.86326905 19.65730864 10.40084709 3.75548276 0.82350651 0.10791016 0.00870814 0.00039229 0.00001003 0.00000015
5 52 0.89679379 0.06259109 0.03552487 0.00410516 0.00088797 0.00008839 0.00000825 0.00000047 0.00000002 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000
9 678 561 408 105.58490110 1.00000000
535
Compared to 5 X 5 2 , 5 X 5 1 is much more dispersed, reaching its maximum in the 10y4 –10y5 class and going very slowly to zero. This seems to indicate that a very large number of small coefficients is required to obtain an accurate description of the FCI vector. It is impossible to judge from our results the amount of information contained in such a large number of small coefficients corresponding to highorder excitations. In particular, it can be expected that a good approximation to high-order excitations could simply be obtained as the product of coefficients of lower excitations, as done in coupled-cluster formalism. To investigate this point, we are considering for performing a ‘cluster analysis’ of the FCI vector, with the same formalism described in Ref. w7x. However, because of the huge computational effort required by such an analysis, this will be discussed in a future work. 5. Conclusions
between twice the energy of an isolated atom and the energy of two atoms at infinite distance, 2 EŽN. – E`ŽN2 .. The results are reported in Table 4. We see that the ŽSC. 2 method has the smallest error, less than one-sixth of the corresponding MRCI value. On the other hand, both AQCC and ACPF tend to considerably overcorrect the MRCI result. We note that even if the ŽSC. 2-CI formalism does not guarantee an upper bound to the FCI energy, the ŽSC. 2-CI results always lie above FCI. This is not surprising, since many linked contributions, due to the more excited configurations of the CI matrix, are missing. We report in Table 5 the dispersion of the FCI vector at R s 2.10 bohr. As expected, the number of small coefficients is extremely large. It reaches its maximum in the 10y9 –10y1 0 class. For each class of coefficients, we computed the contribution to the norm of the vector, defined as 5 X 5 2 ' Ý < CI < 2 . I
The largest contribution to the 5 X 5 2 is concentrated in the highest class and goes below 10y8 for thresholds smaller than 10y8 . We also computed the sum 5 X 5 1 of the absolute values of the coefficients, 5 X 5 1 ' Ý < CI < . I
We have shown that FCI benchmark calculations involving up to 10 10 Slater determinants, and taking into account all the CI coefficients, are now possible. This permits the comparison of approximate versus FCI results using reasonable basis sets for small molecules. In the present case, FCI results on the N2 molecule have been used to test the accuracy of a number of approximate methods. Our conclusions are derived from a single molecule, so one must be cautious about the generality of these results. As far as this system is concerned, some conclusions can be drawn. At the equilibrium distance, the good performance of MP2 compared to CIŽSD. is striking. A similar behaviour is not shown by the MR perturbative approaches when compared to MRCI. A important second point concerns the size-consistency correction. It is relatively easy to improve the MR-CI results to include in an approximate way the effect of unconnected higher excitations. However, a much more difficult task is to do that in a balanced way for different values of the internuclear distance, so as to improve the dissociation energy. With respect to this point, it appears that the ŽSC. 2 method gives excellent results. While in the present implementation it involves a SR dressing, generalization to a MR dressing is in progress.
536
E. Rossi et al.r Chemical Physics Letters 310 (1999) 530–536
Acknowledgements This work was partly supported by EEC within the ‘Human Capital Mobility Program’ ŽContract No. ERBCHRXCT96r0086 ., by the Italian ‘Ministero per l’Universita` e la Ricerca Scientifica e Tecnologica’ ŽMURST., and by the University of Bologna. We wish to thank CINECA for supplying computer time on CRAY T3E. References w1x P.J. Knowles, N.C. Handy, Chem. Phys. Lett. 111 Ž1984. 315. w2x J. Olsen, B.O. Roos, P. Jørgensen, H.J.Aa. Jensen, J. Chem. Phys. 89 Ž1988. 2185. w3x S. Zarrabian, C.R. Sarma, J. Paldus, Chem. Phys. Lett. 155 Ž1989. 183. w4x R.J. Harrison, S. Zarrabian, Chem. Phys. Lett. 158 Ž1989. 393. w5x J. Olsen, P. Jørgensen, J. Simons, Chem. Phys. Lett. 169 Ž1990. 463. w6x G.L. Bendazzoli, S. Evangelisti, J. Chem. Phys. 98 Ž1993. 3141. w7x G.L. Bendazzoli, S. Evangelisti, Int. J. Quantum Chem. Symp. 27 Ž1993. 287. w8x S. Evangelisti, G.L. Bendazzoli, R. Ansaloni, E. Rossi, Chem. Phys. Lett. 233 Ž1995. 353. w9x S. Evangelisti, G.L. Bendazzoli, R. Ansaloni, F. Durı, ` E. Rossi, Chem. Phys. Lett. 252 Ž1996. 437. w10x N. Ben-Amor, S. Evangelisti, D. Maynau, E. Rossi, Chem. Phys. Lett. 288 Ž1998. 348. w11x E. Rossi, G.L. Bendazzoli, S. Evangelisti, J. Comp. Chem. 19 Ž1998. 658.
w12x L. Gagliardi, G.L. Bendazzoli, S. Evangelisti, J. Comp. Chem. 18 Ž1997. 1329. w13x E.R. Davidson, J. Comp. Phys. 17 Ž1975. 87. ˚ Malmqvist, B.O. Roos, Theor. Chim. w14x P.-O. Widmark, P.-A. Acta 77 Ž1990. 291. w15x MOLCAS, version 4, K. Andersson, M.R.A. Blomberg, M.P. ˚ Malmqvist, P. Fulscher, G. Karlstrom, R. Lindh, P.-A. ¨ ¨ Neogrady, ´ J. Olsen, B.O. Roos, A.J. Sadlej, M. Schutz, ¨ L. Seijo, L. Serrano-Andres, ´ P.E.M. Siegbahn, P.-O. Widmark, Lund University, Sweden, 1997. w16x MOLPRO is a package of ab initio programs written by H.-J. Werner, P.J. Knowles, with contributions from J. Almlof, ¨ R. Amos, A. Berning, D.L. Cooper, M.J.O. Deegan, A.J. Dobbyn, F. Eckert, S. Elbert, K. Hampel, R. Lindh, A.W. Lloyd, W. Meyer, A. Nicklass, K. Peterson, R. Pitzer, A.J. Stone, P.R. Taylor, M.E. Mura, P. Pulay, M. Schutz, ¨ H. Stoll, T. Thorsteinsson. w17x K. Andersson, Theor. Chim. Acta 91 Ž1995. 31. w18x H.-J. Werner, P.J. Knowles, J. Chem. Phys. 89 Ž1988. 5803. P.J. Knowles, H.-J. Werner, Chem. Phys. Lett. 145 Ž1988. 514. w19x R. Gdanitz, R. Ahlrichs, Chem. Phys. Lett. 143 Ž1988. 413. w20x P.G. Szalay, R.J. Bartlett, Chem. Phys. Lett. 214 Ž1993. 481. w21x P.G. Szalay, R.J. Bartlett, J. Chem. Phys. 103 Ž1995. 3600. w22x L. Fusti-Molnar, ¨ ´ P.G. Szalay, J. Phys. Chem. 100 Ž1996. 6288. w23x J.L. Heully, J.P. Malrieu, Chem. Phys. Lett. 199 Ž1992. 545. w24x J.P. Malrieu, I. Nebot-Gil, J. Sanchez-Marin, J. Chem. Phys. 100 Ž1993. 1440. w25x I. Nebot-Gil, J. Sanchez-Marin, J.P. Malrieu, J.L. Heully, D. Maynau, J. Chem. Phys. 103 Ž1995. 2576. w26x J. Sanchez-Marin, I. Nebot-Gil, D. Maynau, J.P. Malrieu, Theor. Chem Acta 92 Ž1995. 241. w27x I. Nebot-Gil, J. Sanchez-Marin, J.P. Malrieu, J.L. Heully, D. Maynau, Theor. Chem Acta 95 Ž1997. 215. w28x I. Nebot-Gil, J. Sanchez-Marin, J.P. Malrieu, J.L. Heully, D. Maynau, Chem. Phys. Lett. 234 Ž1995. 45.