Journal of Molecular Structure (Theochem) 709 (2004) 15–23 www.elsevier.com/locate/theochem
Sturmian orbitals and molecular structure Vincenzo Aquilanti*, Andrea Caligiana Dipartimento di Chimica, Lab. Di Chimica, Universita` di Perugia, Via Elce di Sotto No. 8, 06123 Perugia, Italy
Abstract In this paper, we introduce the basic steps of a computational scheme for molecular structure problems, where Sturmian orbitals appear as building blocks in a valence-bond type of approach. The Sturmian basis set approach to the solution of quantum chemical problems is the configuration space counterpart of hyperspherical analysis in momentum space and leads to promising algorithms alternative to the usual SCF-plus-configuration interaction hierarchies. The mathematical ingredients are considered, particularly the possibility of exploiting relationships with Slater-type orbitals in order to circumvent the bottleneck of calculation of matrix elements. We also propose a procedure to solve two-center overlap integrals between Sturmians based on a momentum-space approach. The example of the electronic structure of the hydrogen molecule is worked out and results of computations are given. The optimal choice of scaling parameters for Sturmian sets is also discussed. q 2004 Published by Elsevier B.V. Keywords: Sturmian sets; Hydrogonoid orbitals; Slater-typeorbitals; Momentum-space; Diatomic molecules
1. Introduction The Schro¨dinger equation (SE) for the non-relativistic hydrogen atom is analytically soluble in a closed form, and therefore the use of hydrogenoid orbitals as expansion basis sets to calculate the electronic wave functions of atoms and molecules has been extensively considered since the early days of quantum mechanics. Sturmian orbitals of the hydrogenoid type have been classified [1] and parabolic [2] and elliptic [3,4] variants studied, but the most usual ones are those corresponding to the familiar separation of variables in polar coordinates. Since the early days of quantum chemistry, the closely related Slater-type orbitals (STO) and other variants have been proposed and amply tested as alternatives. However, in order to tackle the obstacles encountered for the efficient computation of matrix elements (specifically the bielectronic integrals), preference has been given over the years to the Gaussian orbitals (GTO), particularly convenient for such purposes in spite of their slower convergency rate. * Corresponding author. Fax: C39-075-585-5606. E-mail address:
[email protected] (V. Aquilanti). 0166-1280/$ - see front matter q 2004 Published by Elsevier B.V. doi:10.1016/j.theochem.2003.10.070
The current software packages available for ab initio calculations of atomic and molecular properties, typically exploit Gaussian basis set expansions: the procedures almost invariably initiate with Self Consistent Field (SCF) or Multiconfiguration Self Consistent Field (MC-SCF) calculations yielding molecular orbitals as linear combinations of GTO centered on the various nuclei of the molecule, and then proceed with post-Hartree–Fock techniques (Configuration Interaction, Coupled Clusters, etc.). Among the disadvantages of such procedures, we mention the difficulties encountered in the calculation of long-range interactions and of the properties of excited states, and the slow convergence in reproducing the cusps in the electronic densities at nuclei. The Sturmian orbitals can be straightforwardly expressed as linear combinations of Slater-type orbitals. Although, as we have seen, the use of the latter has been largely superseded by the GTO, for which the slower convergence is often amply balanced by the much simpler algorithm for the calculation of integrals, it is important to note that the employed techniques are practically identical. On the contrary, the approach involving Sturmian orbitals will be shown here to proceed along a route fully alternative to modern ab initio quantum chemistry: it will be seen that
16
V. Aquilanti, A. Caligiana / Journal of Molecular Structure (Theochem) 709 (2004) 15–23
there is no need of iterative Hartree–Fock calculations, because the introduction of the Sturmian expansion of the electronic wave function leads after algebraic manipulations to a generalized eigenvalue problem (GEP) where the matrix elements contain mono- and bi-electronic integrals among Sturmian orbitals and the eigenvalue is not the electronic energy E but the momentum wave number p0 such that EZKp20 =2 (atomic units are implied). The typical integrals of quantum chemistry depend explicitly on the nuclear coordinates Ri, while parameters like siZp0Ri appear in the Sturmian approach. The basic feature is, therefore, that the Born–Oppenheimer separation is not settled before the calculations. The si parameters are instead fixed and only after the solution of the GEP from the relationships RiZsi/p0 one finds out what is the point in the potential energy surface which corresponds to an energy Kp20 =2: Another unusual property of this approach is that the eigenvalues obtained by solving the GEP come out as correct solutions not only for ground states, but also for the excited states corresponding to the same symmetry. The approach to atomic structure as a many-electron onenucleus problem has been extensively developed in recent years. In Appendix A, the reader will find the procedure of building a polyelectronic basis set as a product of ‘monoelectronic’ Sturmian orbitals [5–11]: an alternative method, which defines N-electronic Sturmians as solutions of the 3N-dimensional hydrogen atom can be elegantly formulated and also has a momentum space counterpart in terms of hyperspherical harmonics [12–14], but still lacks of developments in the algebraic and computational implementation in order to be effectively competitive. Applications so far have been limited to atomic properties and to the one-electron molecule HC 2 (see Refs. [2,15] and references therein). We next consider the application of the Sturmian orbital procedure to the ab initio quantum chemistry of molecules as many-electron many-nuclear systems. The recent proposal, as reviewed in Avery’s book [16] (see also Ref. [17]), involves a MOLCAO approach where one-electron molecular orbitals are computed first, and then their products are employed to build-up N-electron configurations. Such an approach is currently being implemented, one difficulty being the recovery of the Born–Oppenheimer separation, because, as indicated above, the matching of si parameters and Ri coordinates occurs through scaling by the eigenvalues p0, and a fixed si calculation provides basis sets for different physical Ri for different orbitals. This paper presents an account of a promising method for the application of the Sturmian orbital approach to molecular structure calculations. The new proposal (Section 2) for the extension of generalized Sturmian sets to molecular structure is inspired by concepts and formalisms previously elaborated. In the framework of quantum chemistry jargon, it will be recognized as a valence-bond type of approach. Section 3 is devoted to a discussion of the practical obstacles to the implementation of this approach,
the calculation of integrals, which can be in part circumvented through the angular momentum algebra applicable in the momentum space perspective, and in part tackled by exploiting the connections with STO and the existence of programs dealing with integrals involving the latter. Results and conclusions are in Section 4, where it is shown that the calculation of the potential energy curve for the H2 molecule provides encouraging results already at the level of minimal basis set. Appendix A further contains some considerations about the choice of the ‘energy’ in the equation that defines the Sturmians basis sets.
2. Many-center many-electron systems In this section, we propose an extension of the Sturmian formalism for the study of systems having many nuclei and electrons. Background information on the atomic case is relegated to Appendix A. The SE for a N-electronic nonrelativistic molecule is 1 p2 K V2 C V0 ðxÞ C V 0 ðxÞ C 0 JðxÞ Z 0 (1) 2 2 where V2 Z
N X
V2j
(2)
jZ1
V0 ðxÞ Z
X
V0;a ðxÞ Z
a
N XX a
iZ1
Za jxi K Xa j
(3)
and V 0 ðxÞ Z
N X 1 r iOj ij
(4)
The energy of the system is indicated as Kp20 =2; and x is the set of all the electronic coordinates {x1,x2,.,xN}. The index a is referred to the nuclei of the molecule. Corresponding formulas for the atomic case are given in Appendix A. We Pexpand the electronic wavefunction as JðxÞZ n Cn;p Fn;p ðxÞ; where basis functions Fn,p(x) are defined by the differential equation: 1 p2 K V2 C bn;p Vp‡ ðxÞ C 0 Fn;p ðxÞ Z 0 (5) 2 2 The index n labels the different solution of the equation; in the following the meaning of p is explained. The potential Vp‡ ðxÞ is the sum of N contributions, each of them representing the interaction of one electron with one of the nuclei: Vp‡ ðxÞ Z K
N X jZ1
Zp;j jxj K Xp;j j
(6)
V. Aquilanti, A. Caligiana / Journal of Molecular Structure (Theochem) 709 (2004) 15–23
This means that actually there is not only one equation defining the basis set, because there are many choices for Vp‡ ðxÞ : the basis set can be divided into several subsets, each of them associated to an equation like (5) having a different Vp‡ ðxÞ: The index p plays the role of label for the various subsets. As expansion basis sets alternative Coulomb Sturmians are known and have been recently classified by us (see Section 4). The commonly used are the polar Sturmian orbitals, which are defined as: ðn K l K 1Þ! 1=2 unlm ðr; q; fÞ Z ð2km Þ3=2 2n½ðn C lÞ!3 2lC1 !eKkm r ð2km rÞl LnKlK1 ð2km rÞYlm ðq; fÞ
(7)
The L2lC1 nKlK1 are Laguerre polynomials as in Ref. [18] and Ylm is a spherical harmonic. A peculiarity of such a basis set is that the ‘energy’ Kkm2 =2 cannot be positive, and so the set is complete without any continuum states. The Sturmian orbitals obey a potential weighted orthonormality relation: ð 1 1 1 dx un 0 l0 m 0 ðxÞ unlm ðxÞ Z dn0 n dl 0 l dm 0 m (8) km r n and are normalized to unity: Ð dx unlm ðxÞunlm ðxÞ Z 1
(9)
The functions Fn,p(x) are determinants made of polar Sturmian orbitals, provided that two subsidiary conditions are satisfied: km2 C km2 0 C km2 00 C/C Z p20
(10)
In this case the basis set is constituted by two bielectronic functions, F1(x)Z1sA(1)1sB(2) corresponds to the potential V1‡ ðxÞ Z K
V2‡ ðxÞ Z K
1 J1Sg ðx1 ; x2 Þ Z pffiffiffi ½1sA ð1Þ2sB ð2Þ C 1sB ð1Þ1sA ð2Þ 2 1 ! pffiffiffi ½að1Þbð2Þ K bð1Það2Þ 2
km3 1sð1Þ Z p
(16)
1=2
eKkm r1
(17)
The exponentpkffiffimffi can be computed combining Eqs. (11) and (12): it is 1= 2: In general, the exponent of an orbital will depend on the principal quantum numbers of all the orbitals forming the configuration Fn(x) in which it appears. Now we will calculate an approximation to H2 electronic energy exploiting the wave function (13). Entering Eq. (13) inside Eq. (1) and writing the potential V0(x) as in Eq. (16) one can obtain the 1!1 secular problem: ð
½1sA ð1Þ1sB ð2Þ C 1sB ð1Þ1sA ð2Þ 1 1 p2 ! K V2 C V1‡ ðxÞ C V2‡ ðxÞ C C 0 2 r12 2 !1sA ð1Þ1sB ð2Þd1d2 Z 0
(18)
The latter expression can be simplified with the aid of Eq. (5): ½1sA ð1Þ1sB ð2Þ C 1sB ð1Þ1sA ð2Þ !
ð1 K b1 ÞV1‡ ðxÞ C V2‡ ðxÞ C
1 1sA ð1Þ1sB ð2Þd1d2 Z 0 r12 (19)
Taking into account that the basis set functions associated to the same operator Vp‡ ðxÞ enjoy potential weighted orthonormality relationships ð
(13)
(15)
According to Eq. (7), 1s Sturmian orbitals are defined as
ð
This approach can be interpreted as a Valence Bond method, since we perform CI over configurations which are products of non-orthogonal atomic orbitals. In order to clarify the concepts exposed and show how this method works in practice, we will show an application to the simplest molecule H2. The crudest approximation to the ground state 1 Sg wavefunction is:
1 1 K jx1 K XB j jx2 K XA j
V0 ðxÞ Z V1‡ ðxÞ C V2‡ ðxÞ
(11)
(12)
(14)
We can express the overall interaction between electrons and nuclei as:
A combination of Eqs. (10) and (11) leads to the relationship: p0 bn;p Z rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
2 Z 2 Za 0 2 Za 00 a C C C/ 0 n n n 00
1 1 K jx1 K XA j jx2 K XB j
in Eq. (5), and F2(x)Z1sB(1)1sA(2) corresponds to
and nkm n 0 km 0 n 00 km 00 Z Z Z/Z bn;p Za 0 Za 0 Za 00
17
dxFn0 ðxÞVp‡ ðxÞFn ðxÞ Z dn0 ;n
2E bn
(20)
18
V. Aquilanti, A. Caligiana / Journal of Molecular Structure (Theochem) 709 (2004) 15–23
we can modify Eq. (19) as follows: ð 2E 1 ‡ K 2E C 1sA ð1Þ1sB ð2Þ V2 ðxÞ C 1sA ð1Þ1sB ð2Þd1d2 b1 r12 ð 1 C 1sB ð1Þ1sA ð2Þ ð1 K bÞV1‡ ðxÞ C V2‡ ðxÞ C r12
elements and then we can get the eigenvalue p0; eventually, pffiffiffi the physical internuclear distance is obtained as RZ 2s=p0 :
3. The integrals. Relationships with Slater-type orbitals
!1sA ð1Þ1sB ð2Þdx Z 0
The evaluation of integrals is always an important issue of quantum chemical methods, because it can involve large (21) efforts in terms of computer times and therefore may limit the concrete possibility of extending the calculations to If we write down all the integrals explicitly and exploit cover complex systems. This discussion will pertain Eq. (12), we can rewrite Eq. (21) as: ð ð 2 ð ð ð pffiffiffi 1s2A ð1Þ 1s2B ð2Þ 1sA ð1Þ1s2B ð2Þ d1 K 1s2A ð1Þd1 d2 C K 2p0 C p0 K 1s2B ð2Þd2 d1d2 jx1 K XB j jx2 K XA j r12
ð ð ð ð p0 1sA ð1Þ1sB ð1Þ 1sA ð2Þ1sB ð2Þ p ffiffi ffi d1 K 1sA ð1Þ1sB ð1Þd1 d2 C 1K K 1sA ð2Þ1sB ð2Þd2 jx1 K XA j jx2 K XB j 2 ð ð ð ð ð 1sA ð1Þ1sB ð1Þ 1sA ð2Þ1sB ð2Þ 1sA ð1Þ1sB ð1Þ1sA ð2Þ1sB ð2Þ C K 1sA ð2Þ1sB ð2Þd2 d1 K 1sA ð1Þ1sB ð1Þd1 d2 C d1d2 jx1 K XB j jx2 K XA j r12 Z0
(22)
Eventually, we divide the latter by p0 and write the NAI’s in a more compact form, exploiting the definitions ð 1 1 0 0 0 u ðx K RÞ dx un 0 l0 m0 ðxÞ Snlm ðRÞ h (23) nlm km jx K Rj nlm Wnlm n 0 l 0 m 0 ðRÞ h
1 km
ð dx
un0 l 0 m0 ðxÞ
1 u ðxÞ jx K Rj nlm
Ð 0 0 0 Mnnlml m h un 0 l0 m0 ðxÞunlm ðxÞdx Then the expression for the eigenvalue is: hpffiffiffi pffiffiffi p0 Z 2 C 2W100 100 ðsÞ ð 1 1sA ð1Þ1sB ð1Þ1sA ð2Þ1sB ð2Þ C d1d2 p0 r12 ð 2 1 1sA ð1Þ1s2B ð2Þ C d1d2 p0 r12 i pffiffiffi 100 100 100 K2 2M100 100 ðsÞS100 ðsÞ =½1 C M100 ðsÞS100 ðsÞ
specifically to the mono- and bi-electronic integrals encountered in Section 2. 3.1. Analytical formula for M
(24)
Let us consider the overlap integral between two polar Sturmians having the same exponent and centered on two different nuclei a and b, separated by a distance R:
(25)
Ð 0 0 0 lm Mb;n a;nlm h ua;n 0 l 0 m 0 ðxÞub;nlm ðxÞdx
(27)
Working in momentum space, the latter is equivalent to: Ð
Ð ja;n0 l 0 m 0 ðpÞjb;nlm ðpÞdp Z eKip,R jnlm ðpÞjn0 l 0 m0 ðpÞdp
(28)
We are using the symbol jnlm(p) for Sturmians in momentum space. Inserting in Eq. (28) the plane wave expansion [19] ð26Þ
The M’s are two-center overlap integrals between Sturmians: The integrals S and W were introduced by Shibuya and Wulfman [19], who solved the former in a momentum space framework. In Section 3 we show that M integrals can be evaluated in momentum space exploiting angular momentum algebra, analogously to the case of the S and W integrals, which we have discussed recently from this same perspective [15]. Mono- and bi-electronic integrals in Eq. (26) actually do not depend on the value of p0: rather, they depend on the product sZkm$R (R being the internuclear distance). Therefore, we pick up a value for s, compute the matrix
expðip,RÞ Z ð2pÞ3=2
X n;l;m
un;l;m ðRÞ
p20 C p2 jn;l;m ðpÞ 2p20
(29)
the relationship between jnlm(p) and hyperspherical harmonics [14,16,20], jn;l;m ðpÞ Z
4p5=2 0 Yn;l;m ðc; w; 4Þ C p2 Þ 2
ðp20
and the volume element in p-space:
2 3 k m C p2 dp Z dU 2km
(30)
(31)
V. Aquilanti, A. Caligiana / Journal of Molecular Structure (Theochem) 709 (2004) 15–23
we get an integral involving hyperspherical harmonics: X 0 0 0 lm 3=2 Mb;n un00 l 00 m00 ðRÞ a;nlm Z ð2pÞ n 00 l 00 m 00
ð
4km5=2 ! dU 2 Ynlm ðUÞYn0 l 0 m0 ðUÞYn00 l 00 m00 ðUÞ ðp C km2 Þ2 (32) One can easily transform the fourth relationship in Eq. (3) of Ref. [4] into 4km5=2 ð1 C cos cÞ Z 2 km3=2 ðkm C p2 Þ2 2
0 0
lm Mb;n a;nlm
0
19
From a three-term recurrence relationship involving Gegenbauer polynomials [22], we have obtained the following expression: cos cYNLM ðUÞ 1 ZK 2
rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðN K LÞðN C L C 1Þ 1 YNC1;LM ðUÞ K NðN C 1Þ 2
rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðN C LÞðN K L K 1Þ YNK1;LM ðUÞ ! NðN K 1Þ
(36)
(33)
Inserting Eqs. (35) and (36) into Eq. (34), we get a solution for the integral: 9 8 0 n K1 n K1 N K1 > > > > > > >
3=2 X 0 = < 2 2 2 > 2p nn Nð2l C 1Þð2l 0 C 1Þ 1=2 0 0 0 Z hl ; m ; l; mjL; Mi n K 1 n K 1 N K 1 > > km 2p2 > NLM > > 2 2 > > > ; : 20 l l L rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðN K LÞðN C L C 1Þ ðN C LÞðN K L K 1Þ ðN K LÞðN C L C 1Þ ! 1C ð1 K dN;1 Þ uNLM ðRÞ K uNC1;LM ðRÞ C 4NðN C 1Þ 4NðN K 1Þ NðN C 1Þ rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðN C LÞðN K L K 1Þ 1 uNK1;LM ðRÞð1 K dN;1 Þ C K NðN K 1Þ 4ðN C 1Þ sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðN K LÞðN C L C 2Þ½ðN C 1Þ2 K L2 1 ! uNC2;LM ðRÞ C NðN C 2Þ 4ðN K 1Þ 9 sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi = 2 2 ðN C LÞðN K L K 2Þ½ðN K 1Þ K L ! uNK2;LM ðRÞð1 K dN;1 Þð1 K dN;2 Þ ; NðN K 2Þ (37)
Then the integral becomes:
3=2 X ð 0 0 0 2p lm 00 00 00 Mb;n Z u ðRÞ dUð1 n l m a;nlm km n 00 l 00 m 00 C 2cos c C cos2 cÞYnlm ðUÞYn0 l 0 m0 ðUÞYn00 l 00 m00 ðUÞ
(34)
The product of two hyperspherical harmonics can be expressed as (Eq. (2.32) in Ref. [21]): X nn 0 ð2l C 1Þð2l 0 C 1Þ 1=2 Ynlm ðUÞYn0 l 0 m 0 ðUÞ Z 2p2 N;L;M 9 8 0 n K1 n K1 N K1 > > > > > > > > > 2 2 > = < 2 0 0 0 !hl ; m ; l; mjL; Mi n K 1 n K 1 N K 1 YNLM ðUÞ > > > > > 2 2 2 > > > > > ; : 0 l l L ð35Þ
Generally, the use of a basis set as defined in Eq. (5) leads to monoelectronic integrals S, W and M containing two orbitals which can differ in the exponents appearing in Eq. (7). This represents a limit for the application of momentum space harmonic analysis, at least according to its present developments. To circumvent this bottleneck and yet be able to test the perspective power of the method, we explored the possibility of making recourse to the relationship between Sturmian orbitals and STO. 3.2. Numerical evaluation Since Ferna´ndez Rico and collaborators [23,24] have proposed, implemented and made available an efficient program for the calculation of two-center monoelectronic integrals involving STO, based on the recurrence relationships which they enjoy, and since the Sturmian orbitals can be given as simple linear combinations of the STO (explicitly, indicating the latter as jnlm:
20
V. Aquilanti, A. Caligiana / Journal of Molecular Structure (Theochem) 709 (2004) 15–23
u1;0;0 Z j1;0;0 pffiffiffi u2;0;0 Z j1;0;0 K 3j2;0;0 u2;1;m Z j2;1;m pffiffiffi pffiffiffiffiffi u3;0;0 Z j1;0;0 K 2 3j2;0;0 C 10j3;0;0 rffiffiffi pffiffiffi 2 u3;1;m Z 2 j2;1;m K 5j3;1;m 3
(38)
u3;2;m Z j3;2;m pffiffiffiffiffi pffiffiffi pffiffiffiffiffi u4;0;0 Z j1;0;0 K 3 3j2;0;0 C 3 10j3;0;0 K 35j4;0;0 rffiffiffi pffiffiffi pffiffiffiffiffi 3 u4;1;m Z 5j2;1;m K 5 j3;1;m C 21j4;1;m 2 pffiffiffi 3 u4;2;m Z pffiffiffi j3;2;m K 7j4;2;m 2 etc.) we found it expedient to use the program in order to calculate the bicenter monoelectronic integrals S and M. The W integrals can instead be calculated analytically even though the exponents differ. However, this program does not calculate the threecenter monoelectronic integrals ð 1 dx unlm ðxÞ u 0 0 0 ðx K Xa0 Þ (39) jx K Xa j n l m which can occur in molecules with three or more nuclei. The most demanding effort for our valence-bond Sturmian approach appears, as usual in quantum chemistry, when the calculation of the bielectronic multicenter integrals is needed: ðð 1 dxdx 0 utA ðxÞutB ðxÞ u ðx 0 ÞutD ðx 0 Þ (40) jx K x 0 j tC Our test cases presented in Section 2 were tackled in a rather tortuous manner, making recourse to another program due to Ferna´ndez Rico and coworkers [31] which in turn in order to obtain integrals for STO expands them as linear combinations of GTO. The resulting Gaussian integrals are finally calculated using available techniques. The ensuing disadvantage of this inelegant procedure is that, since a satisfactory fit of an STO requires about 10 Gaussians, an integral involving four STO’s requires the calculation of 104 integrals among GTO. Moreover our integrals such as Eq. (40) are linear combinations of integrals among STO: for example an integral involving four 3s-Sturmians requires a sum over 34 integrals of the STO program, and therefore 34! 104 integrals over Gaussians. Therefore, this program, which recommends itself for the accuracy which can be adjusted due to flexibility in the choice of the number Gaussians required to approximate a STO and satisfactorily covers the range of permitted quantum numbers (nmaxZ6, lmaxZ3), although appropriate for test calculations, must be superseded by more direct approaches to tackle integrals such as
Fig. 1. Comparison between the benchmark ab initio potential energy curve for H2 denoted as (a) [25], and the Sturmian approach. Curve (c) corresponds to the single configuration function in Eq. (13) and curve (b) has been obtained including in the expansion all configurations from all Sturmian orbitals up to 3d. In the latter case, the variance in the energy minimum with the one given in Ref. [25] is C0.3%. Numerical values are given in Tables 1 and 2.
Eq. (40). Alternatives should also be explored, involving, for example, elliptic or other suitable coordinates [32]. 4. Exemplary results, further remarks and conclusions Results on the potential energy curve for the ground state of H2, calculated from the procedure described in Section 2, are reported in Fig. 1. Even if really minimal, the basis set of Eq. (13) gives a good qualitative picture of the chemical Table 1 Potential energy curve for the H2 molecule, obtained by the procedure described in Section 2 with a basis set made of only one function (Eq. (13)) R (a.u.)
V(R) (a.u.)
1.0821753 1.2020791 1.3254692 1.4521964 1.5820863 1.7149417 1.8505432 1.9886483 2.1289914 2.2712856 2.4152266 2.560499 2.7067849 2.8537741 3.0011745 3.1487233 3.2961955 3.4434107 3.590235 3.7365797
K1.1423625 K1.1611956 K1.1694278 K1.1701995 K1.1657679 K1.1577794 K1.1474508 K1.1356935 K1.1231972 K1.1104875 K1.0979637 K1.0859263 K1.0745945 K1.0641206 K1.0545992 K1.046076 K1.0385553 K1.0320072 K1.0263758 K1.0215866
A graphical representation is given in Fig. 1(c).
V. Aquilanti, A. Caligiana / Journal of Molecular Structure (Theochem) 709 (2004) 15–23
bond: the equilibrium distance is close to the exact one, and the asymptotic behavior is correct. A more elaborate calculation is also shown in Fig. 1 and Tables 1 and 2. We have tested that even modest size calculations employing larger basis sets yield accurate results also for excited states. The account presented in this paper on developments regarding the use of Sturmians orbitals in quantum chemistry demonstrates how fruitful was the idea put forward long ago by Shull and Lo¨wdin of ‘natural spin orbitals’ [26]. The first steps of the proposed approach show perspectives of advantages over conventional methods, particularly for excited states. Not considered in this work have been the alternatives which open by exploiting the various coordinates systems to express Sturmian functions: they have been recently classified both in Table 2 Potential energy curve for H2 molecule, obtained by the procedure described in Section 2 with a basis set made of Sturmian orbitals up to 3d R (a.u.)
V(R) (a.u.)
0.7290924 0.80695213 0.88703841 0.96934499 1.0538422 1.1404765 1.2291706 1.3198246 1.412318 1.5065128 1.6022571 1.6993889 1.7977411 1.8971455 1.9974369 2.0984569 2.2000568 2.3021 2.4044637 2.5070395 2.6097341 2.7124688 2.8151786 2.9178118 3.020328 3.1226974 3.2248993 3.3269203 3.4287536 3.5303976 3.6318546 3.7331303 3.8342327 3.9351713 4.0359568 4.1366007 4.2371146 4.33751 4.4377982 4.5379902
K0.90468312 K0.97216731 K1.0204927 K1.0543036 K1.0770535 K1.0913646 K1.0992606 K1.102322 K1.1017926 K1.0986553 K1.0936868 K1.0874984 K1.0805681 K1.0732652 K1.0658709 K1.0585944 K1.0515868 K1.0449529 K1.0387603 K1.0330478 K1.0278317 K1.0231117 K1.0188752 K1.0151008 K1.0117613 K1.0088259 K1.0062621 K1.0040368 K1.0021176 K1.000473 K0.99907337 K0.99789094 K0.9969001 K0.99607738 K0.99540153 K0.99485338 K0.99441577 K0.99407342 K0.99381281 K0.99362199
A graphical representation is given in Fig. 1(b).
21
configuration [3] and in momentum space [4]. Since the computation of integrals is a feature which will eventually decide the crucial potentialities and scope of the Sturmian approach, the alternative basis sets and the relationships among them are an important goal to further study; some progress have been made for the parabolic set [2,13]. Needless to say, these orbitals might find a specific role also whenever in quantum chemistry there is a need for an expansion in a complete basis set. Further extensions which are of relevance for the application to heavier systems are those involving the inclusion of relativistic effects: see Chapter 8 in Ref. [16] and recent work by Szmytkowski [27–30]. For the particular scheme introduced in Section 2 the numerical results on H2 indicate that the development of a direct algorithm for the calculation of bielectronic multicenter integrals would greatly contribute to practical implementations to systems of increasing complexity. Our current effort on this problem shows the importance of using momentum space techniques, which allow the exploit of powerful angular and hyperangular momentum algebra. The work presented here can be appreciated in a wider context in a recent review article [33]. See also [34] for further applications and [35] for a review on the STO integral problem. Acknowledgements We thank Simonetta Cavalli, Cecilia Coletti and Gaia Grossi for their cooperation on various aspects of this work, and John Avery and Francesco Tarantelli for several inspiring discussions. This research is supported by the Ministero dell’ Istruzione, dell’ Universita` e della Ricerca (MIUR), the Ente Nazionale per l’Energia Alternativa (ENEA), the Agenzia Spaziale Italiana (ASI) and by European Contract ‘Theoretical Studies of Electronic and Dynamical Processes in Molecules and Clusters’ [Contract no. 1-1999-00121], and COST Actions D9 and D23. Appendix A. Atomic Sturmians and optimization of the ‘energy’ parameters In this appendix, we will provide basic formulas for the atomic case and will also discuss a technical issue, the optimal choice of the scaling parameter in Sturmian sets, which is of interest also in the context of the molecular applications considered in the main text. The SE for a N-electron free atom is 1 p2 K V2 C V0 ðxÞ C V 0 ðxÞ C 0 JðxÞ Z 0 (A1) 2 2 where V2 Z
N X jZ1
V2j
(A2)
22
V. Aquilanti, A. Caligiana / Journal of Molecular Structure (Theochem) 709 (2004) 15–23
V0 ðxÞ Z K
N X Z jZ1
(A3)
rj
and V 0 ðxÞ Z
N X 1 r iOj ij
(A4)
The electronic energy has been indicated as Kp20 =2; and the nuclear charge as Z. The wavefunction is expanded by a linear combination of Slater determinants X JðxÞ Z Cn Fn ðxÞ (A5) n
Our basis set is not completely defined yet, because we have to establish a rule for assigning a value to the exponents kj. We demand our basis functions Fn(x) to be a solution to the differential Equation [5,7,10,11,16] 1 p2 K V2 C bn V0 ðxÞ C 0 Fn ðxÞ Z 0 (A6) 2 2 which is the counterpart in the atomic case of Eq. (5) in the main text. The following is also similar. This is true if we impose two subsidiary conditions for the exponents k12 C k22 C .kN2 Z p20
show that actually using the basis set defined in Eq. (A6) does not imply any iteration: one has to calculate just one matrix, and diagonalize it. Substitution of expansion (A5) into Eq. (A1) yields: X 1 2 p2 K V C V0 ðxÞ C V 0 ðxÞ C 0 Cn Fn ðxÞ Z 0 (A11) 2 2 n Then, taking into account that the functions Fn(x) obey Eq. (A6), we can eliminate the kinetic energy operator from Eq. (A11): X ½Kbn V0 ðxÞ C V0 ðxÞ C V 0 ðxÞCn Fn ðxÞ Z 0 (A12) n
Multiplying Eq. (A12) by a function Fn0 ðxÞ and integrating over the coordinates, exploiting Eq. (A10) we obtain: XÐ dðxÞFn0 ðxÞ½V0 ðxÞ C V 0 ðxÞFn ðxÞ K 2Edn0 ;n Cn Z 0 n
ðA13Þ Then if we divide the latter equation by Kp0 and introduce the definitions ð 1 dðxÞFn0 ðxÞV0 ðxÞFn ðxÞ (A14) Tv00 ;n hK p0
(A7) Tv0 0 ;n hK
and n1 k1 Z n2 k2 Z . Z nN kN Z Zbn
2
ð
dðxÞFn0 ðxÞV 0 ðxÞFn ðxÞ
(A15)
(A8)
Eqs. (A7) and (A8) imply that the exponents’ value is: p0 ffi kj Z qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi (A9) nj n12 C n12 C/C n12 1
1 p0
N
Tv0 ;n hTv00 ;n C Tv0 0 ;n
(A16)
we arrive to a concise expression for our secular problem: X ½Tv0 ;n K p0 dv0 ;n Cn Z 0 (A17) n
The solutions of Eq. (A6) enjoy a potential weighted orthonormality relationship [11,16]: ð 2E dxFn0 ðxÞV0 ðxÞFn ðxÞ Z dn0 ;n (A10) bn The definition of the generalized N-electron Sturmian basis set Fn(x) is analogue to the definition of one-electron hydrogenic Sturmian functions: the energy ðKp20 =2Þ is constant for all the functions, while the ‘charge’ bn plays the role of an eigenvalue. We notice that the same Sturmian orbital may have different exponents if it appears in two different pconfigurffiffiffi ations. For example, in FnZ1sp2ffiffiffiffiffiffi one ffi has k1s Z p0 = 2; while in FnZ1s2s one has k1s Z p0 = 5=4: We have employed the unknown energy Kp20 =2 in the definition of our N-electron basis set. The exponents kj depend on the pvalue ffiffiffiffiffiffiffiffiffi of p0, i.e. on the value of the electronic energy ðp0 Z K2EÞ; consequently it looks like we are following an iterative approach, in which we have to guess an initial value for p0 and then iterate to get improved values for p0 until convergency is achieved. Now we are going to
The integral (A14) can be easily expressed taking into account Eqs. (A7), (A8) and (A10): 2E p 1 d0 Z 0d0 Z p0 bn v ;n bn v ;n bn ! ! X 2 1=2 X 1 1=2 Z kj dv 0 ;n Z Z dv 0 ;n n2 j3n j3n j
Tv00 ;n Z K
(A18)
The integral Tv00 ;n does not depend on p0’s value, and this is true also for the electronic repulsion integral T 0 n 0 ,n: in the radial part of the integral, one can perform a change of integration variable from r to p0r, so that p0’s value is not needed to solve the integral. In order to calculate energies and wavefunction for the electronic states of an isolated atom, one has just to compute and diagonalize the matrix Tn 0 ,n. We have chosen to use the energy of the system under investigation, Kp20 =2; in the definition of the basis set (Eq. (A6)). It would be interesting to understand whether this is the best possible choice for the definition of the basis set; perhaps a different value for the energy in Eq. (A6)
V. Aquilanti, A. Caligiana / Journal of Molecular Structure (Theochem) 709 (2004) 15–23
would give a basis set that has better convergency properties. Let us try and define a basis set by the equation: 1 2 ~ p~ 20 ~ F ðxÞ Z 0 K V C bn V0 ðxÞ C (A19) 2 2 n The orbital exponents be defined according to Eqs. (A7) and (A8). In principle p~0 can take any value. If p~ 0 Z p0 ; the basis set F~ n ðxÞ coincides with Fn(x). Expanding the wavefunction by the basis set defined in Eq. (A19) and eliminating the kinetic energy operator as in Eq. (A12), we can transform the SE into: X p~ 20 K K b~ n V0 ðxÞ C V0 ðxÞ C V 0 ðxÞ K E Cn F~ n ðxÞ Z 0 2 n ðA20Þ Multiplying on the left by F~ n0 ðxÞ; integrating over the electronic coordinates, dividing by p~0 and considering Eq. (A10), we get: X S 0 p~ Tn;n0 K p~ 0 dn;n0 C 0 Sn;n0 C n;n E Cn Z 0 (A21) 2 p~ 0 n
The symbol Sn,n 0 stands for the overlap between basis functions: Ð Sn;n0 h dxF~ n0 ðxÞF~ n ðxÞ (A22) Eq. (A21) is a generalized eigenvalue problem, where the eigenvalue is the energy E. We have to investigate which is the best value for p~0 ; and whether it coincides with p0. In the particular case of a single basis function F~ n ðxÞ (a trivial 1!1 matrix) our problem becomes: Tn;n K
E p~0 C Z0 p~0 2
(A23)
So we get the energy as a function of the parameter p~0 : EZ
p~ 20 K p~ 0 Tn;n 2
(A24)
2 =2; when The energy assumes the minimum value, KTn;n p~ 0 Z Tn;n : Exploiting Eq. (A17) would give the same value for the energy, so we conclude that in the particular case of a basis set made of one function, the best value for p~ 0 is p0. For the more interesting case of a larger basis set, we are not able to give any general proof, but considering a specific case we demonstrated that the best value for p~ 0 is not always p0. Let us approximate the ground state of Helium atom by a linear combination of two functions: F~ 1 ðxÞZ 1s2 ð1 SÞ and F~ 2 ðxÞZ 1s2sð1 SÞ: According to Eq. (A21), the approximated ground state energy is the lowest solution of the second order equation: p~ E T1;1 K 0 C Symm: p~ 0 2 (A25) S1;2 E Z 0 p~ 0 p~ 0 E T2;2 K C T1;2 C S1;2 C p~ 0 2 2 p~0
23
If we assign p~ 0 a value of 2.47 a.u., Eq. (A25) gives EZK2.8557 a.u., while exploiting Eq. (A17) one obtains p0Z2.389 a.u. and thus EZK2.8533 atomic units. This means that in general p0 is not the best choice for the definition of the energy of the basis set.
References [1] V. Aquilanti, S. Cavalli, C. Coletti, Phys. Rev. Lett. 1998; 809. [2] V. Aquilanti, S. Cavalli, C. Coletti, G. Grossi, Chem. Phys. 209 (1996) 405. [3] V. Aquilanti, A. Caligiana, S. Cavalli, Int. J. Quantum Chem. 92 (2003) 99. [4] V. Aquilanti, A. Caligiana, S. Cavalli, C. Coletti, Int. J. Quantum Chem. 92 (2003) 212. [5] O. Goscinski, Adv. Quantum Chem. 41 (2002) 57. This article contains historical remarks and a reprint of the Preliminary Research Report No. 217, Quantum Chemistry group, Uppsala University, 1968. [6] J. Avery, J. Mol. Struct. (Theochem) 458 (1999) 1. [7] J. Avery, Adv. Quantum Chem. 31 (1999) 201. [8] J.P. Gazeau, A. Maquet, J. Chem. Phys. 73 (1980) 5147. [9] J. Avery, J. Math. Chem. 21 (1997) 285. [10] V. Aquilanti, J. Avery, Chem. Phys. Lett. 267 (1997) 1. [11] V. Aquilanti, J. Avery, Adv. Quantum Chem. 39 (2001) 71. [12] J. Avery, D.R. Herschbach, Int. J. Quantum Chem. 41 (1992) 673. [13] V. Aquilanti, S. Cavalli, C. Coletti, Chem. Phys. 214 (1997) 1. [14] V. Aquilanti, S. Cavalli, C. Coletti, D. Di Domenico, G. Grossi, Int. Rev. Phys. Chem. 20 (2001) 673. [15] V. Aquilanti, A. Caligiana, Chem. Phys. Lett. 366 (2002) 157. [16] J. Avery, Hyperspherical Harmonics and Generalized Sturmians, Kluwer, Dordrecht, 2000. [17] J. Avery, S. Sauer, Quantum Systems in Chemistry and Physics, vol. 1, Kluwer, Dordrecht, 1998. [18] P.M. Morse, H. Feshbach, Methods of Theoretical Physics, McGrawHill, New York, 1953. [19] T. Shibuya, C.E. Wulfman, Proc. R. Soc. A286 (1965) 376. [20] J. Avery, Hyperspherical Harmonics, Applications in Quantum Theory, Kluwer, Dordrecht, 1989. [21] B.R. Judd, Angular Momentum Theory for Diatomic Molecules, Academic Press, New York, 1975. [22] M. Abramowitz, I.A. Stegun, Handbook of Mathematical Function, Dover, New York, 1964. [23] J.F. Rico, R. Lo´pez, G. Ramı´rez, J. Comput. Chem. 9 (1988) 790. [24] J.F. Rico, R. Lo´pez, M. Paniagua, G. Ramı´rez, Comput. Phys. Commun. 64 (1991) 329. [25] W. Kolos, L. Wolniewicz, J. Chem. Phys. 49 (1968) 404. [26] H. Shull, P.O. Lo¨wdin, J. Chem. Phys. 30 (1959) 617. [27] R. Szmytkowski, J. Phys. B 30 (1997) 825. [28] R. Szmytkowski, J. Phys. A 31 (1998) 4963. [29] R. Szmytkowski, Phys. Rev. A 65 (2001) 012503. [30] R. Szmytkowski, J. Phys. B 35 (2002) 1379. [31] J.F. Rico, R. Lo´pez, I. Ema, G. Ramı´rez, Comput. Phys. Commun. 105 (1997) 216. [32] J.F. Rico, R. Lo´pez, G. Ramı´rez, C. Tablero, Phys. Rev. A 49 (1994) 3381. [33] V. Aquilanti, A. Caligiana, in: E.J. Bra¨ndas, E.S. Kryachko (Eds.), Sturmian Orbitals in Quantum Chemistry: An Introduction Fundamental World of Quantum Chemistry, vol. I, Kluwer, Dordrecht, 2003, pp. 297–316. [34] J. Avery, V. Aquilanti, A. Caligiana, Adv. Quantum Chem. 47 (2004) in press. [35] C. Guidotti, O. Salvetti, N. Durante, V.T. Lamanna, G.P. Arrighini, Int. J. Quantum Chem. 93 (2003) 59.