Use of Integral Approximations in Non-Empirical LCAO MO SCF Calculations on XCN Systems BY G. DOGGETT AND A. MCKENDRICK Chemistry Department University of Glasgow Glasgow W.2. Received 30th October 1968 The accuracy of three multicentre integral approximations has been evaluated for a number of LCAO MO SCF calculations on the ground states of HCN and FCN. Mcban’s results for HCN using a basis set of best atom atomic orbitals are best reproduced by two integral approximations suggested by Lowdin while the Mulliken approximation is poor. All three integral approximations yield total energies lying below the ab initiu value. But double-zeta calculations on HCN and FCN using a fixed CN internuclear separation both yield minimum total energies which straddle the ab initiu values. The calculated equilibrium CH and FC internuclear separations straddle the experi- mental values of 0.1063 nm and 0.1262 nm respectively RCH = 0-1041 nm and RFC = 0.1 196 nm when using the Mulliken approximation and RCH = 0.1 104 nm and RFC = 0-1 34 nm when using the best Lowdin approximation.Dipole moments and one-electron energies are also calculated for each approximation. The charge distributions are discussed in terms of effective atomic charges and one-electron difference density functions. 1. INTRODUCTION There have been few detailed investigations into the accuracy of invoking integral approximations in calculations on polyatomic molecules. Too often integral approximations are made in molecular calculations without first evaluating their accuracy. This could be carried out for e.g.by attempting to reproduce the results of an ab initio calculation. Unless such a procedure is adopted the results obtained from the approximate calculation may be of little value. But detailed evaluation of approximate methods has only become possible recently as a consequence of the rapidly increasing number of ab iizitio calculations. The possibility of finding a reliable approximate method offers some hope for making reasonably accurate calculations on molecules which are at present beyond the scope of a completely ab initio calculation. The linear nitrile systems XCN are of interest and the present work was initiated in an attempt to elucidate the ground state electronic structures of HCN and FCN. The effect of X on the electron distribution in the CN fragment is of particular interest.The most convenient method of calculation is the single configuration LCAO MO SCF approach using a basis of either simple Slater or double-zeta atomic orbitals. The basic problem in this kind of calculation is in the evaluation of the three-centre one- and two-electron integrals particularly when dealing with extended basis sets. The search for a suitable integral approximation is made by first calculating the ground state electronic structure of HCN using the same molecular parameters as in the accurate ab initio calculation of McLean.’ It is then possible to assess the accur- acy of a particular approximate ab iizitio calculation by a direct comparison of the calculated total energy dipole moment one-electron energies and LCAO coefficients with the corresponding values as determined by McLean.l 32 G.DOGGETT AND A . MCKENDRICK 33 The basic problem in approximating a multicentre integral is in finding the opti- mum expansion for a given overlap density in terms of suitably chosen atomic densities. The most widely used approximation is the one suggested by Mulliken,’ hereafter referred to as the MA d ) a ( l ) d ) b ( l ) = *sab{d)a(l)d)a(l) f d ) b ( l ) d ) b ( l ) ] * (1 * 1) A slightly modified version of the MA which is made invariant to orthogonal trans- formations of the atomic orbital basis set-the IMA has also been investigated in view of the work by Pople Santry and SegaL3 Although (1.1) preserves charge it is basically unsatisfactory because the centroid of the overlap density is not maintained. Lowdin suggested that an asymmetric partitioning of the overlap density was more satisfactory d ) a ( l ) + b ( l ) = sab{Al#a(l)d)a(l) + A2d)b(1)d)b(1)} (1.2) where A + A2 = 1.This approximation is referred to as the partial Lowdin approxi- mation or PLA and A2 is chosen so that the dipole moments of the densities of the left- and right-hand sides of (1.2) are identical. A more complete approximation has also been suggested by L o ~ d i n ~ hereafter called the full Lowdin approximation or FLA + a ( l ) + b ( l ) = ~ ~ ~ d ) a ( l ) + a ~ ( l ) s a ~ b + A ~ ~ d ) b ( l ) ~ b ~ ( l ) s b ~ a a’ b’ where again A + A2 = 1 and A2 is chosen to reproduce the dipole moment of the overlap density. The summations over a’ and 6’ include all orbitals on a or b having principal quantum number equal to or less than that of The calculations are performed within the framework of SCF theory as outlined by McWeeny.6 This analysis assumes the basis functions are mutually orthogonal.Thus the f and G matrices initially calculated using a basis of ordinary atomic orbitals must be first transformed by the matrix S-* which symmetrically ortho- gonalizes the set of basis functions or +b respectively. - x = s-+ x s-4 where X stands for f or G. No difficulties were experienced with the convergence of the iterative procedure and the occasional oscillatory behaviour could always be traced to errors in the data. All one- and two-centre molecular integrals are evaluated numerically except for the one-centre two electron integrals which are evaluated analytically. The remaining three-centre integrals are evaluated according to the MA PLA or FLA as required.2. THE ELECTRONIC STRUCTURE OF HCN The molecular geometry for the pilot calculation on HCN is taken from the work of McLean RCH = 0.1058 nm (2-000 a.u.) RCN = 0.1157 nm (2.187 a.u.). The calculated energies and dipole moments are compared with McLean’s results in table 1 for the integral approximations under discussion. The first point to be noted from table 1 is that the use of a multi-centre integral approximation in this study of HCN lowers the total energy relative to the ab initio value. There is also little difference between the MA and IMA calculations. There is no reason to expect the approximate energy to lie above or below the ab initio value. This follows since the use of an integral approximation induces an arbitrary change in the effective Hamiltonian matrix thereby causing the calculated energy to suffer an arbitrary displacement from the ab initio value.If the total 2 34 I N T E G R A L APPROXIMATIONS IN SCF C A L C U L A T I O N S energy is used as a criterion of accuracy then the FLA is the best integral approxima- tion as it gives a total energy within 0.05 % of the McLean value. However the PLA calculation reproduces the dipole moment more satisfactorily. The accuracy of each calculation incorporating an integral approximation is further checked by calculating the one-electron energies. These energies listed in table 2 again show the superiority of the Lowdin approximations over the Mulliken approximation. The occupied molecular orbital energies are apparently better reproduced by the PLA calculation.TABLE 1 .-CALCULATED MOLECULAR PROPERTIES FOR HCN USING BEST ATOM ATOMIC ORBITALS IMA MA PLA FLA McLean 1 electronic energy - 116,6159 - 116.6079 - 11 6.5314 - 116.4692 - 116.4236 total energy E - 92.7397 - 92.7317 - 92.6552 - 92.5930 - 923474 kinetic energy T 91.2547 91.2412 91.4292 91.4905 91.3209 potential energy Y - 183.9944 - 183.9728 - 184.0843 - 184-0834 - 183.8683 -2T/V 0-9919 0.99 19 0.9933 0.9940 0.9933 dipole moment 1.309 1.304 2.033 1.803 2.100 All energies in a.u. (1 a.u. = 27.20963 eV). Dipole moments are in Debyes (1 D = 0.33356 x lo-'' C m) TABLE 2.-cOMPARISON OF ONE-ELECTRON ENERGIES FOR HCN (BEST ATOM ATOMIC ORBITAL BASIS SET) IMA MA PLA FLA McLean * 10 - 15.7427 - 15.7469 - 15.7287 - 15.7479 - 15.7402 20 -11.3972 -11.3918 -11.4035 -11.4142 -11.4277 30 - 1.2404 - 1.2430 - 1.2526 - 1.2575 - 1.2522 40 - 0.8263 - 0.8158 -0.7934 - 0.781 1 - 0.7965 50 -0.5691 - 0.5646 - 0.5644 - 0.5735 - 0.5582 In - 0.4960 - 0.4974 - 0.5030 - 0.51 37 - 0.5074 60 0.8400 0.6860 0.5215 0.3535 0.3648 7a 1.0923 1.1502 1.0799 0.961 3 1.0860 2rG 0.2592 0-2603 0.2557 0.2453 0-25 16 5o-ln gap 0.073 1 0.0672 0.0615 0.0598 0.0508 1n-271 gap 0.7551 0.7577 0.7586 0.7590 0.7590 1st I.P.(eV) 13.49 13.53 13.69 13.98 13.81 All energies are in a.u. But overall since the 60 level is particularly sensitive to the choice of integral ap- proximation the FLA is considered as the best choice. Two orbital energy gaps are also given in table 2 as well as the estimated value for the molecular ionization potential deduced by application of Koopmans' theorem. A Mulliken population analysis gives the results listed in table 3.The Lowdin approximations again yield results which are more in accord with those of McLean. The comparison between the dipole moments or effective atomic charges although useful does not reveal in which regions of space the wave function is inaccurate. G . DOGGETT AND A . MCKENDRICK 35 Some information on this problem is forthcoming from an examination of the differ- ence density function AP obtained by subtracting the approximate density function from the one given by McLean.l This Ahp function is shown in fig. 1 for the FLA calculation. The inaccuracies are mainly confined to regions close to the nuclei. A similar situation obtains for the other approximations but as expected the regions of inaccuracy are more extended in the MA calculation.Fig. 1 also shows that the overall effect of the multicentre integral approximation is to remove charge from the bonding regions and deposit it in the vicinity of each nucleus this behaviour is also apparent in the AP plots for the other integral approximations. The PLA or FLA gives results which are very close to the accurate ab irzitio values the MA on the other hand is poor. TABLE 3.-EFFECTIVE ATOMIC CHARGES DEDUCED FROM A MULLIKEN POPULATION ANALYSIS OF THE HCN MOLECULAR WAVE FUNCTION (BEST ATOMIC ORBITAL BASIS SET) N C H IMA + 0-083 - 0.261 +0*178 MA + 0.09 1 - 0.300 + 0.209 PLA - 0.100 -0.181 + 0.281 FLA -0.113 - 0.1 57 + 0.270 McLean - 0.08 1 -0.155 + 0.236 \ '\O.O i0.0 \ 0.0 N C H - 0.6k FIG. 1.-The upper part of the figure shows contours of the difference density function AP obtained by subtracting the FLA one-electron density function for HCN from that of McLean the basis functions for both the FLA and McLean calculations are best atom atomic orbitals.The lower part of the figure shows the variation of AP along the internuclear axis. The units of AP are electrons (a. u .)-" The success of the FLA and PLA indicates that a further study of the cyanide systems would be profitable. The results of some IMA and FLA calculations on HCN using a double-zeta basis a CN internuclear separation of 0.1159 nm (2.190 a.u.) and five different CH internuclear separations are shown in fig. 2. The IMA and FLA calculations yield minimum total energies of - 92.8970 and - 92.7449 a.u. respectively at CH internuclear separations of 0.1041 nm (1.968 a.u.) and 0.1 104 nm(2.086 a.u.) respectively so straddling theexperimental value of 0.1063 nm (2.009 a.u.).The molecular dipole moments are 1.75 and 2.59D respectively as against the experimental value * of 3.00D. But the assumed CN internuclear 36 INTEGRAL APPROXIMATIONS I N SCF CALCULATIONS separation is 0.0004 nm longer than the experimental value,8 to facilitate comparison with the FCN calculations in the next section. have given the results of a double-zeta calculation on HCN using RCH = 0.1066 nm (2.014 a.u.) and RCN = 0.1153 nm (2.179 a.u.). Since this molecular configuration is close to the IMA and FLA minimum energy configurations it is interesting to find the IMA and FLA total energies straddling the McLean and Yoshimine value of -92.8369 a.u. McLean and Yoshimine -92.710- -92.730- - 92.750 ? ? - I G - 9 2 .8 6 0 - -92.880. -92.900- 1 I . I 1-6 1-8 2.0 2.2 2'4 RCH(a.u-) FIG. 2.-A plot of the total energy E(a.u.) against RcH(a.U.) for the IMA (lower curve) and FLA (upper curve) calculations on HCN using a basis of best atom atomic orbitals. The one-electron energies follow the same general trends as before. In both approximate ab initio calculations all energy levels except 60 and 7a become more strongly binding as RCH decreases. Also the 1n-27-c gap is virtually constant over the range of RCH values. The predicted ionization potentials are now 15.43 and 15.67 eV respectively for the IMA and FLA calculations. A Mulliken population analysis of the molecular wave function gives rather similar results to those found earlier. In the IMA calculation the effective charge on H remains approximately constant at about +0-24 and as RCH decreases n charge is transferred preferentially to C.This 7-c transfer is accompanied by a smaller transfer of a charge from C to N. The results for the FLA calculation cannot be generalized so readily. As RCH decreases charge is transferred from H and redistri- buted predominantly on C while N increases its share of charge to a smaller extent. A more satisfactory representation of the way in which the electronic charge is redistributed on molecule formation is obtained by plotting contours of the difference density function : G . DOGGETT AND A. MCKENDRICK 37 The neutral atoms placed at the appropriate nuclear sites in the molecule are assumed to have the following configurations H Is C ls2 2s2 2px 2py N ls2 2s2 2px 2py 2pz where the 2p atomic orbitals are directed along the internuclear axis.Contours of the AP function are shown in fig. 3 for the FLA calculation at the experimental value of RCH. Electronic charge is displaced from the vicinity of the nuclei and accumu- lated in the bonding regions as well as the region usually ascribed to the N lone pair. Furthermore as the CH internuclear separation is decreased less charge is removed from the vicinity of the proton in contradiction to the results obtained from the Mulliken population analysis. N C H -0.6 FIG. 3.-The upper part of the figure shows contours of the difference density function AP for HCN obtained from the difference between the FLA and free atom one-electron density functions situated at their respective nuclei.The basis functions for the calculation are double-zeta atomic orbitals. The lower part of the figure shows the variation of AP along the internuclear axis. The units of AP are electrons (a.~.)-~. 3. THE ELECTRONIC STRUCTURE OF FCN The calculations on FCN are performed in a similar manner to those described for HCN. A double-zeta ' basis set is chosen and RCN is again fixed at 0.1 159 nm which is now the experimental value.' IMA and FLA calculations are performed using five different FC internuclear separations. In addition one PLA calculation is performed at the experimental geometry (see later). The IMA calculations yield a minimum energy of - 192-2624 a.u. at RFC = 0.1 196 nm (2.26 a.u.). The corresponding value for the dipole moment is 1-66 D (expt.,' 2.17 D).The apparent anomaly in the FLA calculation at the experimental value of RFC requires further investigation. The results of a PLA calculation at the same molecular geometry yield a total energy of - 191.6657 a.u. which is reasonable in view of the earlier cal- culations on HCN. McLean and Yoshimine find a total energy of - 191-6275 a.u. (no other molecular properties are quoted) for a double-zeta calculation on FCN using a slightly different molecular geometry RFC = 0.1260 nm (2.381 a.u.) and RCN = 0.1 165 nm (2.202 a.u.). If the anomalous FLA calculation is disregarded the estimated total energy and equilibrium value for RFC are - 191.566 a.u. and 0.134 nm (2.54 a.u.) respectively. Hence the IMA and FLA energies probably straddle the ab initio energy as was found in the earlier double-zeta calculation on The variation of total energy with RFC is shown in fig.4. 38 INTEGRAL APPROXIMATIONS I N SCF CALCULATIONS HCN. The IMA and FLA calculations also predict equilibrium values for RFC which straddle the experimental value of 0-1262 nm (2.385 a.u.). However the PLA may in fact be superior to the FLA as not only is the total energy better but so is the molecular dipole moment 1-68D as against an estimated 1.54D for the FLA calculation. The question as to why the FLA calculation breaks down should be answered. The determination of A2 involves the calculation of a ratio of two numbers of which the denominator is given as the difference of two summations. This difference is usually not small but for the particular choice of molecular parameters under con- sideration it is small enough to cause appreciable instability in the value of A2 e.g.A2 = 5.97 for the 2pz 2p carbon-fluorine overlap density. That this instabilityis the I 191.650 I I 1 - 192.1 5 0 - 192.250 ia I I 0 1 I I I I -T- I I I I . 2.0 2.2 2.4 2.6 2.8 RFC(a.u.) FIG. &A. plot of the total energy E(a.u.) against RFc(a.U.) for the IMA (lower curve) and FLA (upper curve) calculations on FCN using a basis set of double-zeta atomic orbitals. exception rather than the rule is confirmed to some degree by the reasonable values obtained for A2 at the other FC internuclear separations. So far it has been impossible to find any rules to predict when instability will occur. The best method of avoiding the instability is to examine the values of A at various molecular geometries before performing the SCF calculations.By this means any configuration leading to instability can be avoided. A comparison of the one-electron energies for the IMA and PLA calculations at the experimental molecular configuration is given in table 4. Also included are the results of McLean and Yoshimine lo for an extended basis set calculation incor- porating 3d and 4fatomic orbitals on each centre. Although their molecular con- figuration differs slightly from the one used here it is satisfying to find the PLA results for the occupied molecular orbital energies are in such good agreement with those of McLean and Yoshimine. However the IMA calculation predicts a different sequence of occupied energy levels with the 70 and 2n transposed with respect to the other two calculations. G .DOGGETT AND A . MCKENDRICK 39 The charge distribution in FCN is examined by means of Mulliken population analyses and maps of appropriate difference density functions. The population analyses give the following results for the effective atomic charges -0.183 + 0.358 -0.175 F C N -0.123 + 0.263 - 0.140 where the IMA and PLA values are given above and below each atom respectively. These charges are conventionally interpreted in terms of a displacement of charge TABLE 4.-A COMPARISON OF ONE-ELECTRON ENERGIES FOR THE IMA AND PLA CALCULATIONS OF FCN USING A BASIS SET OF DOUBLE-ZETA-ATOMIC ORBITALS 1 0 2 0 30 40 50 60 70 171 271 IMA - 26.413 1 - 15.7950 - 11 *6723 - 1.7721 - 1.41 12 - 0.9192 -0.6610 - 0.971 6 - 0.723 1 PLA - 26.4593 - 15.7749 - 11.6230 - 1.7679 - 1.3652 - 0.9391 - 0.6581 -0.8158 - 0.5618 McLean and Yoshimine 10 - 26.4375 - 15.6115 - 11.4018 - 1.7673 - 1.2523 - 0.9308 -0.6012 -0.8141 - 0.4984 80 0.5473 0.473 1 0,1922 90 1.4127 1.3190 - 3lT 0.5089 0,1619 0.2359 All energies are in a.u.The McLean and Yoshimine results have been obtained from a calculation using an extended basis set and a slightly different molecular geometry (see text). f 0.0 FIG. 5.-The upper part of the figure shows contours of the difference density function AP for FCN obtained from the difference between the FLA and free atom one-electron density functions situated at their respective nuclei. The basis functions for the calculation are double-zeta atomic orbitals. The lower part of the figure shows the variation of AP along the internuclear axis. The units of AP are electrons (a.~.)-~.40 INTEGRAL APPROXIMATIONS I N SCF CALCULATIONS from C to each end of the molecule. The difference density function obtained by subtracting the appropriate atomic densities from the FLA molecular density function is shown in fig. 5. Again the pattern is similar to that found for HCN charge is displaced from the regions around the nuclei and accumulated in the bonds or in regions normally considered as being occupied by N or F lone pairs. The plot of AP along the molecular axis shows this effect very clearly. 4. CONCLUSION The systematic use of three integral approximations in LCAO MO SCF cal- culations on HCN and FCN has been examined. Preliminary calculations on HCN using the same molecular parameters as McLean show that the two Lowdin approxi- mations are superior to the Mulliken approximation.Further calculations on HCN and FCN using a double-zeta basis set of atomic orbitals yield total energies which straddle the ab initio value. The FLA appears to be breaking down in one of the calculations on FCN. The PLA devoid of such instability provides a simpler yet reliable approximation for reproducing the results of an ab initio calculation with a surprising degree of accuracy. One of us A. McK. thanks the S.R.C. for the award of a Research Studentship. A. D. McLean J. Chem. Physics 1962,37 627. R. S. Mulliken J. Chim. Physique 1949 46,497 675. J. A. Pople D. P. Santry and G. A. Segal J. Chem. Physics 1965,43 S129. P. 0. Lowdin J. Chem. Physics 1953 21 374. P. 0. Lowdin Adv. Physics 1956 5 1. R. McWeeny Rev. Mod. Physics 1960 32 335. E. Clementi Tables ofAtomic Wave Functions (I.B.M. J. Res. Develop. Suppl. 1965). J. Sheridan and J. K. Tyler Trans. Faraduy Soc. 1963 59 2661. A. D. McLean and M. Yoshimine Int. J. Quantum Chem. 1967 1% 313. Develop. Suppl. 1968). lo A. D. McLen and M. Yoshimine Tables of Linear Molecule Wave Functions (I.N.M. J . Res.a