Ellinger Y., Defranceschi M. (eds.) Strategies and applications in quantum chemistry (Kluwer, 200

146	M. DEFRANCESCHI ET AL.

Nuclear attraction, electron-electron repulsion and exchange terms. Using the Eqs. 6 and

7, these contributions are respectively written in terms of the quantity Wij(q) previously

defined in Eq. 11 :

and

The detailed expression for is then :

The momentum space equivalent of the total energy E, is :

3. Principles for numerical resolutions

Because of the terms

and

explicit solutions to Eq. 3 cannot be obtained in

position space. In such cases approximate solutions are usually expressed as truncated

linear combinations of basis functions (LCAO expressions). In spite of its successes, the

LCAO approximation experiences various difficulties (truncation limits, nature of the

basis functions, etc.) hard to estimate and which are not entirely controllable [51].

QUANTUM CHEMISTRY COMPUTATIONS IN MOMENTUM SPACE

147

Furthermore their incidence are very dependent upon the nature of the properties [52,53]. Due to computer limitations, basis sets cannot be extended indefinitely and direct numerical evaluations seem the ultimate solution for molecules [54]. In position space this is a viable alternative for diatomic molecules [55,56], but it cannot be extended easily

to polyatomic systems. Formulated in momentum space, the HF equations have not

explicit solutions and the difficulties to express them in terms of basis functions are analogous to those encountered in the r-space. However the momentum space HF

equations give way to numerical approaches in which Coulombic interactions become

tractable even for polyatomic molecules [7] ; among other advantages, these equations,

Eqs. 13 and 20, do not require coordinate systems adapted to the geometry of the molecules to remove Coulombic singularities. In both equations the only singular

contribution comes from the factor.

3 . 1 . VARIATION-ITERATION PROCEDURE

In both position and momentum spaces, iterative procedures are necessary to solve the HF equations. Starting from a trial orbital an approximate orbital, is obtained after k+1 iterations from Eq. 13 rewritten as:

The procedure is repeated until convergence is reached. Since we are interested in bound

states where no problem of divergence or cusps conditions is raised. But the

method can be adapted to more general situations by introducing a translation of the energy origin in Eq. 21.

Numerical and computational problems associated with the implementation of the approach for routine use fall in two main categories : (a) numerical integration and (b)

enforcement of the orthogonality and renormalization of the numerical orbitals during the

iteration steps. Many different integration schemes have been considered in the past,

some of which will be detailed in the section 3.2. As concerns orthonormalization, at

148	M. DEFRANCESCHI ET AL.

each step the new iterates

even if initially orthonormal, need to be

renormalized and orthogonalized to form true canonical HF orbitals. Great care must be exercised in selecting orthogonalization procedures, for instance the so-called Löwdin's symmetric orthogonalisation procedure [57], often used in Quantum Chemistry, mixes all

the orbitals simultaneously, tends to contaminate all the iterates, and impairs the

convergence of the iterative steps. Schmidt orthogonalization does better (since it allows to choose the sequence of orthogonalization) but looses track of the symmetry of these

orbitals. Finally, the canonical orthogonalization performs a maximum in mixing of all states. We have shown [31] for the ground state of Be and B+ that the Gram-Schmidt procedure turns out to be more appropriate in most cases, but with very good trial functions, Löwdin's symmetric procedure yields equivalent results. In all cases reported in Table 1, Gram-Schmidt orthonormalization has been used.

3.2.NUMERICAL TECHNIQUES

Different integration schemes have been considered. To cancel the singularity factor in

Eq. 13 by the integration volume element, Navaza and Tsoucaris have proposed the use of spherical polar coordinates. However, because of the convolution integrals, interpolation schemes are needed in these coordinates since arguments (p-q) do not necessarily belong

to the grid points. The computation time increases as the square, , of the number of

points of the integration grid, and for large systems, this time becomes prohibitive.

Another point of view has been to focus on these convolution integrals and treat them via

a more economical fast Fourier transform procedure. In this case, the computation time

increases only as , but at the expense of an approximate treatment of the

singular factor [58,59]. Variants [60,61] based on the Fock transformation have also been proposed to deal with the infinite limits of integration resorting to a one-to-one

correspondence between intervals and . At the present time,

none of the approaches has been satisfactory enough to bring the fully numerical momentum quantum chemistry calculations beyond a stage of prematurity. Furthermore, computational tests [25] on helium atom have shown the importance of accuracy and convergence of the integrals. It seems that straightforward numerical calculations are not readily applicable and our work is now directed toward mixed numerical and analytical procedures.

QUANTUM CHEMISTRY COMPUTATIONS IN MOMENTUM SPACE

149

3.3.SEMI-ANALYTICAL TECHNIQUE

The method presented here allows, starting with trial gaussian functions, a partial analytical treatment which we have used to improve the LCAO-GTO orbitals (trial functions) essentially obtained from all ab initio quantum chemistry programs. As in r-

representation,	trial	functions	(Eq. 21) are conveniently expressed as linear
combinations	of	functions	themselves written as linear combinations of

gaussian functions (LCAO-GTO approximation)

and

where	is a normalized gaussian function expressed in momentum space. As
	belong to the Sobolev space		direct Fourier transformation leads to a
set	that fulfills the criterion about the convergence of the energy and wave
function (the completeness of the orbital			bases		is not sufficient
to guarantee the convergence of the energy and wave function in the norm of						; to
ensure this convergence the set			must be complete in		. The expression
for the first iterate		based on trial functions		expressed	as LCAO-GTO
expansions is thus :

The various quantities entering Eq. 24 are deduced when the trial orbitals

are

expressed as linear combinations of Gaussian functions, they are expressible in terms of

150	M. DEFRANCESCHI ET AL.

known transcendental functions. In the case of s-functions, two basic integrals,

respectively denoted by and

have to be solved to obtain

The details of the calculations can be found in Ref. 35, the final expressions are :

and Daw	dx is the so-called Dawson function. The individual terms
appearing in Eq. 24 are :
Nuclear attraction term:

QUANTUM CHEMISTRY COMPUTATIONS IN MOMENTUM SPACE

151

Fourier transform of orbital products :

Electronic repulsion term :

Exchange term :

So the first iteration transforms the trial wave functions expressed as linear combinations

of gaussian functions into an expression which involves Dawson functions [62,63]. We

have not been able to find a tabular entry to perform explicitly the normalization of the first

iterate, accordingly this is carried out numerically by the Gauss-Legendre method [64].

One of the drawbacks of the first iteration, however, is that computation of energy

quantities, e.g. orbital and total energies, requires to evaluate the integrals occurring in

Eq. 3 on the basis of the		Unfortunately, the transcendental functions in terms of
which the	are expressed at the end of the first iteration do not lead to closed form

expressions for these integrals and a numerical procedure is therefore needed. This constitutes a barrier to carry out further iterations to improve the orbitals by approaching the HF limit. A compromise has been proposed between a fully numerical scheme and the

simple first iteration approach based on the fact that at the end of each iteration the

entail the main qualitative characteristics of the exact solution and most

152 M. DEFRANCESCHI ET AL.

importantly the right asymptotic decay. The idea is thus to fit the iterated analytical

functions	obtained at the	step on a finite set of gaussian functions and then use
these fitted functions as a new set of trial functions			The advantage is twofold.

First, with exponents and linear coefficients specific for each orbital, energies and functions are quickly improved. Second, the problematic convolution products and integrals are efficiently computed in terms of the gaussian functions obtained to represent

the	The analytical functions		are represented as linear combinations of
gaussian functions,		This fit	is carried out using a modified version of the

Gausfit package [65] developed by Stewart [66] for gaussian fits of Slater functions. The resulting functions are analytically orthonormalized.

For atoms, the radial part of is expressed as a linear combination of spherical

gaussians, which, in the case of 2p orbitals writes as :

Given a radial function

to fit, one minimizes the variance,

where

is a function which weights the contributions to the integral according their

expected importance [28]. From several tests on Be and Ne we have found that the following weight functions are quite efficient:

Gaussian functions do not have the right asymptotic decay due to too low amplitudes in regions of large p values, therefore representations in terms of gaussians are of much

slower convergence than Slater functions. Since contributions from high momenta are essential to the energy, a second degree polynomial, Eq. 35, is used to enforce them in the

valence orbitals.

A set of nine gaussians allows a satisfactory fit with low variance, Eq. 34, values : about

for the 1s and 2p orbitals and	for the 2s orbital. The valence orbitals having
node(s) are slightly more difficult	to fit. Under these conditions, the iterative scheme

QUANTUM CHEMISTRY COMPUTATIONS IN MOMENTUM SPACE

153

converges to results close to HF limit, but obviously it cannot approach it completely because the fit is based on a limited number of gaussian functions.

4. Results and advantages

In Tab. 1 are given the various results obtained in our group ; the precision of the method

used for the resolution as well as the main interest of these results are summarized. When

results have been obtained numerically (section 3.1) the method is denoted num-SCF or num-MCSCF according to the level of theory used ; when an analytical treatment (section 3.2) has been performed, the denotation is analyt-gauss if the trial functions were expressed as linear combinations of gaussian functions or analyt-Slater if the trial functions were expressed as linear combinations of Slater functions. Finally when a semi-analytical treatment (section 3.3) has been done the method is called analyt-gauss*.

Results fall into three categories : the first one corresponds to pure numerical results on which have demonstrated the feasibility of numerical calculations. They have also provided momentum wavefunctions for physical quantities

such as Compton profiles [17], (e,2e) cross-sections [26]. In the second category we have investigated the possibilities of using a variation-iteration procedure defined in

momentum space to improve the one-electron states for various chemical systems

expressed as linear combinations of gaussian functions. Significant improvements in

energy quantities and properties sensitive to the shape of the wave function (Compton

profile, momentum distribution, etc.) were indeed noted. In particular, the first iteration transforms the trial wave function expressed as linear combinations of gaussian functions in an expression which involves Dawson functions. An asymptotic analysis carried on

the first iterate discloses a behavior quite close to the exact one. In the third category, the semi-numerical approach is used to provide physical quantities. Similarly to the position

space approach it is based on the variation principle which guides the changes of the

wavefunction	: the closer the energy E to	the nearer the trial wave function	the
ground state	In LCAO-SCF-MO schemes however, the function obtained		by

minimizing the total energy does not necessarily give a good description of properties such as multipole moments, while in momentum space due to the capacity of the method to improve the quality of a wavefunction significant improvements have been obtained e.g. for the dipole moment of the hydrogen fluoride [38].

154	M. DEFRANCESCHI ET AL.

QUANTUM CHEMISTRY COMPUTATIONS IN MOMENTUM SPACE

155