J . Phys. Chem. 1990, 94, 5496-5498
5496
problems, the well-understood nonrigid bender model of ammonia,” and to use a Hamiltonian of the Hougen-Bunker-Johns (HBJ) type.’* This explains why we always express the nonhindered motion by means of an angular coordinate. Next it is necessary to select the computational method. Basically, our choice is limited to the variational method or the adiabatic approach. If we decide upon a variational calculation, it is best to split the Hamiltonian as H = Ho + HI,where all terms coupling the “large amplitude” and “harmonic” motions are collected in H,.Homay be constructed as a HBJ-type Hamiltonian: R
which gives us the eigenfunctions that are then used as a basis set for the variational calculation. Since in practical applications to polyatomic molecules such variational calculations are hardly feasible for more than four modes, it is normally necessary to look for a less rigorous approach. The “adiabatic” approach seems to be well suited for this purpose. In this approach we also distinguish between “large amplitude” and “harmonic” motions, but in addition we introduce the condition of “adiabaticity”, which means that “harmonic” (“fast”) motions must be considerably higher in energy than the “large amplitude” (“slow”) motions. The idea is of course the same as that of slow nuclei and fast electrons in electronic structure theory, and we make use of this analogy to give the adiabaticity condition a semiquantitative expression. In ammonia, for example, the highest vibrational frequencies are of about 3500 cm-I, compared to the energy of the lowest excited singlet state of about 50000 cm-I. The same ratio of energies is satisfied for a pair of vibrational frequencies of 200 and 2900 cm-I. For some types of vibrational motion the energy separation may be considerably smaller and we may-on the basis of experience accumulated in the literature-still expect a weak interaction. When the adiabatic approximation is applied, the “harmonic” motions may be treated separately, and we may obtain for their vibrational energy an explicit expression U(Sl, ...,Smvrrcl, ...,u,) as a function of coordinates SI,..., S, for large amplitude motions and quantum numbers of the “harmonic” modes, v , + ~ , ..., v,. In the simplest case, U(Sl, ..., S,,,V,+~,..., u,) may be evaluated in the harmonic approximation. If a potential energy fit also contains cubic (and higher) force constants for coordinates (17) Spirko, V. J . Mol. Specrrosc. 1983, 101, 30. (18) Hougen, J. T.; Bunker, P. R.; Johns, J. W. C. J. Mol. Specrrosc. 1970, 34, 136.
SH1, ..., S,, the anharmonicity may be accounted for by per-
turbation theory. The energy of “harmonic” modes for the preselected quantum numbers ulrcI, ..., v, may be then added to a pure potential V(Sl, ..., SJ, and the resulting effective potential may be used for a rigorous treatment of large amplitude motions. We use this potential for the construction of the HBJ-type Hamiltonian, and the respective Schrijdinger equation is solved numerically to give a manifold of vibrational levels for a set of preselected quantum numbers vrrcl,..., v,. The condition of adiabaticity was well satisfied for propargylene. The 0 1 frequency for the large amplitude CCH bending mode was about 100 cm-I, and that for the near harmonic CC stretching was about 1400 cm-I. We considered this difference large enough for a safe use6 of the adiabatic approach. The situation with cyclobutadiene was different. The 0 1 frequencies corresponding to coordinates (4)-(6) were 1630, 1420, and 1220 cm-I. The small differences between these frequencies preclude the use of the adiabatic approach, and we have therefore treated’ cyclobutadiene by a variational calculation. To summarize this section, we list a sequence of steps that should be followed in a theoretical approach to a nonrigid molecule: (i) Find all extrema on the potential surface and calculate analytically for them second and, if possible, higher energy derivatives. (ii) Calculate the energy at additional points, and get a fit for the whole potential surface. (iii) Inspect all the minima, and estimate the number of vibrational levels they can accommodate. Identify the nonhindered large amplitude motions that are most associated with the anharmonicity of the system. (iv) Design a suitable set of coordinates, and express the potential function in these coordinates. Identify in the potential function the terms representing the large amplitude motions, the harmonic potential for the other modes, and the coupling between the kinds of motions. (v) Drop the coupling terms from the potential function and add the corresponding kinetic energy terms to get Ho. Obtain the basis set functions. These are typically numerical solutions for the large amplitude motions and harmonic oscillator eigenfunctions for the other modes. (vi) Perform a variational calculation to mix these harmonic basis functions. If it is not feasible because of too large a number of modes, apply the adiabatic approximation. Registry No. Propargylene, 2008-19-7; cyclobutadiene, 1 120-53-2.
-
-
Nonnuclear Attractors in the Li, Molecule J. Cioslowskif Department of Chemistry and Supercomputer Computations Research institute, Florida State University, Tallahassee, Florida 32306-3006 (Received: November 1, 1989; In Final Form: February 2, 1990) Topological features of the electron density in the Li2 molecule, calculated at the HF/6-311+G* and MP2/6-311+G* levels, are described and rationalized within the catastrophe theory. Depending on the internuclear distance, the density exhibits zero, one, or two nonnuclear maxima. The phase diagrams for the electron density are preented and the bifurcation points are characterized by the respective critical exponents for the differences in the positions of extrema, electron densities, and the electron density Laplacians at the extrema. Introduction Since the electron density of an isolated atom falls off expo. nentially at large distan=,l* one a n easily demonstrate that the electron density of any homonuclear diatomic molecule is expected to possess no maxima other than the cusps at nuclei, provided the ‘Camille and Henry Dreyfus Foundation New Faculty Awardee.
internuclear distance R is sufficiently large. Moreover, it can be proven rigorouslylb that the electron density in the H2+molecular ion has no maxima other than the nuclear CUSPS for all values of (1) (a) Alrichs, R.; Hoffmann-Ostenhof, M.; Hoffmann-Ostenhof, T.; Morgan 111, J. D. Phys. Rev. 1981, ,423,2106. (b) Hoffmann-Ostenhof, T.; Morgan 111, J. D. J . Chem. Phys. 1981, 75, 843.
0022-3654/90/2094-5496$02.50/0 0 1990 American Chemical Society
The Journal of Physical Chemistry, Vol. 94, No. 14, 1990 5497
Nonnuclear Attractors in the Li Molecule
R. In most (homonuclear or heteronuclear) diatomic molecules, this property is observed empirically. However, there are some notable exceptions to this rule, the most known example being the Li2 molecule in which the presence of the midpoint maximum in the HartretFock electron density was recognized in as early as 1 9 S 2 Subsequently, it has been conjectured3 that inclusion of electron correlation would result in disappearance of the nonnuclear maximum. However, the recent studies demonstrated that the maximum persists even if the electron density is calculated within the MRD CI approach.' Similar observations are pertinent to the Na2 mole~ule.~ As pointed out by the authors of ref 4, the nonnuclear maximum is converted to an ordinary minimum when the internuclear distance is decreased to R = 4 au. In the language of the theory of atoms in molecules3 this means that the corresponding nonnuclear attractor disappears as R becomes smaller than some critical value. At the critical value of R, the molecular structure undergoes a bifurcation catastrophe.6 We believe that the unusual features of the electron density in the Li2 molecule deserve closer investigation. In the present paper, the results of calculations on both the Hartree-Fock (HF) and correlated electron densities, carried out for various values of R, are reported. The bifurcation points are located and classified by using the formalism of the catastrophe theory. Overview of the Catastrophe Theory The catastrophe theory, popularized by Thom,' has gained a widespread application in describing natural processes. It uses the concept of potentials V(c,z),parametrized by the control spuce C, on the behaoior space I. In the particular case of I being one-dimensional, the ensuing catastrophes are called cuspoids* and can be classified by considering the generic potentials of the form V(c,z) =
A(z
-
ZO)"+'
+ u,(c)(z - 20)" + ... + u~(c)(z- 10) + B (1)
In the above equation, A and B are constants, whereas the coefficients ulr ...,a, depend on c. The beauty of the catastrophe theory stems from the fact that the same formalism can be applied to very diverse phenomena, such as biological processes,oscillatory chemical reactions, and molecular structures, to name a few. There are only a few different kinds of elementary catastrophes. In particular, in this paper we deal with the fold and cusp catastrophes. The use of catastrophe theory to describe the electron density in molecules has been pioneered by Collard and Hall? In the description of topological features of the electron density in diatomic molecules, the potential Vcomsponds to the electron density along the chemical bond, and z is the distance from the center of mass; and the control space C is one-dimensional, the control variable c being the internuclear distance R. Two mast elementary catastrophes correspond to n = 3 and n = 4. The fold catastrophe occurs when the electron density behaves locally as p(R,z) = A ( z - 20)3
+ al(R)(z - zo) + B
(2)
with A > 0. For a l ( R )> 0, the density possesses no minima or maxima, whereas for a l ( R ) C 0 there is a pair of one minimum and one maximum at (2) Besnainou, S.; Roux, M.; Daudel, R. Compt. Rend. 1955, 241, 31 1. ( 3 ) Bader, R. F. W.; Nguyen-Dang, T. T.; Tal, Y. Rep. Prog. Phys. 1981, 44, 893. ( 4 ) Gatti, C.; Fantucci, P.; Pacchioni, G. Theof. Chim. Acta 1987,72.433. ( 5 ) Cao, W. L.; Gatti, C.; MacDougall. P. J.; Bader, R. F. W. Chem. Phys. Lett. 1987, 141, 380. (6) Bader, R. F. W.; Tal, Y.; Anderson, S. G.; Nguyen-Dang, T. T. Isr. J . Chem. 1980. 1 9 , 8 . ( 7 ) Thom, R. Stabilitl Structurelle et Morphog8nEse; Benjamin-Addison-Wdey: New York, 1973. (8) Woodcock, A. E. R.; Poston, T. A Geometrical Study of the Elementary Catastrophes. In Lecture Notes in Mathematics; Springer-Verlag: New York, 1974; Vol. 373. ( 9 ) Collard, K.; Hall, G. G. In?. J . Quant. Chem. 1977, 12, 623.
z = 20 f [-a1(R)/3A]"2
(3)
The catastrophe is triggered at the bifurcation point, R = %, for which a l ( R ) = 0. Since, in the vicinity of &, al(R) depends linearly on IR - 4,this results in the following critical exponents: tz= 1/2, for the difference between the positions of the extrema of p: k,(R - Ro('l2
AZ = 2[-~l(R)/3A]'/~
(4)
{, = 3/2, for the difference between electron densities at the extrema:
Ap = (4/3)[-~1'(R)/3A]'/~
N
k,lR - Ro13/2
(5)
and {' = 1/2, for the difference between the values of the electron density Laplacian, K = V2p,at the extrema: AK N k,lR
- Rol'12
(6)
The cusp catastrophe occurs when the electron density behaves locally as p(R,z) = A(z
- z0)4 + a2(R)(z- zo)2+al(R)(z - zo) + B (7)
with A > 0. For homonuclear molecules zo = 0 and a l ( R ) 0, because of symmetry. In this case eq 7 becomes
+
+
~ ( R , zP;): A Z ~ a 2 ( ~ ) z 2 B
(8) The catastrophe is triggered at R = R,,for which a2(R)= 0. For a2(R)> 0 the density possesses a minimum at z = 0, whereas for a2(R)C 0 there is a maximum at z = 0 and a pair of minima at
z = f[-a2(R)/2A]'/2
(9)
Because
AZ = [-a2(R)/2A]'/2 = k,lR - R o y Ap = -Uz2(R)/4A AK
k,)R - R0I2
k,)R - Rol
(10)
(1 1) (12)
the critical exponents for the bifurcation point associated with a cusp catastrophe are tz= 112, {, = 2, and t, = 1. A further classification of the bifurcation points is made possible by taking into account the sign of the derivative D = -(a/aR)%2(R)IR-&
(13)
where n = 3 for the fold catastrophe and n = 4 for the cusp catastrophe. For D > 0 one deals with either a (+)-fold or a (+)-cusp bifurcation point. The opposite sign of D corresponds to either a (-)-fold or a (-)-cusp bifurcation point. One easily realizes that the bifurcation points of the (+)-type generate two additional extrema as R increases and passes through &. The bifurcation points of the (-)-type result in coalescence of two extrema at R = R,, and therefore decrease the number of extrema by two for R > Ro.As a corollary of this observation one finds out that if a molecule possesses no nonnuclear maxima at small values of R, then the total number of the (+) and (-) bifurcation points, counted over the entire range of internuclear distances, must be equal.
Results and Discussion The calculations on the Liz molecule have been carried out within the 6-31 1+G* basis set.I0 Both the H F density and the MP2 relaxed density" have been computed. In the latter a p proach, the electron density is calculated as the energy derivative with respect to a perturbation being Dirac's delta function. At the HF level, the equilibrium Li-Li distance is found at R = 5.262 au, and the corresponding energy is -14.870 352 5 au. At the MP2 (IO) Krishnan, R.; Frisch, M. J.; Pople, J. A. J. Chem. Phys. 1980, 72, 4244. Clark, T.; Chandrasekhar, J.; Spitznagel, G. W.; Schleyer, P. v. R. 1. Comput. Chem. 1983.4, 294. (1 1) Salter, E. A.; Trucks, G. W.; Fitzgerald, G.; Bartlett, R. J. Chem. Phys. Lett. 1987, 141, 61. Trucks, G. W.; Salter, E. A.; Sosa, C.; Bartlett, R. J. Chem. Phys. Lett. 1988, 147, 359.
5498
The Journal of Physical Chemistry, Vol. 94, No. 14, 1990
4.
Cioslowski TABLE I: Bifurcation Points in tbe Liz Molecule (HF Electron Density)
I
(+)-cusp bifurcn pt RO Zn
5; k. {P
kP
i-‘
k,
A
B
(-)-fold C , , Ct
3.898 0.0 1 12 f 6 . 0 3 X 10-1 2 -3.33 x 10-3 1 9.27 X
6.189 0.0 1 12 f1.19 2 4.63 x io-‘ 1 -4.09 X IO-’
6.853 f1.193 112 f1.27 x 10-1 312 2.51 x 104 112 -5.06 X 10”
TABLE II: Bifurcation Points in the Li2 Molecule (MP2 Electron Density)
(+)-cusp I.
I
I
I
2.
3.
4.
J
5.
Figure 1. The phase diagram for the Liz molecule (the HF density). See text for further explanation.
bifurcn ut RO Zn
r;
k.
i-; k,
.rc
k,
I.
I.
2.
3.
4.
5.
Figure 2. The phase diagram for the Li, molecule (the MP2 relaxed
density). See text for further explanation.
level, the energy minimum of -14.915 5294 au is attained at R = 5.172 au. The dependence of the topological features of the electron density in a diatomic molecule, XY, on the distance R is easily visualized with the help of a phase diagram. The positions of the density extrema are plotted by using the distances to the nuclei X and Y as the x and y variables, respectively. The minima are denoted by thick solid lines and the maxima by thick broken lines. The topological features of the electron density a t the internuclear distance R can be found by drawing a straight line given by the equation y = R - x. The phase diagram of a “normal” homonuclear molecule consists of a straight solid line given by x = y . The phase diagram of a normal heteronuclear molecule has a curved solid line. The phase diagrams of the Liz molecule are displayed in Figures 1 and 2. For both the HF and MP2 electron densities, there are two (+)-cusp bifurcation points, A and B, and a pair of the (-)-fold ones, C1and C1. The domain between the solid line CIACzand the broken line C1BC2is a basin of a nonnuclear attractor arising from the midpoint maximum in the electron density. For larger values of R, this attractor is split into a pair of nonnuclear attractors with their basins located inside boundaries given by the points BCID and BC2D, respectively. The remaining areas are basins of the nuclear attractors.
A
B
(-)-fold Ci, C,
f6.01 X IO-1 2 -3.22 X 1 9.13 X IO-*
6.038 0.0 1/2 f1.23 2 6.26 X IO-’ 1 -4.98 X lo-’
6.510 f 1.034 112 f8.17 X IO-’ 3f2 3.63 X lo-’ 112 -5.88 X lo-’
3.904 0.0 If2
Within the 6-31 1+G* basis set, the electron density possesses zero, one, or two nonnuclear maxima. The positions of the bifurcation points and the corresponding values of k and t are presented in Tables I and 11. Inclusion of electron correlation affects the position of the bifurcation point A only slightly. The onset of both the (+)-cusp (the point B) and the (-)-fold (the points C I and C,) catastrophes is shifted toward shorter internuclear distances by the electron correlation effects. Therefore, the domain of existence of the nonnuclear attractors in the Liz molecule contracts upon inclusion of electron correlation. One should comment on the possible basis set dependence of the present results. Test calculations show that adding two more d functions to the basis set on each lithium atom does not alter the shape of the phase diagram. However, since the midpoint electron densities are quite small, it is possible that for even larger basis sets the curve ClBC2 would become convex rather than concave. This would result in disappearance of the bifurcation points CI and C2. As a consequence, the point B would become a (-)-cusp bifurcation point. The electron density would exhibit zero or one nonnuclear maxima depending on R. Conclusions
The phase diagrams aid in visualizing the dependence of the topological features of the electron densities in diatomic molecules on the internuclear distance. In case of the molecules with nonnuclear attractors, such as the Liz system, the bifurcation points can be easily classified with the help of the catastrophe theory. The critical exponents provide further insight into the properties of the bifurcation points. Finally, the influence of electron correlation on the topological features of the electron density can be studied by employing the HF, CI, or relaxed MBPT of CC densities. Acknowledgment. This research was partially supported by the Camille and Henry Dreyfus Foundation New Faculty Award Program and by the US.D.O.E. under the contract DE-FCO585ER250000. The calculations were performed on the ETA-IOG supercomputer located at the Florida State University. The author thanks Prof. R. J. Bartlett and Dr. J. Watts (University of Florida, Gainesville) for their help with the MP2 relaxed densities program and Dr. J. P. Ritchie (Los Alamos National Laboratory) for interesting discussions. Registry No. Li,, 14452-59-6.
