Content-type: text/HTML
Molecular mechanics (MM) [,] is currently a versatile and very popular tool in laboratory and industrial practice. It remains the most practical way of analyzing potential energy surfaces (PES) for large molecules even though linear scaling O(N) methods of quantum chemistry are becoming available [,]. The nature of the MM as a combination of totally empirical classical force fields allows to realize its main drawbacks (for example, inapplicability to highly correlated systems) and advantages (fast evaluation of the total energy and its gradient and high accuracy of molecular geometries obtained) and to characterize them as the most preferable according to the "quality/cost" criterion but with the field of application limited to certain combination of properties, processes, and classes of molecules.
Despite its long history and a wealth of successful applications the MM approach remains not substantiated theoretically. It is clear that the molecular PES can be expanded up to the second order in nuclear displacements near the equilibrium geometry and this fact is often considered as a general substantiation of the possibility to use the MM-type expansions for the total energy []. However, it is very questionable since (i) the transferability of the elements of the dynamic (second derivatives) matrix between molecules is not proven, and (ii) the geometry variables (bond lengths and valence angles) used in the definitions of the force fields do not diagonalize the dynamic matrix even approximately, whereas the standard MM PESs are by construction diagonal in the bond-length/valence-angle representation.
Substantiation of MM is not merely an academic question. The last years have demonstrated a growing interest [,,,] to hybrid quantum mechanical/molecular mechanical (QM/MM) schemes where different parts of the entire molecule are treated by different methods. The relatively small part of the molecule (reactive center) is treated by an adequate quantum mechanical (QM) method while the rest (inert environment) is described by classical force fields. Thus constructed QM/MM methods combine the universality of the QM description and the low cost of the MM that makes them an optimal tool for analysis of the structure of large molecules including those of biological importance. Though the QM/MM methods shift the limits for numerical applications they require more insight into the electronic structure underlying the MM schemes since it is not clear how to construct the junction between quantically and classically treated subsystems. Though many ad hoc recipes for the junction construction are proposed in the literature (see, for example, reviews [,,]), the problem has not been yet solved. Meanwhile, the absence of sequential derivation for the intersubsystem junction leads to numerous artefacts in QM/MM schemes reviewed in Ref. []. The general formal solution of this problem proposed in Refs. [,] is based on the assumption that a generic MM scheme can be obtained from a suitable QM method by a series of approximations. Thus the prospect of development of effective QM/MM schemes constitutes the practical aspect of the general task of the MM substantiation. An important prerequisite for the announced derivation of MM from QM is the choice of the quantum chemical method underlying the MM description. By taking a proper form of the trial wave function having some common features (or resemblance) with the MM schemes we can significantly simplify the task. The representation of the total energy of the covalently bound molecule in the current MM schemes is tightly connected to the concept of two-center bonds. The energy is taken as a sum of intrabond, bending, torsion and non-bonding contributions. The specific forms of force fields depend on the implementation and quite a lot of schemes with systems of parameters fitted to reproduce various characteristics of different classes of molecules have been proposed [,,,,,,]. In the course of the evolvement of the MM approach increasingly sophisticated contributions are added to the original simplistic picture [] which allow to extend the MM approach to more and more complex and not transparently tractable classes of the molecules [,,]. Thus an acceptable substantiation of the MM should not be reduced to deriving or validating any specific MM scheme. By contrast, a generic mechanistic picture must be obtained on the basis of the adequate QM description.
The standard QM methods in that or another manner based on the delocalized picture of molecular electronic structure provided by the SCF approximation are not quite suitable for the purpose of the MM substantiation. First, the general SCF energy expression contains contributions of the same form for all pairs of atoms and separation of the energy into bonding and non-bonding contributions is not built in to the SCF methods. Second, even applying different localization techniques intended to reconcile the SCF picture of the electronic structure with the chemist's view relying upon the bonds concept leads not to very much progress due to "tails" of the localized orbitals which are difficult for handling and absolutely non-transferable (it is also impossible to prove theoretically the transferability of the MO LCAO coefficients or, even, of the SCF density matrix elements).
The above discussion allows to specify some criteria for selecting a QM method suitable to serve as a starting point for the MM derivation: discrimination and different treatment of bonding and non-bonding contributions to the energy, strict locality of one-electron states, variational determination of the electronic structure parameters (ESPs). Failure to fulfil any of these criteria leads to a necessity to take whichever of the mentioned features as assumptions (as for example, in Ref. [], where the analysis was based on the PCILO method []).
The unsuitabilty of the SCF approximation as a basis for the additivity concept in general (and thus for the MM in particular) was clear even on very early stage of the theory development. Yet then (see []) it had been proposed to use a geminal based description for substantiating the additive methods and the bond energy concept. However, the one-electron carrier spaces to be used for geminal consctuction had not been sufficiently specified and the entire geminal-based scheme had not been explored to a due extent by adequate numerical studies. Only recently we proposed [,] a semiempirical method which satisfies the suitability criteria formulated above. It is based on the trial wave function in the form of antisymmetrized product of strictly local geminals (SLG). Variational determination of strictly local one-electron states and of the density matrix elements provides the flexibility of the wave function necessary to describe the electronic structure and properties of organic compounds. Numerical estimates based on the MINDO/3 parameterization of the semiempirical Hamiltonian [,] and three (MNDO, AM1, and PM3) parameterizations of the NDDO family [] have shown that the method proposed supercedes the SCF one in description of the heats of formation, molecular geometries, and ionization potentials of organic molecules. Incidentally, it provides a good basis for description of molecular electronic structure in chemical terms. The method is well suited for description of large molecules since it belongs to a class of O(N) methods and can be used for constructing hybrid schemes []. Also the approximate transferability of the ESPs of the SLG approach was numerically demonstrated in Ref. [], where it was shown that the energy calculated with some ESPs "fixed" and characteristic for a particular class of atoms and for bonds is close to that obtained by direct minimization of the SLG energy. More refined treatment was proposed in Refs. [,], where linear response relations for the form of one-electron states were used to derive explicit expressions for the angular dependence of the energy in the case of sp3 hybridized carbon and nitrogen atoms. The numerical estimates of the MM force field parameters obtained in [] on the basis of the analytical expressions for the constants of force fields are close to those accepted in the MM. Special attention was paid to off-diagonal force fields and to possible sources of the angular dependence of the energy [] as well as those of piramidalization potential in nitrogen-containing compounds []. The development of these ideas has led to the formulation of the "deductive molecular mechanics" [], which represents atoms by their hybridization tetrahedra (see below) of different shapes adjusted according to geometry and valence state variations rather than by harmonically interacting point masses ("balls-and-springs") used in standard MM. All these results are obtained within the assumptions of (i) the perfect transferability of ESPs characterizing chemical bonds and (ii) of the validity of linear response relations for hybridization with respect to geometry variations. In the present paper we address the problem of transferability of the density related ESPs, explore the precision and the validity limits for the linear response formulae for the shapes of hybridization tetrahedra and consider the possibility of recovering the standard MM description from the deductive molecular mechanics.
The paper is organized as follows: in Section we briefly review the main features of the SLG method (with the MINDO/3 semiempirical parameterization) relevant to constructing the MM description. In Section the ESPs characterizing the density distribution in chemical bonds and lone pairs as they appear in the SLG method are analyzed and the features assuring their transferability are singled out. In Section we briefly discuss the structure of the hybridization manifold and give expressions for the variations of the hybridization related ESPs in response to geometry changes. Results of numerical experiments are given in these Sections to support our theoretical derivation of transferability and of linear response relations. In Section we discuss different possible approximate descriptions of the ESPs and energy leading to various versions of the deductive MM theory. In Section we consider the possibility of deriving the standard MM picture from the DMM by projecting out excessive (from the MM point of view) variables of the latter and give theoretical expressions for the force fields' parameters of the standard MM. Finally we discuss the relation between the transferability of the density related ESPs and that of the MM force fields.
Constructing the SLG trial wave function according to [,] requires the following moves. First, the one-electron basis of the strictly local hybrid orbitals (HOs) must be contructed. These orbitals are obtained by an orthogonal transformation of the s and p AOs for each "heavy" (non-hydrogen) atom. These transformations are represented by 4×4 orthogonal matrices hA Î O(4) for each heavy atom A. All the HOs are assigned either to respective two-electron chemical bonds or to electron lone pairs. Each chemical bond refers to two such HOs - | r ñ and | l ñ (right- and left-end ones, respectively). Each lone pair is formed by one HO only (a right one for the sake of definiteness).
Chemical bonds and lone pairs are described by singlet two-electron functions
- geminals
[] taken in the form originally proposed by Weinbaum [].
With use of the second quantization notation they are written as:
| (0) |
| (0) |
| (0) |
The wave function of electrons in the molecule is then taken as the antisymmetrized
product of the geminals given by Eqs. (1), (2):
| (0) |
The SLG energy is a function of the intrabond matrix
elements of spinless one- and two-electron density matrices:
| (0) |
| (0) |
If the Hamiltonian of the MINDO/3 form [] in the HO basis is used
with the SLG trial wave function the total energy can be written in a form
somewhat close to the MM energy with interactions between bonded and non-bonded
atoms treated in different ways and closely relates it to that given in Ref.
[] in the context of analysis
of a variety of additive schemes of molecular energy:
| (0) |
| (0) |
| (0) |
Molecular integrals entering the above expressions depend on molecular geometry and on the orthogonal matrices hA for all non-hydrogen atoms A. The expressions for the matrix elements in the HO basis are given in Ref. [].
The energy expression Eq. (7) depends on the amplitudes of bond geminals through the values of Pmtt¢ and Gmtt¢, and on the form of HOs through the molecular integrals. The amplitudes of two ionic (um, vm) and one covalent (Heitler-London type, Ö2wm) configurations in Eq. (1) for each geminal are determined with use of the variational principle as well as the matrices hA of transformation of the AO basis to the HO one.
This comprises the essence of the semiempirical SLG method.
In the previous Section we reviewed briefly the semiempirical implementation of the SLG method for analysis of electronic structure and expressed the total molecular energy in the form Eq. (7) which allows the representation of the molecular PES as a sum of local increments. These increments depend on the ESPs of two classes (i) those defining the hybridization of atomic basis sets and (ii) the intrabond density matrix elements. In this Section we concentrate upon the proof of transferability of the electron density matrix elements, related to the geminal amplitudes, and on the structure of the hybridization manifold as it appears in the SLG approximation.
As it is reported in Section 1 the energy in the SLG approximation is a function of one- and two-electron density matrices. Their matrix elements are in turn expressed through the geminal amplitudes, appearing while diagonalizing the effective bond Hamiltonians. Thus any analysis of the properties of the density ESPs starts from description of the latter.
Within the original SLG approach [,] the geminals are
characterized by the amplitudes (see Eq.
(1) um, vm, and Ö2wm = zm, which
simplifies the normalization condition Eq. (3) for the amplitudes
to: um2+vm2+zm2 = 1. To determine them the effective Hamiltonians for each bond
geminal are constructed. The optimal values of these amplitudes are the
solutions of the eigenvector problem (see also []):
| (0) |
The matrix elements of the effective bond Hamiltonians are defined as (with the MINDO/3
Hamiltonian):
| (0) |
| (0) |
The calculations of Refs. [,] performed on organic compounds of different classes (alkanes, alcohols, amines etc.) have demonstrated a remarkable stability of all the geminal related ESPs. The values of the polarity Prrm - Pllm do not exceed 0.07 by absolute value for the compounds containing carbon, nitrogen, and hydrogen atoms (for the situation with oxygen and fluorine see below). Also the ionicity u2m + v2m for a rich variety of bonds has a stable value about 0.4. The bond orders 2Prlm all acquire values between 0.92 and 1.0. These features though not completely unexpected, since the transferability of the parameters of the single bonds in organic compounds is well known experimentally, require a theoretical explanation.
In order to provide the required explanation we notice that the effective Hamiltonians for the bond geminals can be represented as a sum of the unperturbed part which when diagonalized yields an invariant, i.e. exactly transferable, values of the ESPs and of a perturbation responsible for specificity of different chemical compositions and environments.
Pseudospin operator of the bond geminal
Let us introduce a pseudospin [^(t)]m operator corresponding to the
pseudospin value tm = 1. The matrices of its components
in the basis of the configurations defining the geminal are given by:
| (0) |
| (0) |
| (0) |
| (0) |
Perturbative estimate of ESPs with respect to noncorrelated bare Hamiltonian
One can try to estimate the optimal values of the ESPs specific for each bond and
molecule perturbatively by using the linear response approximation []. According
to the latter the response d
á A
ñ of a quantity described
by the operator A to the
perturbation lB of the
Hamiltonian (where l is the parameter characterizing the intensity of
the perturbation) has the form:
| (0) |
| (0) |
| (0) |
Perturbation of the density matrix elements for correlated ground state
In order to overcome the above failure of the perturbative estimation of the
two-electron density and of the bond orders let us consider
a symmetric bond. This would correspond to a different
decomposition of the effective bond Hamiltonian than that of Eq. (15).
We assume that the contribution to the effective bond Hamiltonian which is
proportional to [^(t)]zm2 is included into the unperturbed (zero
order) Hamiltonian. The problem then reduces to a 2×2 matrix
diagonalization. The ESPs, as they appear from solution of this
problem, are:
| (0) |
| (0) |
| (0) |
Now, when the total ionic contribution to the geminal is calculated exactly
(variationally), the bond polarity can be estimated perturbatively in the
linear response approximation, but with the correlated ground state of the
symmetric effective bond Hamiltonian taken for evaluating the Green's function. It
can be conveniently done with use of a dimensionless asymmetry parameter:
| (0) |
| (0) |
The expression for the bond polarity coincides with that of Eq. (24)
even if the second-order perturbation correction to the wave function is used
(i.e., the contribution to the bond polarity proportional to mm2
is absent). At the same time the second-order corrected bond ionicity and bond
order have the following form:
| (0) |
Another archetypical form of two-electron group is the lone pair. As it is
mentioned above the lone pair is described by a degenerate geminal containing
the contribution of only one ionic configuration. For the sake of definiteness
we set it to be the right-end ionic configuration of the corresponding
degenerate bond (the amplitude um becomes equal to unity, see Eq.
(2)). The ESPs related to the lone pair can be readily evaluated:
| (0) |
The above analytical results must be controlled by numerical estimates in order to get a feeling of the real sense of the "first" and ßecond" orders. Table represents the results of calculations on the ESPs á [^(t)]zm ñ , á [^(t)]zm2 ñ , and á [^(t)]+m ñ by the SLG method (Eq. (10)) and by the approximate formulae Eqs. (19), (20), (22), (24), and (25) for some characteristic bonds in small molecules. The results show that in the case of bonds with small polarity all the formulae perform very well. The most precise approximations Eqs. (24) and (25) give results which perfectly coincide with the exact (SLG-MINDO/3) ones even for very polar O-H and F-H bonds. Also estimates according to the asymptotic (zm >> 1) formulae Eq. (20) give reasonable results for the ESPs of the bonds in not too polar molecules at their equilibrium geometries. The main source of stability of the bond order values is the validity of the above limit which in its turn takes place due to the fact that the difference between one- and two-center electron-electron repulsion integrals (Dgm) at interatomic separations characteristic for chemical bonding is much smaller than the resonance interaction at the same distance. The data of Table illustrate the difference between that which may be called MM atom types. For example, the primary C-H bonds in the ethane and propane molecules have very similar ESPs at the same time somewhat differing from those for the secondary C-H bonds in the propane molecule.
Further analysis of the quantities mm allows to single out two types of
factors loaded upon this parameter: those related to the bond itself (which are
again hybridization dependent) and the rest describing the environment of the
bond. These factors contribute additively:
| (0) |
| (0) |
The contribution
to the bond asymmetry
coming from the environment of the bond is:
| (0) |
In the framework of the SLG scheme the structure of one-electron basis states is defined by orthogonal transformations of AOs for each atom with an sp-valence shell. The energy expression Eq. (7) is the function of the parameters defining these transformations. The 4×4 O(4) matrix hA of transformation from the AO to the HO basis set on the atom A depends on six angular variables. Three of them (pseudorotation angles [(w)\vec] b = (wsx,wsy,wsz) with subscripts indicating pairs of basis AOs mixed by the corresponding 2×2 Jacobi rotations) define the structure of the HOs (s-/p-mixing and relative directions of the HOs) while other three (quasirotation angles [(w)\vec] l = (wyz,-wxz,wyz)) define the SO(3) matrix performing rotation of the set of four HOs as a whole (prefix quasi refers to the fact that no physical body rotates under its action, only the system of HO's). Generally, the transformation of orbitals caused by a pseudorotation forms a set of HOs which is known as hybridization pattern (like sp3, sp2 etc.) which is more or less stable, while the set of quasirotation angles is totally non-transferable, depends on the relative placement of bonded atoms and, obviously, is governed by the resonance contribution to the energy since only the latter depends on the directions of the HOs.
The mathematical description of hybridization is based on employing the algebraic group structure of the hybrids'
manifold. Due to the latter any small variation of HOs in a vicinity of a given
set of HOs represented by a 4×4 orthogonal matrix h can be expressed
with use of the SO(4) matrix H close to the unity matrix:
| (0) |
| (0) |
The analysis of the properties of the ESPs pertinent to the SLG approximation
performed in Sections 0.1 and 0.1 allows to
rewrite the energy Eq.
(7) as follows:
| (0) |
Further
components of the description are those related to the HOs. The latter enter into
the theory through the Hamiltonian matrix elements in the HO basis.
The matrix elements entering Eqs. (7), (32) are either
invariant with respect to basis transformations (the interatomic Coulomb
interaction gAB) or can be uniquely expressed through contributions
of s-AO to the HOs (the one-center matrix elements). The only class of
molecular integrals depending on the whole structure of the HOs (including
directions) is that of the resonance integrals. As we mentioned in Section 0.1
each spx-HO can be considered as a normalized quaternion
(s,[(v)\vec]). Following [] we represent the entire system of HOs at
any given atom by four vector parts [(v)\vec] m of the corresponding orthonormal
quaternions. Even this representation is superfluous since only six Jacobi angles
suffice to describe the system of HOs of each given atom completely. Nevertheless,
usage of the vector parts is visual. If the latter are assumed to have the
corresponding nucleus as their common origin the tetrahedral shape thus obtained
contains (with an excess) all necessary information about the system of HOs
of the given atom. In [] such a construct was called the hybridization
tetrahedron of the heavy atom at hand. Using the hybridization tetrahedra as elements
of the theoretical construct allows further discrimination of possible approximations.
Due to the mentioned dependencies of the molecular integrals on the Jacobi angles both
the FA and TA approximations to the energy Eq. (32) depend
on the relative orientation of the hybridization tetrahedra through the bond
resonance integrals brmlmRmLm. The resonance integrals depend also on
the shapes (relative weights of the s- and p-contributions to the HOs which
ultimately define the interhybrid angles) of the hybridization tetrahedra. All other
terms in Eq. (32) in the FA and TA approximations depend only on the shapes
of the hybridization tetrahedra. This leads to the possibility to either fix the
relative weights of the s- and p-orbitals (FO i.e. fixed orbitals
approximation) at spn (n = 1 ¸3) or any other allowable values and by
this fix the shapes of the hybridization tetrahedra which thus become interacting rigid
bodies or to allow the relative weights of the s- and p-orbitals to be tuned
thus leading to the TO - tuned orbitals - picture of the flexible
hybridization tetrahedra. Whichever combination of the FA or TA treatments for the
density matrix elements on one hand with the FO or TO treatments for the HOs on the
other hand results in a representation of the molecular energy Eq. (32)
as such of the system of tetrahedral bodies (rigid or flexible) whose interactions and
self energies depend on distances between their centers, their shapes and relative
orientations. For example, the energy variation due to small pseudo- and quasirotations
of the hybridization tetrahedron in the vicinity of the
equilibrium for sp3-hybridized carbon atom in the FA approximation is given
by a diagonal quadratic form:
| (0) |
The content of the deductive molecular mechanics as formulated in Ref.
[] and above is a description of the
molecular energy in the form of Eq. (32) as a function
of shapes and mutual orientations of the hybridization tetrahedra and of
geometry parameters. On the other hand the standard MM can be qualified as a
scheme directly parameterizing the molecular energy as a function of molecular geometry
only. From this point of view the Jacobi angles variables [(w)\vec] b, [(w)\vec] l
describing the shapes and orientations of hybridization tetrahedra
are superfluous
and must be excluded. This can be done by finding the
response of the corresponding ESPs to the variations of bond lengths and valence angles
with use of linear response relations between different
subsets of variables pertinent to the DMM picture. To do so let us consider
a minimum of the energy with respect to both geometry and the ESPs.
In the vicinity of a minimum
the energy can be expanded upto second order with respect to
nuclear displacements q and variations of the ESPs x:
| (0) |
| (0) |
The main use of the formulae Eq. (35) is for exclusion of the angular variables describing the hybridization tetrahedra from the mechanistic picture. Now we estimate the precision of the linear response relations (Eq. (35)) between geometry and hybridization variations themselves by numerical study of elongation of one C-H bond and deformations of valence angles. We consider tetrahedral methane molecule as a reference (its parameters then correspond to subscript 0 in Eqs. (34), (35)). First of all, we notice that the ÑxÑxE matrix further simplifies for methane since sLmm = 1 and, therefore, simple analytical expressions become possible. Also we remark that the FA approximation is adequate here since, for example, even very large elongation of one C-H bond by 0.1 Å leads to changes of the bond geminal amplitudes u,v, and w not exceeding 0.003. The same applies to the averages of the pseudospin ([^(t)]) operators.
Linear response of hybridization to bond elongation
Let us consider first the relation between hybridization
and elongation of the C-H bond. For this end we need the mixed second order
derivatives coupling the bond
stretching with the hybridization ESPs. For every C-H bond in methane we can
introduce diatomic coordinate frame with the z axis directed along the
bond and express the resonance integral as:
| (0) |
| (0) |
Formula Eq. (37) gives the
analytical expression for the coupling coefficient between the bond elongation
(in Å ) multiplied by unit vector of this bond direction and changes of
pseudorotation angles d[(w)\vec] b in methane (in radians). Its numerical
value C1 is 0.2764 rad·Å-1. This distortion corresponds
to the following form for the matrix of small transformation of the whole set of HOs
(matrix H in Eq. (30)):
| (0) |
Linear response of hybridization to valence angle deformation
The linear response relations between the molecular shape and the shape of hybridization
tetrahedron are rather tricky due to complex structure of the
hybridization manifold.
The molecular shape can be formally characterized by unit vectors with origin
at an atom considered, taken as a center, and pointing to those bonded to the
central one.
In the case of methane
the deformations
of thus defined coordination polyhedron are small rotations of unit vectors [(e)\vec]m
directed from the carbon atom to hydrogen atoms.
Their small rotations d[(j)\vec] m form an 8-dimensional space which decomposes to
a direct sum of two subspaces: one 3-dimensional corresponding to rotations
of the molecule as a whole and another 5-dimensional corresponding to independent
variations of valence angles. The former one is precisely mapped on the 3-dimensional
space of quasirotations d[(w)\vec] l
while the latter (5-dimensional) one is mapped on the
3-dimensional space of pseudorotations d[(w)\vec] b corresponding to
changes of the shape of
the hybridization tetrahedron [].
Due to very general theorems of linear algebra [] there exists at
least a two-dimensional kernel in the space of deformations of molecular shape which
maps to zero deformation of hybridization tetrahedron. In [] the term
"hybridization incompatible" has been coined for the deformations from this kernel. The
structure of
deformations laying in the kernel of the mapping is quite simple: they
are produced by equal variations of opposite (spiro) valence angles. In
contrast, the variations which correspond to increase of one valence angle by
dc and decrease of its spiro counterpart by the same value fall into "coimage"
of this mapping i.e. to the subspace which one-to-one maps to the space of
pseudorotations d[(w)\vec] b. The deformations in the coimage
can be
called "hybridization compatible". It is clear that only these latter
variations should be considered. It is also clear that any variation of the valence
angle is a sum of equal amounts of
hybridization compatible and hybridization incompatible deformations.
The denominator in the linear
response relation Eq. (35) is the same as for Eq.
(37) while the relevant block of the ÑxÑqE matrix
(with q taken as a difference of two opposite valence angles) is proportional
to bzsCH. Applying the linear response technique to the
"hybridization compatible" variation of two spiro valence angles we obtain:
| (0) |
The coupling coefficient between the change of the pseudorotation vector and
totally hybridization compatible deformation of valence angles can be easily
found by Eq. (39). Its numerical value C2 is -0.20734 for the equilibrium
interatomic distance. The considered distortion produces the following HO
transformation matrix (matrix H in Eq. (30)):
| (0) |
The smallness of the coupling coefficient C2 even for the totally hybridization compatible deformations allows to qualitatively understand certain features of the electronic structure of cyclopropane as it appears in the SLG approach: a very large distortion of the C-C-C valence angle from the tetrahedral to 60° one leads only to a relatively small distortion of the corresponding interhybrid angle. We model this process by strongly deforming the methane molecule. Simple estimate is based on Eqs. (39), (40) and runs as follows. The valence angle deformation when going from methane to cyclopropane is of 49.5° (=109.5°-60°); only one half of it is hybridization compatible; after multiplying by C2 this yields the value of the interhybrid angle between the HOs corresponding to the üntouched" C-H bonds of 114.5° (i.e. the angle variation amounts only 5°). From the energy minimum condition for hydrides it follows that the C-H bonds indeed must follow the directions of the HOs. Numerical experiments performed with use of the SLG-MINDO/3 method show that if one of the H-C-H valence angles is fixed at the cyclopropane value of 60° the energy minimum corresponds to its spiro counterpart of 115°. These results can be directly compared with the experimental H-C-H valence angle in the cyclopropane molecule which equals to 115.1°.
Analogous estimate can be applied to cyclobutane. In this case we consider the distorted methane molecule with angle 90°. The response of the HOs to the deformation is proportional in our model to deviation of the valence angle from the tetrahedral one. The deviation of the C-C-C angle from the tetrahedral one in cyclobutane (19.5°) amounts 40% of that in cyclopropane. Therefore, we can expect that about the same ratio will be observed for the deviations of the H-C-H valence angle from the tetrahedral one in the cyclobutane and cyclopropane molecules. In fact, this ratio in the SLG-MINDO/3 numerical experiment is about 39%.
Announced transition from the DMM model of molecular PES to a model dependent on molecular
geometry is formally obtained by inserting Eq. (35) to Eq. (34) which yields
| (0) |
| (0) |
| (0) |
We consider in more details the energy curve for the C-H bond. The curve
corresponding to the sp3 hybridization of carbon atom and to the symmetric
TA picture (Eq. (20)) is given by Fig. 3. It has correct qualitative
behaviour for all interatomic separations. The minimum depth on this curve is
approximately -0.23 a.u. and can be considered as the "pure" energy of the
C-H bond. It is generally accepted in the literature that the energy of
C-H bond is approximately 0.15 a.u. At the same time the latter value is
thermodynamic one while our value is obtained by extracting the contributions
to the energy intrinsic to this bond and excluding the interaction between the
bonds. The difference between the thermodynamic value for the bond energy and
that obtained from the SLG energy in the FA picture can be explicitly written
in quite simple form:
| (0) |
It is interesting to compare the form of the bond energy curve Fig. 3 with the Morse potential. To this purpose we tried to approximate the curve of Fig. 3 by the Morse function D0[1-exp(-a(r-re)/re)]2 by minimizing the area between two curves in the interval from 0.72 Å to 2.50 Å . With the parameters D0 and re fixed at the values equal to the minimum depth and position on the curve (0.2295 a.u. and 1.078 Å ) the optimal value of parameter a is then 2.306 but with these parameters two curves are in fact quite different (the area between curves is almost 11% of area between the bond energy curve and the abscissa). If we optimize all three parameters of the Morse curve they become slightly modified D0=0.2333 a.u., re=1.045 Å , and a=2.295. This reduces the area between the curves by 30%. It should be concluded that the energy profile in the TA approximation is not partcularly well reproduced by whatever Morse curve.
In order to estimate the parameters of harmonic
force fields we consider the symmetric correlated single bond, where the
energy can be obtained without any reference to its
environment. In our case the derivative of the bond energy with respect to
a geometry parameter q has the form:
| (0) |
The same concepts can be used to determine the elasticity constant for the
bond stretching by taking the second derivative of the energy with respect
to the bond length. In the FA picture we get:
| (0) |
| (0) |
Analogous treatment of the energy terms quadratic in
valence angles' deformations yields the bare estimate for the harmonic
bending constant in the form:
| (0) |
In the previous Section we provided the exclusion of the angular variables characterizing the shapes and orientations of the hybridization tetrahedra from the mechanistic DMM model of molecular PES. This results in a model announced in the Introduction, which is similar to the standard MM models but is obtained by the sequential derivation from the QM (SLG) model of molecular electronic structure. As it is mentioned the transferability of the ESPs characterizing chemical bonds in molecules and linear response relations for hybridization ESPs are main components of deriving MM theory of molecular PESs from corresponding QM theory. Both these features have been mathematically derived and numerically checked in Section 1.
Despite its long history the very term "transferability" remains somewhat vaguely defined synonym of äll the best" in parameterization schemes, referring largely to their capacity to be used without change for any molecule in a sufficiently wide class of similar ones. From quantitative point of view this concept have got some attention in two related areas. First we mention the estimates of transferability given in Ref. [] where that of the semiempirical quantum chemical parameters has been related to the fact that the corresponding quantities remain the same for all molecules of similar structure upto the second order with respect to overlap integrals between AOs residing at neighbour atoms. That allowed to define the transferability for the quantum chemical parametrs (ultimately, for the Hamiltonian matrix elements) as invariance of some quantity to a given order of precision with respect to a small parameter. Analogously in Ref. [] the problem of constructing transferable dynamical matrices in relation to analysis of vibrational spectra has been considered. The stability of the dynamical matrix was analyzed with respect to small parameter of relative mass variation under the isotope substitution in a series of related molecules.
The importance of the transferability of the geminals has been pointed out yet in []. It was stated that the assumption of the transferability of the geminal amplitudes is a prerequisite for that of the bond energy. However, in [] the geminal transferability had not been shown and the authors concentrated on the statements equivalent to the transferability of the MM bond stretching force fields. The transferability of the MM force fields must be considered as an important chracteristic of this approach. The reasons to treat that or another force field as a transferable between two specific molecules or classes of molecules are either purely pragmatic or this question is solved on the school-wise grounds []. As far as we know in the literature there were no attempts reported to prove this property of the force fields from any general point of view. Here we present a step towards quantitative analysis of the transferability of the MM force fields by proving the transferability of the density matrix elements - equivalent to that of the geminal amplitudes, but more directly related to the energy. Under the assumptions given by Eqs. (15), (16) the averages of the pseudospin operators (and thus all the bond ESPs) are invariant in that sense that they do not depend on environment of the bond under consideration and even on particular composition of the bond, i.e. on the nature and the hybridization of the atoms the bond connects. This corresponds to the FA picture (see above and Ref. []). It is important that the invariance (at the established level of precision) of the density matrix elements can be proven only for the basis of the variationally determined HOs - a specific characteristic of the SLG approach [,]. In the basis of AOs the density matrix elements are not invariant even approximately. Though the approximation sufficient to obtain formally these invariant results (the SCF approximation) is very crude it, nevertheless, breakes only at large interatomic separations which normally are not covered by any MM-like approximation. This result allows to pose further questions: to what extent the density ESPs' invariance may stand further improvements of the description and whether it is possible to relate the invariance of the density matrix elements with the transferability of the MM force fields. To answer these questions we notice that the invariant values of ESPs can be improved by perturbative corrections (the TA picture) reflecting all diversity of chemical compositions and environments the bond may occur in. Nevertheless, all the variety of perturbations is characterized by two small dimensionless parameters: zm-1 Eq. (21) and mm Eq. (23). Both parameters depend on the atoms connected by the bond, their separation, and their hybridization. The perturbative treatment allows to estimate the precision of transferability. For example, using Eqs. (22), (24), and (25) we conclude that the bond order is the quantity transferable upto second order with respect to both zm-1 and mm; the ionicity (the total weight of the ionic configurations) is transferabile upto second order with respect to mm and upto first order with respect to zm-1; the bond polarity is transferable upto first order with respect to both zm-1 and mm. The second order transferability of bond orders explains to certain extent the success of the concept of ßingle bond" suitable for a large variety of chemical bonds. Note that the second order transferability takes place for the bond orders also in the case when we employ the SLG bond wave function with the correct asymptotic behavior despite the fact that the transferable numerical value itself is obtained from the SCF wave function which does not have the correct asymptotic behavior. Within this picture all specific characteristics of the force field are loaded into parameters of the (effective) Hamiltonian, which are numbers specific either for a given atom in certain hybridization state or for a pair of such hybrid states of atoms - ends of the bond. Although the density ESPs can be considered as constants independent on any details of molecular composition or geometry, the force fileds which are basically sums of products of ESPs by matrix elements of molecular Hamiltonian are geometry dependent and composition specific. The force fileds thus obtained are expected to be the same for the same composition of the bond and to depend weakly (to the extent of the variance of the mm1 parameters) on the environment. These properties are basically much more than necessary for substantiating whichever MM-like description.
In the present paper we discussed the problem of deriving the MM representation of the molecular PES from a relevant QM description. Using the SLG wave function we analyzed the ESPs related to bond geminals and to hybridization tetrahedra. It was shown that the bond-related parameters can be represented as functions of parameters of the MINDO/3 Hamiltonian in the HO basis, transferable from one molecule to another. The functional form of the ESPs found is valid at arbitrary interatomic separations. At the interatomic separations close to the equilibrium bond lengths characteristic for the MM-like treatments two approximations, both suitable for substantiation of the ESPs transferability were considered. One is the fixed geminal amplitudes approximation which results in perfectly transferable numbers referring to the ESPs in question. Another, more exact is the tuned geminal amplitudes approximation which takes into account small corrections to the invariant ESPs. Two small parameters characterizing specificity of the bond and effects of its environment were introduced. By this the whole manifold of quantum chemical parameters defining the effective bond Hamiltonian boils down to only two relevant parameters zm -1 and mm. The presence of such only two-dimensional manifold and smallness of the parameters for a wide range of bonds in quite different environments essentially explains the transferability of the density related ESPs. This allows for a family of mechanistic models describing molecular PESs in terms of hybridization tetrahedra with interactions dependent on distances between their centers and on mutual orientations. Linear response relations for variation of hybridization parameters due to elongation of chemical bonds or specific changes of valence angles are considered. They allow to exclude the angular variables describing the shapes and orientations of hybridization tetrahedra and to represent the molecular energy in both the FA and TA approximations as that of the system of interacting point masses [] ("balls-and-springs" picture) depending on the molecular geometry only. This energy has a form of a sum of local (bond) increments corresponding to the force fields of the standard MM. The estimates of the paramemters of these force fields coming from the analytical expressions are compared with those obtained in numerical experiments showing the high accuracy of analytical estimates. The reasons for this possibility are both the transferability of the bond-related ESPs and the linear response relations for the hybridization tetrahedra are numerically tested for their precision.
This work was performed with financial support of the RFBR through the grants 02-03-32087, 04-03-32146, and 04-03-32206. It has been completed during the stay of A.M.T. in the RWTH, Aachen in the frame of the Alexander von Humboldt Postdoctoral Fellowship which is gratefully acknowledged as is the kind hospitality of Prof. R. Dronskowski. A.L.T. gratefully acknowledges valuable discussions with Profs. N.F. Stepanov, A.A. Levin, and I. Mayer.
Molecule | Bond | z-1 | m | á [^(t)]zm ñ | á [^(t)]zm2 ñ | á [^(t)]+m ñ | ||||||||
(10) | (19) | (24) | (10) | (19) | (20) | (25) | (10) | (20) | (22) | (25) | ||||
H2 | H-H | 0.127 | 0.000 | 0.000 | 0.000 | 0.000 | 0.437 | 0.436 | 0.437 | 0.437 | 0.992 | 0.992 | 0.992 | 0.992 |
CH4 | C-H | 0.181 | 0.093 | 0.065 | 0.095 | 0.065 | 0.414 | 0.410 | 0.411 | 0.414 | 0.982 | 0.984 | 0.984 | 0.982 |
NH3 | N-H | 0.172 | 0.097 | 0.069 | 0.098 | 0.069 | 0.419 | 0.414 | 0.415 | 0.419 | 0.983 | 0.986 | 0.985 | 0.983 |
H2O | O-H | 0.161 | 0.341 | 0.244 | 0.345 | 0.247 | 0.464 | 0.419 | 0.420 | 0.466 | 0.959 | 0.987 | 0.987 | 0.959 |
HF | F-H | 0.285 | 0.540 | 0.309 | 0.562 | 0.308 | 0.451 | 0.357 | 0.363 | 0.457 | 0.925 | 0.962 | 0.959 | 0.929 |
C2H6 | C-C | 0.193 | 0.000 | 0.000 | 0.000 | 0.000 | 0.405 | 0.404 | 0.405 | 0.405 | 0.982 | 0.982 | 0.981 | 0.982 |
C-H | 0.181 | 0.071 | 0.050 | 0.072 | 0.050 | 0.413 | 0.410 | 0.411 | 0.413 | 0.983 | 0.984 | 0.984 | 0.983 | |
C-C | 0.193 | 0.016 | 0.011 | 0.016 | 0.011 | 0.405 | 0.404 | 0.405 | 0.405 | 0.982 | 0.982 | 0.981 | 0.982 | |
C3H8 | C1-H | 0.181 | 0.074 | 0.052 | 0.075 | 0.052 | 0.413 | 0.410 | 0.411 | 0.413 | 0.983 | 0.984 | 0.984 | 0.983 |
C2-H | 0.180 | 0.051 | 0.035 | 0.051 | 0.035 | 0.412 | 0.410 | 0.411 | 0.412 | 0.984 | 0.984 | 0.984 | 0.984 | |
Cyclo- | C-C | 0.199 | 0.000 | 0.000 | 0.000 | 0.000 | 0.402 | 0.400 | 0.402 | 0.402 | 0.981 | 0.981 | 0.980 | 0.981 |
propane | C-H | 0.180 | 0.084 | 0.059 | 0.086 | 0.059 | 0.414 | 0.410 | 0.411 | 0.414 | 0.983 | 0.984 | 0.984 | 0.983 |
N2H4 | N-N | 0.290 | 0.000 | 0.000 | 0.000 | 0.000 | 0.361 | 0.355 | 0.361 | 0.361 | 0.960 | 0.960 | 0.958 | 0.960 |
N-H | 0.176 | 0.092 | 0.065 | 0.093 | 0.065 | 0.416 | 0.412 | 0.413 | 0.416 | 0.983 | 0.985 | 0.984 | 0.983 | |
N-C | 0.219 | 0.025 | 0.016 | 0.026 | 0.016 | 0.393 | 0.390 | 0.393 | 0.393 | 0.977 | 0.977 | 0.976 | 0.977 | |
CH3NH2 | N-H | 0.173 | 0.098 | 0.069 | 0.099 | 0.069 | 0.419 | 0.412 | 0.415 | 0.419 | 0.983 | 0.985 | 0.985 | 0.983 |
C-H | 0.187 | 0.078 | 0.053 | 0.079 | 0.053 | 0.410 | 0.406 | 0.408 | 0.410 | 0.982 | 0.983 | 0.982 | 0.982 | |
O-C | 0.218 | 0.276 | 0.178 | 0.282 | 0.179 | 0.420 | 0.391 | 0.394 | 0.421 | 0.964 | 0.977 | 0.976 | 0.964 | |
CH3OH | O-H | 0.158 | 0.328 | 0.236 | 0.332 | 0.240 | 0.462 | 0.421 | 0.422 | 0.464 | 0.962 | 0.988 | 0.988 | 0.961 |
C-H | 0.182 | 0.077 | 0.053 | 0.078 | 0.053 | 0.413 | 0.409 | 0.411 | 0.413 | 0.983 | 0.984 | 0.984 | 0.983 | |
CH3F | F-C | 0.355 | 0.417 | 0.212 | 0.442 | 0.208 | 0.382 | 0.322 | 0.333 | 0.383 | 0.930 | 0.942 | 0.937 | 0.932 |
C-H | 0.180 | 0.089 | 0.062 | 0.091 | 0.062 | 0.414 | 0.410 | 0.411 | 0.414 | 0.982 | 0.984 | 0.984 | 0.982 |
Molecule | Bond | m0 | m1 | á [^(t)]zm ñ | á [^(t)]zm2 ñ | á [^(t)]+m ñ | |||
(10) | (24) | (10) | (25) | (10) | (25) | ||||
H2 | H-H | 0.000 | 0.000 | 0.000 | 0.000 | 0.437 | 0.437 | 0.992 | 0.992 |
CH4 | C-H | 0.075 | 0.018 | 0.065 | 0.052 | 0.414 | 0.413 | 0.982 | 0.983 |
NH3 | N-H | 0.091 | 0.006 | 0.069 | 0.065 | 0.419 | 0.419 | 0.983 | 0.984 |
H2O | O-H | 0.343 | -0.002 | 0.244 | 0.249 | 0.464 | 0.466 | 0.959 | 0.959 |
HF | F-H | 0.540 | 0.000 | 0.309 | 0.308 | 0.451 | 0.457 | 0.925 | 0.929 |
C2H6 | C-C | 0.000 | 0.000 | 0.000 | 0.000 | 0.405 | 0.405 | 0.982 | 0.982 |
C-H | 0.068 | 0.003 | 0.050 | 0.047 | 0.413 | 0.413 | 0.983 | 0.983 | |
C-C | 0.007 | 0.009 | 0.011 | 0.005 | 0.405 | 0.405 | 0.982 | 0.982 | |
C3H8 | C1-H | 0.068 | 0.006 | 0.052 | 0.047 | 0.413 | 0.413 | 0.983 | 0.983 |
C2-H | 0.061 | -0.010 | 0.035 | 0.042 | 0.412 | 0.413 | 0.984 | 0.983 | |
Cyclo- | C-C | 0.000 | 0.000 | 0.000 | 0.000 | 0.402 | 0.402 | 0.981 | 0.981 |
propane | C-H | 0.092 | -0.007 | 0.059 | 0.064 | 0.414 | 0.415 | 0.983 | 0.982 |
N2H4 | N-N | 0.000 | 0.000 | 0.000 | 0.000 | 0.361 | 0.361 | 0.960 | 0.960 |
N-H | 0.094 | -0.002 | 0.065 | 0.066 | 0.416 | 0.417 | 0.983 | 0.983 | |
N-C | 0.033 | -0.007 | 0.016 | 0.021 | 0.393 | 0.393 | 0.977 | 0.977 | |
CH3NH2 | N-H | 0.102 | -0.005 | 0.069 | 0.073 | 0.419 | 0.419 | 0.983 | 0.983 |
C-H | 0.078 | -0.001 | 0.053 | 0.054 | 0.410 | 0.410 | 0.982 | 0.982 | |
O-C | 0.297 | -0.021 | 0.178 | 0.193 | 0.420 | 0.425 | 0.964 | 0.962 | |
CH3OH | O-H | 0.346 | -0.018 | 0.236 | 0.253 | 0.462 | 0.469 | 0.962 | 0.958 |
C-H | 0.093 | -0.016 | 0.053 | 0.065 | 0.413 | 0.414 | 0.983 | 0.982 | |
CH3F | F-C | 0.431 | -0.015 | 0.212 | 0.215 | 0.382 | 0.387 | 0.930 | 0.932 |
C-H | 0.092 | -0.003 | 0.062 | 0.064 | 0.414 | 0.415 | 0.982 | 0.982 |
r0CH | kCH | kHCH |
Å | mdyn/Å | mdyn/deg |
FA: 1.069 | 8.30 | 0.509 |
TAsymm: 1.078 | 7.77 | |
TApert: 1.096 | 7.17 | |
Standard MM: | ||
[]: 1.113 | []: 4.5 ¸ 4.7 | []: 0.549 |
[]: 1.105 | []: 5.31 | []: 0.508 |
[]: 1.090 | []: 7.90 | []: 0.493 |
1AvH Postdoctoral Fellow, on leave from the Karpov Institute of Physical Chemistry, Moscow, Russia