Towards a possible ab initio molecular mechanics.
I. Stability of density matrix elements¹

May 3, 2007

A.L. Tchougréeff,

A.M. Tokmachev

I. Mayer

Abstract

Our previous studies based on the semiempirical quantum chemical Hamiltonians of MINDO and MNDO types permitted us to draw the conclusion that the success of MM can indeed be understood on the basis of quantum mechanics. The analysis of the geminal-type wave functions constructed by making use of oriented hybrid orbitals showed that the parameters of the wave function and the energy contributions of the individual bonds are indeed well transferable, in good agreement with the simple chemical picture of the systems studied.

It seems to be desirable to develop a similar analysis also at the ab initio level of the theory, because that would connect together the chemical and physical description of molecules and explain the observed transferability of molecular interactions more rigorously than semiempirical theories can do it. On the other hand this treatment might be also useful for developing QM/MM junctions in the cases when the QM part of the systems are described at the ab initio level. We report here results of the first step of this analysis: the geminal parameters expressed in a symmetrically orthogonalized optimized minimal basis set exhibit a degree of stability similar to that observed in the semiempirical case.

Introduction

Chemistry is dominated by the idea of transferability either explicitly formulated or implicitly assumed. The concept of isolated chemical bonds - which are to a great extent independent of their surrounding - allowed the chemists to formulate a consistent approach which was remarkably successful already in the XIX century in interpreting and predicting different structural and reactivity features of various compounds, especially organic ones. All these achievements were reached not only without any relation to the quantum mechanical basics of the molecular electronic structure but even without the knowledge of electrons themselves, not to speak about the detailed understanding of the laws governing their motion. The modern incarnation of that more than a century old approach is represented by the various molecular mechanics" (MM) methods providing direct parameterization schemes for molecular potential energy as function of the coordinates of nuclei of the atoms composing the molecule Dashevskii,BurkertAllinger1982. Although the standard paradigm in the MM area is that there is no direct relationship between the quantum theory of molecular structure and the MM force fields, the entire MM theory is tacitly based upon the concept of the isolated two-center chemical bonds, which is of course a crude representation of a basically quantum mechanical nature of the ground state many-electron wave function of the molecule (see below).

Meanwhile, modern quantum chemistry gives very little for understanding of the fundamental transferability properties mentioned above. At best, they are reproduced in the calculations (numerical experiments) perfomed at different levels of sophistication, so that these important features remain experimental facts still missing whatever theoretical explanation. However, as stressed by Coulson, any explanation why must be given in terms of concepts which are regarded as adequate or suitable. So the explanation must not be that the electronic computer shows that D(H-F) >> D(H-H), since this is not an explanation at all, but merely a confirmation of experiment" []. Also a theoretical understanding of the origin of transferability observed in experiments of both real and numerical types could be of a siginificant practical use in the context of the hybrid QM/MM methods, since it would potentially provide the form of the junction between the parts of a complex molecular system treated on the ab initio QM and the classical MM levels, respectively.

In the recent series of papers [,,] an attempt has been undertaken to reconcile the quantum mechanical description with the idea of transferable bonds. That was done in the framework of a semiempirical quantum chemical model of molecular structure, in which one uses implicitly orthogonalized valence AOs and the Hamitonian is parametrized against experimental data. Technical details of these papers will be referred to throughout the present paper. But here we just want to give some very general notes concerning the overall structure of the theory which may provide an explanation for the bonds' transferability and represent a basis for a further development of a classical picture of the molecular structure, compatible with its quantum mechanical description.

First of all we notice that the transferability of bonds must be a feature pertinent to the wave function, which is built up as an (antisymmetrized) product of well-distinguishable two-electron components. This can be the only basis of the semi-observability of the chemical bonds, as defined by Ruedenberg [], which is tacitly assumed within the entire MM picture and more generally - in chemistry itself, where applies. It is important to mention here that for the molecular systems whose electronic structure does not fit into the picture of two-center two-electron bonds, the models of the MM force field type cannot be easily developed either. It is noteworthy that the molecular Hamiltonian is a rather nontransferable entity: those corresponding to different molecules differ by singular terms of electron-nuclear attraction in the direct physical space; in this respect they are not close in any reasonable mathematical sense. Similar consideration can be applied to the matrix form of the electronic Hamiltonian: those for different molecules act in spaces of different dimensionality and their matrix elements are geometry dependent. Thus there are no a priori signs of the transferabilty of the Hamiltonian as far as the description in the direct physical space is concerned. (By using a refined analysis, some transferable elements of the Hamiltonian could be identified in the framework of a special mixed" second quantization formalism [].) By contrast, formal characteristics of the transferability, which are pertinent to the wave functions, have been established in the framework of the semiempirical theory (see Refs. [,,] and below).

The minimal requirement to a subsystem within a molecular system to be distinguishable or observable is that it has a stable composition. In this respect we may ask that to what extent the numbers of particles residing in different parts of molecular system are well defined, and are, therefore, good quantum numbers. To find this out we have to consider the dispersion (or mean-square deviation) of the numbers of particles in the parts.

Let us introduce operators [^N]_M of the numbers of electrons in some parts M of the molecular system such that:

å M  ^ N   M  = ^ N   ,

(0)

where [^N] is the operator of the number of electrons in the entire molecular system under consideration. For the normalized wave function Y we get:

áY| ^ N   M  |Yñ = nM  ,     "M .

(0)

As we mentioned above, the physical distinguishability or observability of a part of the system is formally expressed as vanishing fluctuation (dispersion) of the number of electrons in it. The dispersion of the number of electrons in the M-th part:

áY|( ^ N   M  -nM) 2|Yñ

(0)

must be as small as possible.

The minimum condition for the overall dispersion

min å M  áY|( ^ N   M  -nM)2|Yñ

(0)

can be trivially satisfied by the function Y taken as the antisymmetrized product of group functions (GF):

Y = ( Õ M  gM+)|0ñ,

(0)

where each of the group functions g_M⁺ is an eigenfunction of the corresponding electron number operator [^N]_M:

^ N   M  gM+| 0 ñ = nMgM+|0 ñ .

(0)

This implies that n_M is an integer number. The above condition can be satisfied if one assumes that the entire space L of one-electron states of the molecular system is divided into disjunct subspaces L_M:

L = Å M  LM; LMÇLN = {0}, M ¹ N.

(0)

which are also assumed mutually orthogonal. Each of the subspaces L_M serves as a carrier space for the M-th electron group. To build up the operators g_M⁺, an orthogonal basis set { c_m} of vectors (functions) is introduced in each of the subspaces L_M. If we define the operator m⁺ as one creating an electron in the one-electron state c_m, then the operators g_M⁺ are constructed as expansions:

gM+ = å { mi}   C{ mi} NM Õ \substack i = 1 mi Î { M}   mi+

(0)

with the amplitudes C_{{ mi}} to be determined variationally. (Each {m_i} represents a selection of n_M elements of the local basis {c_m}.) Within this construct the hitherto vaguely defined parts of the system" become the carrier subspaces. The number of particles operator for each group or carrier subspace can be defined as:

^ N   M  = å m Î { M}   m+m-,

(0)

the average occupancies of the parts n_M become the integer numbers of electrons residing in each of the carrier subspaces, and the dispersions Eq. (3) vanish as each group fuction g_M⁺ is an eigenfunction of the respective operator [^N]_M with the eigenvalue n_M.

While the overall dispersion Eq. (4) vanishes for a wave function represented by a single GF product Eq. (5). This is not the case for the exact molecular wave function, which may be always represented as a linear combination of several functions with a GF structure, but with an arbitrary distribution of electrons between the groups. This physically corresponds to accounting the charge transfer between the subsystems. If the selection of the groups (break up of the entire orbital space into carrier subspaces and assigning integer numbers of electrons to them) was appropriate, then the dispersion Eq. (4) for the exact wave function will be small. This concept can be viewed alternatively as a separation of the Hamiltonian [^H] into its ßignificant" part [^H]₀:

^ H   = ^ H   0  + ^ W   ,

(0)

which commutes with all the number of particles operators introduced above:

é ë ^ N   M  , ^ H   0  ù û = 0,     "M

(0)

and its noncommuting part [^W] playing the rôle of an ïnteraction" or "perturbation" defined with respect to [^H]₀. This perturbation [^W] must be weak. Provided the set of the operators [^N]_M is optimized with respect to this criterion, the entire molecular system naturally falls out into subsystems. Under this condition the magnitude of the dispersions Eq. (4) for the exact wave function is controlled by the same smallness parameters as the convergence of the perturbation series with the perturbation [^W]. The matrix elements of [^W] are non-vanishing only betwen the approximate ground state Y (or other eigenvectors of [^H]₀) and those other states Y^¢ which differ from Y by the the numbers of electrons in some subsystems. Then the distinguishability means that

áY| ^ W   |Y¢ñ DECT << 1,

(0)

where DE_CT stands for the characteristic energy of a charge transfer state.

Our following considerations are based on the trial wave function in the GF form. This form of the wave function assures a simplified structure of the one- and two-electron density matrices and, as a consequence, a simplified expression for the molecular energy, which can be written as a sum of intragroup contributions and Coulomb and exchange interactions between electron groups. Selecting an appropriate form of the electronic wave function for each of the groups allows to formalize the chemical observations on the transferability of certain parts of molecules and to establish transferability of certain types of electronic structure variables (ESVs) as well. The simplest form of the wave function relevant to our purpose of reconstructing a classical bonding picture in organic compounds is based on the antisymmetrized product of two-electron functions (geminals) [,]. This is particularly relevant in the case of usual saturated organic molecules, the electronic structure of which can be described by two-center, two-electron bonds and lone electron pairs - basically following Lewis. In the semiempirical context each carrier space could be specified by a pair of strictly local hybrid orbitals (HOs) centered on the atoms connected by the bond at hand. This form of the trial wave function in combination with the semiempirical Hamiltonians allowed us to derive a mechanistic picture for the molecular energy Tch62,Tch64,Tch73 and to show a path to the junction between QM and MM parts in the hybrid QM/MM methods. It has been demonstrated [] that the matrix elements of one- and two-electron density matrices are stable - and thus transferable for analogous bonds - up to certain degree in two small parameters representing intrabond electron correlation and polarity. Remarkably, the off-diagonal matrix element of the one-electron density matrix - the Coulson bond order - is stable (transferable) within the second order in both these small parameters. It is noteworthy that the transferability of the density matrix elements takes place namely in the basis of the strictly local HOs. Although it is possible that the transferability of the density matrix elements would take place in some other basis as well, this is definitely not the case for an arbitrary basis (e.g., for the original AO basis).

The HOs used for constructing the geminals turn out to be very sensitive to the variations of the molecular Hamiltonian irrespective to the origin of these variations: either the geometry changes or the changes of the chemical composition of the molecular system. The actual form and orientation of the HOs spanning the carrier space for the respective bond geminal takes upon itself a major part of the variations in the Hamiltonian, whereas the density matrix elements in the basis of HOs remain relatively stable as noted above.

In the present paper we start the long way towards an ab initio MM theory. The semiempirical approach to deductive MM serves us as an example and we try to extend it taking into account the complications arising in the ab initio domain. Our first objective is to study the stability (transferability) properties of the density matrix elements. The rest of the paper addresses this issue: we devise a quantum chemical theory underlying the tentative ab initio MM and use it to establish the transferability of the density matrix elements on the solid theoretical basis supported by numerical experiments.

Strictly local geminals formalism for the electronic wave function

The general and formal construction of the classical picture of chemical bonding starts from selecting an adequate form of the trial wave function to be used in the derivations and numerical experiments. The deductive MM method based on the semiempirical Hamiltonians employs a geminal-based wave function. An important restriction on the structure of the bonding geminals is that each of them is constructed in a carrier space spanned by two strictly local HOs | r_m ñ and | l_m ñ (four spin-orbitals):

gm+ = umrma+rmb++vmlma+lmb++wm( rma+lmb++lma+rmb+) ,

(0)

where the amplitudes of configurations u_m, v_m, and w_m are subject to the normalization condition:

um2+vm2+2wm2 = 1

(0)

and should be determined variationally. In the case of a lone electron pair the corresponding carrier space is spanned by one HO and the geminal is trivial:

gm+ = rma+rmb+.

(0)

This type of geminals based on the HOs can be considered as strictly local geminals (SLG) [,]. The overall wave function of electrons in the molecule is then taken as the antisymmetrized product of the geminals, in accordance with Eq. (5).

The operators r_ms⁺ and l_ms⁺ create an electron with the spin projection s in one-electron states centered on the "right" and on the "left" ends of the bond, respectively. These operators are special cases of the general operators m_i⁺ of Eq. (8) for the case of two-electron two-orbital groups (geminals). The structure of these operators should reflect the strictly local character of the geminals. These orbitals are obtained by transformations involving only orbitals of the same atom:

tms+ = å a Î A  htma(A)as+;    t = r,l.

(0)

In this formula a_s⁺ stand for the operators creating electrons in the (orthogonal) basis AOs centered at the atom A; the matrices h(A) have the SO(4) structure and transform valence s- and p-orbitals into a set of four orthogonal HOs. The independent parameters of these matrices produce another subset of the ESVs, which describe the actual hybridization related ones and should also be determined variationally.

This picture is natural for semiempirical settings, but it cannot be directly transferred to the ab initio context where a series of complications arises. First of all, the dimensionality of the atomic subspace is larger than four due to explicit presence of the core shells for the heavy atoms (second row and higher), even if a minimal basi set is used. Another complication comes from the fact that the standard ab initio methods employ non-orthogonal AOs basis sets and one cannot form orthogonal HOs simply by performing intraatomic orthogonal transformations of them. Nevertheless, we can propose a procedure for obtaining approprate one-electron states, which can be naturally assigned to geminals. We perform this procedure in two steps. The first step is an orthogonalization of the initial set of AOs (including orbitals assigned to different atoms). Obviously there is some arbitrariness in this step. We employ the symmetric (Löwdin) orthogonalization to obtain the orthonormal set of orbitals, which is as close to the initial one as only possible. (It also conserves the assignment of the orbitals to atoms.) The second step is an intraatomic unitary transformation of orbitals. It is performed for each atom separately:

tms+ = å a Î A  Htma(A)as+,

(0)

where the operators a_s⁺ correspond to the orthonormalized set of AOs and the transformations H(A) is not limited to four orbitals. The whole transformation is a product of the (atomic-)block-diagonal transformation matrix with blocks given by the matrices H(A) and the symmetrical orthonormalization matrix S^-1/2. It is noteworthy that the order of the two steps can be interchanged, i.e. the whole procedure can be viewed as a formation of non-orthogonal HOs by the transformations H(A) with subsequent orthogonalization of all HOs.

The assignment of HOs obtained by the transformation Eq. (17) to two-electron groups is not yet specified. It is possible to distribute all the orbitals between geminals or to use one orbital as a core and a subset of four orbitals as valence HOs assigned to geminals. At the same time, our ultimate goal of deriving a classical (MM) picture makes preferable to recover hybridization structure Eq. (16) from the general transformation H(A) because it allows to construct the most important valence bond configurations and to relate the structure of one-electron states to the stereochemistry of an atom. We impose a set of restrictions on the structure of the matrix H(A) and consider it as a product of transformations. First we form new shells (sets of basis AOs; an s-shell is given by one function while a p-shell is given by three functions with the same radial part). The optimal sets of s- and p-functions are obtained by two independent orthogonal transformations of the basis s- and p-shells, respectively. Among these shells one s-shell (AO) corresponds to the core, while the remaining s- and p-shells can be combined to produce some sp-shells (basis subsets containing one s-AO and three p-AOs with correct transformation properties). In our subsequent calculations we use only one of these sp-shells. At the second step we apply some hybridization SO(4) transformation mixing one s- and three p-orbitals for each sp-shell. Therefore, we produce the same setting for studying the transferability of the density matrix elements as we had in the semiempirical case. The structure of the HOs and the algorithms for determination of the hybridization-related ESVs are important for the construction of ab initio MM and these questions will be addressed elsewhere. Here we concentrate on the density matrix elements in the basis of these HOs.

The use of only the minimal number of basis orbitals at this stage of development may be justified by the known fact that even the use of STO-NG type minimal basis sets at the RHF level of theory represents an often surprisingly good ``model chemistry". Furthermore, the a posteriori analyses of large scale calculations usually also lead to as many significantly populated effective core and valence orbitals as is the number of functions in the minimal basis sets [,].

Biorthogonal treatment for the geminals with overlap

General biorthogonal quasi-energy formula for geminals

Overlap of the orbitals assigned to the different bonds makes it difficult to calculate expectation values - in particular the total energy - already at the single determinant level of theory, However, the situation there is much simpler because any single determinant wave function can be expressed (up to an insignificant normalization constant) by using orthonormalized one-electron orbitals, which is not the case for geminals.

In the single determinant case one usually calculates the energy after performing an orthogonalization of the orbitals. However, it is easy to show that the energy of a single determinant wave function |Y_SDñ can be given also as

E = áYSD| ^ H   |YSDñ áYSD|YSDñ = á ~ Y   SD  | ^ H   |YSDñ

(0)

where |[(Y)\tilde]_SDñ is the single determinant built up of the orbitals {[(j)\tilde]_i} representing the biorthogonal counterparts of the nonorthogonal orbitals {j_i} contained in the determinant |Y_SDñ:

~ j   i  = å j  (S-1)jijj

(0)

Here s is the overlap matrix of the occupied orbitals having the elements S_ij = áj_i|j_jñ. (When writing equation (18), we have taken into account that á[(Y)\tilde]_SD|Y_SDñ = 1, owing to the biorthogonality of the two sets of orbitals.)

No such strict relationship holds for the geminal wave functions. However, at least at the equilibrium geometries, each geminal is expected to be not too far from its respective Hartree-Fock approximation. Therefore, for qualitative considerations it may be of meaning to consider the biorthogonal quasi-energy of the geminal wave function |Y_Gñ:

~ E   = á ~ Y   G  | ^ H   |YGñ á ~ Y   G  |YGñ

(0)

We note that the equality [E\tilde] = E would hold also if |Y_Gñ were the exact (full CI) wave function.

The idea of using the biorthogonal quasi-energy for qualitative considerations of geminals, describing different chemical bonds is not new; for instance it has been proposed in the book of Surján []. However, there one cannot find all the necessary quantities, and the expressions are in a form not very suitable to be used in our case. Therefore we have decided to present here the result explicitly.

It is easy to see, that if one introduces the biorthogonal transformation of the basis orbitals |l_mñ,|r_mñ then the geminals, |[g\tilde]_mñ formed with the same coefficients u,v and w of the latter will be ßtrongly biorthogonal" to the original ones (provided that the orbitals within same geminal have been orthogonalized). As a conseqence, the equality á[(Y)\tilde]_G|Y_Gñ = 1 holds, one-electron matrix elements appear only between orbitals assigned to the same geminal and two-electron ones involve at most two geminals. A somewhat lengthy but essentially straightforward derivation gave the result:

~ E   = á ~ Y   G  | ^ H   |YGñ á ~ Y   G  |YGñ = á ~ Y   G  | ^ H   |YGñ\notag (0) = n å i = 1  ì í î 2(ui2+wi2)á ~ r   i  | ^ h   |riñ+2wi(ui+vi)(á ~ r   i  | ^ h   |liñ+á ~ l   i  | ^ h   |riñ) (0) +2(vi2+wi2))á ~ l   i  | ^ h   |liñ+ui2[ ~ r   i  ~ r   i  |riri]+2uiwi[ ~ r   i  ~ r   i  |liri]+uivi[ ~ r   i  ~ r   i  |lili] \notag (0) +2uiwi[ ~ r   i  ~ r   i  |riri]+2wi2[ ~ r   i ~ l   i  |liri]+2wi2[ ~ r   i  ~ l   i|rili]+2uiwi[ ~ r   i  ~ l   i  |lili] \notag (0) + +uivi[ ~ l   i  ~ l   i  |riri]+2viwi[ ~ l   i  ~ l   i  |liri]+vi2[ ~ l   i  ~ l   i|lili] ü ý þ +2 å i < j  ~ E   Cx ij  \notag (0)

Here [E\tilde]_ij^Cx describes the Coulomb and exchange interactions between the i-th and j-th geminals in the biorthogonal framework; it is given by the following expression

~ E   Cx ij  = (ui2+wi2)(uj2+wj2) æ è 2[ ~ r   i  ~ r   j  |rirj]-[ ~ r   i  ~ r   j|rjri] ö ø \notag (0) + wi(ui+vi)(uj2+wj2) æ è 2[ ~ r   i  ~ r   j|lirj]-[ ~ r   i  ~ r   j  |rjli] ö ø \notag (0) + (ui2+wi2)wj(uj+vj) æ è 2[ ~ r   i  ~ r   j|rilj]-[ ~ r   i  ~ r   j  |ljri] ö ø \notag (0) + wi(ui+vi)wj(uj+vj) æ è 2[ ~ r   i  ~ r   j|lilj]-[ ~ r   i  ~ r   j  |ljli] ö ø \notag (0) + (ui2+wi2)wj(uj+vj) æ è 2[ ~ r   i  ~ l   j|rirj]-[ ~ r   i  ~ l   j  |rjri] ö ø \notag (0) + wi(ui+vi)wj(uj+vj) æ è 2[ ~ r   i  ~ l   j|lirj]-[ ~ r   i  ~ l   j  |rjli] ö ø \notag (0) + (ui2+wi2)(uj2+wj2) æ è 2[ ~ r   i  ~ l   j|rilj]-[ ~ r   i  ~ l   j  |ljri] ö ø (0) + wi(ui+vi)(vj2+wj2) æ è 2[ ~ r   i  ~ l   j|lilj]-[ ~ r   i  ~ l   j  |ljli] ö ø \notag (0) + wi(ui+vi)(uj2+wj2) æ è 2[ ~ l   i  ~ r   j|rirj]-[ ~ l   i  ~ r   j  |rjri] ö ø \notag (0) + (vi2+wi2)(uj2+wj2) æ è 2[ ~ l   i  ~ r   j|lirj]-[ ~ l   i  ~ r   j  |rjli] ö ø \notag (0) + wi(ui+vi)wj(uj+vj) æ è 2[ ~ l   i  ~ r   j|rilj]-[ ~ l   i  ~ r   j  |ljri] ö ø \notag (0) + (vi2+wi2)wj(uj+vj) æ è 2[ ~ l   i  ~ r   j|lilj]-[ ~ l   i  ~ r   j  |ljli] ö ø \notag (0) + wi(ui+vi)wj(uj+vj) æ è 2[ ~ l   i  ~ l   j|rirj]-[ ~ l   i  ~ l   j  |rjri] ö ø \notag (0) + (vi2+wi2)wj(uj+vj) æ è 2[ ~ l   i  ~ l   j|lirj]-[ ~ l   i  ~ l   j  |rjli] ö ø \notag (0) + wi(ui+vi)(uj2+vj2) æ è 2[ ~ l   i  ~ l   j|rilj]-[ ~ l   i  ~ l   j  |ljri] ö ø \notag (0) + (vi2+wi2)(vj2+wj2) æ è 2[ ~ l   i  ~ l   j|lilj]-[ ~ l   i  ~ l   j  |ljli] ö ø \notag (0)

One may introduce the biorthogonal Coulomb and exchange-type operators by the definitions of their action on an arbitrary orbital j as

^ ~ J   ajbj  j(1) = ó õ ~ a   * j  (2) bj(2) r12 dv2j(1) \notag (0) (0) ^ ~ K   ajbj  j(1) = ó õ ~ a   * j  (2) j(2) r12 dv2 bj(1) \notag (0)

Using them, one may define the effective (biorthogonal) potential due to the electrons in the j-th geminal as

^ ~ U   eff j  = (uj2 + wj2) æ ç è 2 ^ ~ J   rjrj - ^ ~ K   rjrj  ö ÷ ø +(vj2 + wj2) æ ç è 2 ^ ~ J   ljlj  - ^ ~ K   ljlj  ö ÷ ø \notag (0) (0) + wj(uj+vj) æ ç è 2 ^ ~ J   rjlj  - ^ ~ K   rjlj  +2 ^ ~ J   ljrj  - ^ ~ K   ljrj  ö ÷ ø \notag (0)

and then the effective (biorthogonal) core for the i-th geminal will be

^ ~ h   eff i  = ^ h   + n å j = 1  j ¹ i  ^ ~ U   eff j  = - 1 2 D- å a  Za ra + n å j = 1  j ¹ i  ^ ~ U   eff j  (0)

Then the biorthogonal quasi-energy of the whole system becomes

~ E   = n å i = 1  ì í î 2(ui2+wi2)á ~ r   i  | ^ ~ h   eff i  |riñ+ 2wi(ui+vi)(á ~ r   i  | ^ ~ h   effi|liñ+á ~ l   i  | ^ ~ h   eff i  |riñ) \notag (0) + 2(vi2+wi2))á ~ l   i  | ^ ~ h   eff i  |liñ+ ui2[ ~ r   i  ~ r   i  |riri]+ 2uiwi[ ~ r   i  ~ r   i  |liri]+uivi[ ~ r   i  ~ r   i  |lili] \notag (0) + 2uiwi[ ~ r   i  ~ r   i  |riri]+2wi2[ ~ r   i  ~ l   i  |liri]+2wi2[ ~ r   i  ~ l   i  |rili]+ 2uiwi[ ~ r   i  ~ l   i  |lili] \notag (0) + +uivi[ ~ l   i  ~ l   i  |riri]+2viwi[ ~ l   i  ~ l   i  |liri] + vi2[ ~ l   i  ~ l   i  |lili] ü ý þ (0)

Considerations based on excluding effects of BSSE-type

In the present study we are interested in the factors influencing the geminal coefficients and the resulting density matrix elements, rather than in the actual value of the total energy. Assuming that the orbitals r_i, l_i belonging to the same geminal are orthogonalized and that they are optimized by using a sufficiently large local basis, then one may use the conventional (not biorthogonal) integrals when determining the geminal coefficients .(According to our experience gained in the field of intermolecular interactions [], the energy should preferably be calculated without such a substitution.) That is the case because the difference between an integral containing the biorthogonal orbital in the bra" and that containing there the original one can be connected with the finiteness of the basis - these effects cause the known basis set superposition error" (BSSE) in the calculation of intermolecular interactions. They are known to diminish as the basis set improves.

For that reason we are going here to show that the difference between integrals containing the biorthogonal and original orbitals is a basis finiteness effect. As the derivation is not particularly related to the geminal structure of wave functions, we shall use simplified notations m,n, [(n)\tilde] to denote the individual orbitals encountered. Also, the derivation will be pertinent to some ``units" A - which in our case can be a geminal.

Let us consider the difference of ïntra-unit" (m,n Î A) integrals

á ~ n   | ^ h   |mñ-án| ^ h   |mñ

(0)

and introduce the projector [^P]^A on the subspace of the basis orbitals assigned to the given unit. If the orbitals within each unit are orthonormalized, then

^ P   A   = å r Î A  |rñár|

(0)

and we may write

á ~ n   | ^ h   |mñ-án| ^ h   |mñ = \notag (0) = á ~ n   |( ^ P   A   +1- ^ P   A   ) ^ h   |mñ-án| ^ h   |mñ\notag (0) = á ~ n   | ^ P   A   ^ h   |mñ-án| ^ h   |mñ+á ~ n   |(1- ^ P   A   ) ^ h   |mñ (0) = á ~ n   | å r Î A  |rñár| ^ h   |mñ-án| ^ h   |mñ+á ~ n   |(1- ^ P   A   ) ^ h   |mñ\notag (0) = å r Î A  dnrár| ^ h   |mñ-án| ^ h   |mñ+á ~ n   |(1- ^ P   A   ) ^ h   |mñ\notag (0) = á ~ n   |(1- ^ P   A   ) ^ h   |mñ\notag (0)

which would vanish if the function [^h]|mñ could exactly be expanded with the aim of the basis orbitals assigned to group A - thus it is connected with the BSSE effects. Analogous considerations apply for the two-electron integrals, too.

Stability of the density matrix elements

In the present Section we prove that under very nonrestrictive conditions the matrix elements of the one- and two-electron density matrices satisfy the stability (transferability) conditions also for the ab initio effective bond Hamiltonians, which are defined by the chemical composition and geometry of the molecular system. We show that for a broad class of molecules the latter two characteristics do not affect to a significant degree the values of the density matrices in the basis of the symmetrically orthogonalized HOs.

Effective bond Hamiltonians

Within the original SLG approach [,] the geminals are characterized by the amplitudes u_m, v_m, and w_m (see Eq. (13)). We introduce a new notation z_m = Ö2w_m simplifying the normalization condition for the amplitudes to: u_m²+v_m²+z_m² = 1. The effective Hamiltonians for the geminals representing chemical bonds are expressed in terms of the molecular integrals in the HO basis. Each geminal is expanded over three singlet two-electron basis configurations and the optimal values of their amplitudes are given by the solutions of the eigenvector problems (see also ParksParr):

æ ç ç ç è Rm Dm Fm Dm Cm Gm Fm Gm Lm ö ÷ ÷ ÷ ø æ ç ç ç è um zm vm ö ÷ ÷ ÷ ø = em æ ç ç ç è um zm vm ö ÷ ÷ ÷ ø ,

(0)

corresponding to their lowest eigenvalues. The actual expressions for the matrix elements of the effective bond Hamiltonians defined in the basis of orthogonal (and real) HOs are (subscript m is omitted for brevity):

R = 2hrr+[rr|rr], L = 2hll+[ll|ll], C = (hrr+hll)+[rl|rl]+[ll|rr], F = [ll|rr], D = Ö2(hlr+[rl|rr]), G = Ö2(hlr+[rl|ll]),

(0)

where h_tt^¢ stands for the matrix element of the one-electron part of the effective Hamiltonian for the bond under consideration and [tt^¢|t^¢¢t^¢¢¢] are the two-electron Coulomb repulsion integrals written according to the [12|12] convention. The difference with the semiempirical case is that F ¹ 0 and D ¹ G, in general.

It is constructive to perform the transformation to the basis of normalized symmetric and antisymmetric combinations of the ionic configurations [1/(Ö2)](r_ma⁺r_mb⁺± l_ma⁺l_mb⁺). Then the effective Hamiltonian can be written as:

^ H   = 1 2 æ ç ç ç è R+2F+L Ö2( D+G) R-L Ö2( D+G) 2C Ö2( D-G) R-L Ö2( D-G) R-2F+L ö ÷ ÷ ÷ ø ,

(0)

where the first state in the list is the symmetric combination of the ionic configurations, the second state is the covalent configuration, and the last one is the antisymmetric combination of the ionic configurations. This form of the effective Hamiltonian allows to extract the most important composite parameters determining the structure of the bond. A key parameter is the difference of the energies of the symmetric ionic and covalent states:

Dg = 1 2 (R+L)+F-C = 1 2 ( [rr|rr]+[ll|ll]) -[rl|rl] > 0.

(0)

This parameter equals to the difference between the average of the intrahybrid Coulomb repulsions and their interhybrid counterpart. In the simple symmetric case it would correspond to the difference of the one-and two-center Coulomb integrals g₁₁ and g₁₂. Another important parameter is b:

b = - 1 2Ö2 (D+G) = -hlr- 1 2 ([rl|rr]+[rl|ll]) > 0,

describing electron hopping between the ends of the bond. The other two characteristics represent the bond's asymmetry:

R-L = 2(hrr-hll)+[rr|rr]-[ll|ll], Ö2(D-G) = 2([rl|rr]-[rl|ll]).

(0)

The sign of these matrix elements in the case of polar bonds depends on the arbitrary assignment of the "right" and "left" labels to atoms. If the "right" atom is more electronegative than the "left" one then R-L < 0 and, typically, D-G > 0. We use this convention in the following.

It is convenient to use some dimensionless parameters instead of the matrix elements. We use Dg as a measure of all the interactions and introduce a set of dimensionless parameters:

z = 4b Dg , (0) l = L-R Dg , \notag (0) c = Ö2(D-G) Dg , \notag (0) w = 4F Dg . \notag (0)

All the parameters introduced are non-negative and the last two parameters do not appear in the semiempirical case. They allow to re-write the effective Hamiltonian Eq. (63) in a more convenient form:

^ H   = C ^ I   +Dg æ ç ç ç è 1 -z/2 -l/2 -z/2 0 c/2 -l/2 c/2 1-w/2 ö ÷ ÷ ÷ ø .

(0)

It is remarkable that the matrix in the parentheses determines the structure of the geminal irrespective to the average value of energy. In the important case of symmetric bonds (or symmetrized effective Hamiltonians) all the variety of possible solutions (ground state wave functions) is parameterized by a single parameter z depending on the given chemical composition of the bond, the shape and orientation of the HOs and the bond length. Already this observation brings one to the conclusion that some degree of the transferability of geminal amplitudes is can be expected, because the multidimensional manifold of the matrix elements of the molecular electronic Hamiltonians boils down to a single parameter z, which by only this reason can be the same for very different molecules. Below we, however, shall show even more: namely, that the significant components of the electronic structure description in the SLG approximation - the elements of the one- and two-electron density matrices in the HOs basis are described by the formulae which further reduce the possible variance of the the quantities of interest. This result may look trivial because the structure of the symmetric geminal actually depends on only one parameter but the point is that the parameter z has a simple and explicit expression through the parameters of the problem.

Ground state of the effective symmetrized bond Hamiltonian and the stability of the density

The actual solution of the eigenvector problem for the effective symmetrized Hamiltonian (with the antisymmetric ionic state excluded) reduces to diagonalizing 2×2 matrix

æ ç è 1 -z/2 -z/2 0 ö ÷ ø

depending on the parameter z (see. Eq. (70)). We denote the amplitudes of the symmetric ionic and covalent configurations in the ground state as sinj and cosj, respectively. Diagonalizing the Hamiltonian Eq. (70) yields (in the symmetric case) the following solution for the ground state:

tan2j = z

(0)

or in a detailed form,

u  = v = sinj Ö2 = 1 2   æ  ú Ö 1- 1 G(z)   , \notag (0) (0) w = cosj Ö2 = 1 2   æ  ú Ö 1+ 1 G(z)   , \notag (0)

where the notation

G(z) = Ö   1+z2

(0)

is introduced for the hereinafter ubiquitous square root.

The amplitudes of configurations enter the expressions for any observables (including the energy) not by themselves, but only in some combinations: namely, the elements of the one- and two-electron density matrices written in the basis of the HOs spanning the carrier space for the bond geminal under consideration. Those for the m-th bond (geminal) are by definition:

Pmt¢t = å s  á 0| gmtms+tms¢gm+|0 ñ , Gmt¢¢¢t¢¢tt¢ = å s  á 0| gmtms+tm-s¢+tm-s¢¢tms¢¢¢gm+|0 ñ .

(0)

All the necessary matrix elements can be readily determined as simple functions of geminal amplitudes u_m, v_m and w_m taking into account the actual definition of the geminals Eq. (13). This results in the following:

Pmrr = 2( um2+wm2) , Pmll = 2(vm2+wm2) , Pmrl = Pmlr = 2(um+vm)wm, Gmrrrr = 2um2, Gmllll = 2vm2, Gmrrll = Gmlrrl = 2wm2, Gmlrlr = Gmrlrl = 2umvm,Gmrrlr = Gmrlrr = 2umwm,Gmrlll = Gmlrll = 2wmvm

(0)

where we took into account that some of the integrals appear twice in the expression for the energy for different spin projection indices on which the density matrix elements are themselves independent. Only the above quantities each of which is not more complicated than involving one square root irrationality (for the symmetric or symmetrized bond) enter the expression for the energy in terms of the Hamltonian matrix elements in the HOs basis. In the following we shall consider only three density matrix elements P_m^rl = 2w_m(u_m+v_m),   P_m^rr-P_m^ll = 2(u_m²-v_m²),    G_m^rrll = 2w_m²,which naturally reproduce three important characteristics: Coulson-type bond order (off-diagonal element of the first order density matrix), bond polarity and bond covalency.

The solution of the symmetric problem results in the following values of the density matrix elements:

Prl = z G(z) ,    Prr-Pll = 0,    Grrll = 1 2 æ ç è 1+ 1 G(z) ö ÷ ø .

(0)

As we mentioned above, the very possibility to characterize any symmetrical bond by a single parameter explains a lot of the observed transferability: the entire diversity of semiempirical parameters and/or ab initio basis sets boils down to a single number (Anderson's parameter z). The definition of z (Eq. (66)) shows that z®0 in the limit of infinite separation, corresponding to the homolitic bond cleavage. However, the analysis performed in a semiempirical setting Ref. [] indicated that the interatomic separations characteristic for the real bonds correspond to values z >> 1. It means that the density matrix elements have the asymptotic behavior described by:

Prl\asymp 1- 1 2z2 ,    Grrll\asymp 1 2 æ ç è 1+ 1 z ö ÷ ø ,

(0)

where z^-1 serves as a small parameter. The most important feature of the above expression is that the variations of the off-diagonal matrix element of the one-electron density matrix (Coulson-type bond order) with respect to the limiting value of unity are of the second order with respect to the small parameter z^-1. This explains the results of the numerical experiments of Refs. [,] performed on organic compounds of different classes (alkanes, alcohols, amines etc.) at the semiempirical level, where the Coulson bond orders P_m^rl all acquired values between 0.92 and 1.0. (Of course, the asymmetry of bonds had also been properly taken into account.)

The same behaviour is observed in the ab initio case, too. We are going to support the present theoretical conclusions by the numerical results performed in the framework of the fully variational SLG model described in the previous Section employing the 6-31G basis set. We repeated some calculations in larger basis sets 6-311G and 6-311++G and those calculations demonstrated that although the numerical values were affected, all the qualitative conclusions remain valid. Table presents the results of our calculations for a set of simple test molecules with characteristic chemical bonds in different environment. First, the Table presents the numerical values of the parameters determining the structure of the Hamiltonian Eq. (70). We can see that for all molecules at their equilibrium geometry structures the values of z are considerably larger than unity and that the second important is the asymmetry parameter l (see below). The Table also presents the estimates of the density matrix elements Eq. (78), which are compared with their actual values obtained in the ab initio SLG computations. As expected, the two schemes give exactly the same values for non-polar chemical bonds. The agreement is reasonable even for highly polar bonds. At the same time the difference between the two schemes quickly increases when the polarity of the bond becomes larger. It means that further corrections taking into account the bond asymmetry are necessary. But before we perform a more refined analysis of the density matrix elements, we consider another model based on the symmetric structure of the geminal.

Estimates of the density matrix elements from the Heitler-London solution

The non-orthogonality of the one-electron basis functions represents one of the most important differences between the ab initio and semiempirical settings compared in the present paper. In the semiempirical setting as a rule - though nowadays it is not mandatory - the AOs centered on different atoms are assumed orthogonal with all effects of non-orthogonality reloaded to the parameters of the method. On the contrary, in the ab initio setting the overlap of the different one-electron functions is treated explicitly. We have considered above the model based on the single parameter z_m, reflecting the electron correlation within a symmetrized geminal built up of Löwdin-orthogonalized orbitals. It led to establishing the transferability properties and the numerical results similar to those obtained in the semiempircal case []. However, in the ab initio setting the overlap integral S of the non-orthogonalized HOs ascribed to bond geminal can be considered as an alternative parameter, the effect of which on the result may be important. These quantities were not explicitly used in the above estimates since those were performed in the basis of orthogonal two-electron configurations constructed from the orthogonalized HOs. By virtue of this, the overlap integral affected the matrix elements of the effective bond Hamiltonian and the parameter z only indirectly. Nevertheless, we can determine the overlap between the HOs with the knowledge of the transformation matrices H(A) Eq. (17) and use it explicitly for constructing geminal structure models.

It is known that the Heitler-London (HL) wave function gives a good estimate for the bond-length and for the energy of the H₂ molecule. It has the form of the covalent (homopolar) contribution to the bond geminal, but it is built up from the original, i.e. non-orthogonal HOs. Transforming the HL function to the basis of orthogonal HOs results in the appearance of the ionic terms relative to the orthogonal basis. Let us assume that the states | r_m ñ and | l_m ñ employed above are the symmetrically orthonormalized original states | r_m⁰ ñ and | l_m⁰ ñ so that:

æ ç è r l ö ÷ ø = S-1/2 æ ç è r0 l0 ö ÷ ø ,

(0)

where the 2×2 overlap matrix S contains 1 as diagonal elements and the overlap S as off-diagonal ones. The original non-orthonormal states are restored form the orthonormal ones by applying matrix S^1/2. The normalized HL state:

1 Ö 2(1+S2) | rma0+lmb0++lma0+rmb0+ ñ

(0)

can be readily transformed to the form Eq. (13) taking into account that this transformation is a direct product of two S^1/2 transformations for the orbitals. It means that the amplitudes are determned by the elements of the matrix S:

| rma0+lmb0++lma0+rmb0+ ñ = Srma+rmb++rma+lmb++lma+rmb++Slma+lmb+,

(0)

where both sides have the norm-square 2( 1+S²) .

The symmetric form of the geminal representation allows to extract the effective value z_S of the Anderson's parameter z parameter as a function of the intrabond overlap, which is equal to:

zS = 2S 1-S2 .

(0)

This parameter determines all the geminal amplitudes and the density matrices according to the formulae given above. To test this model we performed a series of calculations similar to those described in the previous subsection. The numerical results are presented in Table . This Table demonstrates that the values of z determined on the basis of the overlap are significantly smaller than the values arising from the effective bond Hamiltoniani, although they often have similar trend when going from one bond to another. At the same time, the numerical results within the HL model are relatively poor and this model can be used for the construction of the mechanistic picture for the potential energy surfaces only if some (uniform) scaling of the HL based estimates of density matrix elements is performed.

Bond asymmetry

So far our consideration was limited to the properties of either symmetric or symmetrized effective bond Hamiltonians. At the same time, the difference between the estimates obtained in this manner and the actual values of the density matrices can be significant (although not exceedeing 5% for the Coulson bond orders and 20% for the two-electron density matrix elements) in the case of polar bonds. Qualitatively the symmetrized model obviously does not reproduce at all the difference between the diagonal one-electron density matrix elements P^rr-P^ll. In the semiempirical context Tch73 the bond asymmetry has been taken into account perturbatively. It was natural, since a single unaccounted parameter l measured the interaction between the symmetric and antisymmetric states and the entire effect of the bond antisymmetry could be loaded upon it. Actually, instead of the parameter l another dimensionless parameter m had been introduced []:

m = l G(z) .

(0)

This parameter is small and it naturally arised from the first-order correction to the bond polarity:

Prr-Pll = 2m G(z)-1 G(z)+1 ,

(0)

where the factor 2 is due to the spin summation in Eq. (76) and the rest of the multiplier at m is smaller than 1 and tends to 1 for large z. It was possible to single out two types of contributions to the parameter m. It can be broken down into a sum of hybridization dependent component m₀ corresponding to the bond itself and the rest m₁ describing the environment of the bond:

m = m0+m1,

(0)

where the first part prevails by orders of magnitude.

In the ab initio context the situation is somewhat more complex. First of all, the interaction of the symmetric ground state with the excited antisymmetric one leading to the asymmetry of the bond does not reduce to the difference of the diagonal matrix elements of the original effective Hamiltonian (reflected by l), but contains also a contribution coming from the asymmetry of the off-diagonal Coulomb interaction matrix elements (expressed by c). Next, the energy denominators entering the perturbation theory expansions are modified through the exchange matrix element. Nevertheless, it is possible to generalize the principal stability (transferability) results obtained previously in the semiempirical context also to the ab initio setting. If we number the ground, symmetric excited and antisymmetric ionic states of the symmetrized bond hamiltonian as states 0, 1 and 2, respectively, then the energies and the matrix elements of the perturbation can be written (omitting the constant C in the Hamiltonian Eq. (70) and expressing everything in units of Dg):

e0 = 1-G(z) 2 ,    e1 = 1+G(z) 2 ,    e2 = 1- w 2 , \notag (0) [-2mm] (0) V02 = ccosj-lsinj 2 ,    V12 = - lcosj+csinj 2 . \notag (0)

Using these notations it is easy to write the perturbation theory corrections to the wave function and the density matrix elements. The first-order correction to the bond polarity becomes:

Prr-Pll = 2(m-m¢) G(z)-1 G(z)+1-w     ,

(0)

where the contribution m^¢ comes from c and is defined as:

m¢ = cz G(z)(G(z)-1) .

(0)

The actual values of the parameters m and m^¢ are given in Table . The parameter m^¢ is usually much smaller than m (the only exception is the molecule LiH) where they are of the same order.

The analysis in the semiempirical case [] demonstrated that the off-diagonal density matrix elements P^rl and G^rrll deviate from their ideal" values 1 and [1/2], respectively, in the second order with respect to the bond asymmetry perturbation. A similar result can be obtained in the ab initio casei, too. However, it is not possible to consider (m-m^¢) as a single parameter of the bond asymmetry because l and c enter the matrix elements V₀₂ and V₁₂ in different combinations. Nevertheless, the corrections remain quadratic functions with respect to them (or, equivalently, with respect to m and m^¢), although the expressions become somewhat cumbersome. This result is important because it shows the degree of the transferability for different density matrix elements. Moreover, taking into account that the corrections due to c and w are typically small (and of opposite signs), the equations derived in the semiempirical domain can be good approximations to the problem. We performed the calculations of the density matrix elements employing the second-order corrections with respect to m and m^¢ (In fact, it actually represented the first-order correction as far as bond polarity is conserned.) Our calculations (Table ) show that this approach produces surprizingly accurate results even for strongly polar bonds, where the asymmetrical perturbation is not small.

Non-perturbational estimate of stability (transferability)

The previous analysis of stability (transferability) with respect to the bond polarity is based on the perturbation estimates. At the same time it is possible to propose an alternative approach to the problem of establishing transferability for the density matrix elements. Let us consider the following parametric form for the geminal amplitudes:

u = cosfsinr, \notag (0) v = cosfcosr, (0) z = sinf. \notag (0)

Our previous computations [] show that the values of f and r should be close to p/4 for slightly polar and slightly correlated chemical bonds. Therefore, we use a natural transformation:

f = p 4 +q, \notag (0) (0) r = p 4 +h, \notag (0)

which leads to a new definition of the geminal amplitudes:

u = 1 2 (cosq-sinq)(cosh+sinh), \notag (0) v = 1 2 (cosq-sinq)(cosh-sinh), (0) w = 1 2 (cosq+sinq), \notag (0)

where the parameters q and h are expected to be small. These parameters allow to look on the transferability at a different angle. Namely, the density matrix elements can be written as:

Prl = cos2qcosh\asymp 1-2q2- 1 2 h2,\notag (0) Prr-Pll = (1-sin2q)sin2h\asymp 2h-4qh = 2h(1-4q), (0) Grrll = 1 2 (1+sin2q)\asymp 1 2 +q. \notag (0)

These equations show the orders of the transferability of the density matrix elements with respect to the presumably small parameters q and h. These relations remind the previously established ones using the parameters z^-1 and m (for example, the bond order is transferable up to second order with respect to both parameters).

The transferability properties of the density matrix elements are demonstrated by the above expressions provided the actual values of the parameters q and h are indeed small. Thus, it is necessary to have some reasonable estimates of these characteristics through the parameters of the problem (effective Hamiltonian). To this end we write the energy as an approximate function of q and h omitting all the contributions of the order higher than second with respect to these parameters. This energy function can be explicitly minimized with respect to its variables. As a result, we obtain the following expressions of the sought parameters:

q = z-l2+lc-w 2(z2-wz-l2) , (0) \notag (0) h = zl-zc-l z2-wz-l2 . \notag (0)

We performed the actual calculations of the density matrix elements based on these estimates. The results are presented in Table . The estimates are good in the case of bonds with small polarity (like C-H) but the quality of the estimates deteriorates in the case of polar bonds. It means that the estimates made by making use the perturbation theory are better for the practical calculations. Of course, the somewhat worse results obtained in this case characterize rather the quality of the estimates Eq. (104) than the transferability properties of the density matrix elements.

Discussion

In the present paper we attempted to extend the transferability of density matrix elements on organic compounds of different classes (alkanes, alcohols, amines etc.) established in Refs. [,] in the semiempirical context to the ab initio Hamiltonian. This can be considered as a prerequisite for developing a mechanistric model of molecular potential energy surface which in its turn is consistent with ab initio model of molecular electronic structure. It turns out to be possible provided tha wave function of the molecular system can be presented as an antisymmetrized product of strictlylocal geminals. Under this condition we have demonstrated a remarkable stability of all the geminal related ESVs. The values of the polarity P_m^rr-P_m^ll do not exceed 0.07 by absolute value for the compounds containing carbon, nitrogen, and hydrogen atoms. Also the ionicity u_m²+v_m² for a rich variety of bonds had a stable value about 0.4. These features although were not completely unexpected, since the transferability of the parameters of the single bonds in organic compounds is well known experimentally, required a theoretical explanation. This explanation had been given by the observation that the values of the density matrix elements have to be close to their transferable values:

G0 = 1 2 , P0 = 1.

(0)

which depend neither on molecular geometry nor on the nature nor on hybridization of the atoms bonded within a narrow range. The same as we see now applies to the estimates of the similar quantities in the ab initio context. The small corrections to the transferable values of the density matrix elements can be expressed in terms of the effective bond Hamiltonians characterizing the amount of deviation of correlation and of asymmetry of the considered bond as compared to the transferable values. In the semiempirical context the above observations allowed to develop a series of approximate expressions for the molecular energy functional written for the SLG trial wave function. The approximate formulae for the density matrix elements suggest several options: one may neglect all the bond specific corrections and by this to arrive to the picture where all the density matrix elements in the basis of orthogonal HOs are fixed at their transferable values. This is termed as the FA (fixed amplitudes where the geminal amplitudes u,v,w for all geminals are meant). In this case the energy becomes the function only of molecular geometry and of the rest of the ESVs of the molecular system which describe the shape and the orientation of the HOs spanning carrier subspaces for different geminals (bonds). Alternatives to the the FA picture are numerous. All of them can be uniformly termed as TA (tuned amplitudes) approximations. Only the simplest ones were explored in Refs. [,,,]. They can be characterized as those neglecting any nontransferability corrections of the order higher than one in any of the small parameters. It is reasonable to assume that similar moves can be applied in the ab initio context thus paving the route of a sequential transition from a strict quantum chemical description of molecular electronic structure of now accepted ab initio level to its classical description in terms of bonding and nonbonding force fields.

Conclusion

In the present paper we explored the transferability features of the approximate SLG wave function in the ab initio context and showed that the important elements of molecular electronic structure: elements of the one- and two-electron density matrices are transferable from one molecule to another to the degree controlled by small composition and environment specific parameters. If the mentioned corrections are taken into account the corresponding matrix elements can be calculated as transferable functions of the molecular integrals pertinent to each of the bonds.

Acknowledgments

This work has been generously supported by the TÉT grant No. RUS-12/2004 and the FASI grant No. 02.434.11.7098 as a joint Russian-Hungarian research project, by the Russian Foundation for Basic Research (RFBR) grants Nos. 05-07-90067, 07-03-01128 and by the Hungarian Research Fund (OTKA) grant No. T43558. A.M.T. acknowledges the financial support from a Marie Curie Fellowship (contract MIF1-CT-2005-021198).

References

[]

V.G. Dashevskii. Conformations of organic molecules. 1973 [in Russian].

[]

U. Burkert, N. Allinger. Molecular Mechanics. ACS Publications, 1982.

[]

C.A. Coulson. Rev. Mod. Phys., 32 (1963) 170.

[]

A.L. Tchougréeff, A.M. Tokmachev. Int. J. Quant. Chem. 96 (2004) 175-184.

[]

A.L. Tchougréeff. J. Mol. Struct. THEOCHEM 630 (2003) 243-263.

[]

A.M. Tokmachev, A.L. Tchougréeff. J. Comp. Chem. 26 (2005) 491-505.

[]

K. Ruedenberg. Rev. Mod. Phys., 34 (1962) 326-376.

[]

I. Mayer. Int. J. Quant. Chem. 23 (1983) 341-363.

[]

V.A. Fock, M.G. Veselov, and M.I. Petrashen'. J. Exp. Theor. Phys., 10 (1940) 723 [in Russian].

[]

E. Kapuy. Acta Phys. Acad. Sci. Hung. 11 (1960) 409; ibid. 12 (1960) 185.

[]

A.M. Tokmachev, A.L. Tchougréeff. J. Comp. Chem. 22 (2001) 752-764.

[]

A.M. Tokmachev, A.L. Tchougréeff. J. Phys. Chem. A 107 (2003) 358-365.

[]

I. Mayer. J. Phys. Chem. 100 (1996) 6249-6257

[]

R. Ponec and D.L. Cooper Faraday Discuss. 135 (2007) 31-42.

[]

P. R. Surján. Second Quantized Approach to Quantum Chemistry. Springer-Verlag, Berlin et al. 1989.

[]

I. Mayer Int. J. Quant. Chem. 70, (1998) 41-63.

[]

J.M. Parks and R.G. Parr. J. Chem. Phys. 28 (1958) 335.

[]

A.M. Tokmachev, A.L. Tchougréeff. Zh. Fiz. Khim. 73 (1999) 347 - 357.

[]

A.L. Tchougréeff, A.M. Tokmachev. Int. J. Quant. Chem. 100 (2004) 667-676.

Footnotes:

¹This contribution is dedicated to Prof. I.V. Abarenkov on occasion of his 75-th birthday