5422
J . Phys. Chem. 1994, 98, 5422-5430
Theoretical Studies on the Relative Stability of Neutral and Protonated N,N'-Diarylguanidines in Aqueous Solution Using Continuum Solvent Models Giuliano Alagona' and Caterina Ghio CNR-Istituto di Chimica Quantistica ed Energetica Molecolare. Via Risorgimento 35, 1-561 26 Pisa, Italy
Peter I. Nagyo*+and Graham J. Durant Center for Drug Design and Development, Department of Medicinal and Biological Chemistry, The University of Toledo, Toledo, Ohio 43606-3390, and Cambridge Neuroscience, Cambridge, Massachusetts 021 39 Received: December 16, 1993; In Final Form: March 18, 1994'
The conformational properties in vacuo and in solution of neutral and protonated N,N'-diphenylguanidines have been studied making use of ab initio SCF and MP2 calculations in vacuo and continuum solvent free energy calculations in solution on the STO-3G and 4-31G optimized geometries obtained in vacuo. For the N,N'diarylguanidine series with condensed rings (di-l-naphthyl, di- l-anthracenyl, and di-9-anthracenyl derivatives) the calculations have been carried out on the relevant partial charges located at the molecular mechanical geometries derived from the AMBER force field in MacroModel. The STO-3G solvation free energy is considerably less favorable than the 4-3 1G one for the diphenylguanidines, but only about 2 kcal/mol less favorable for the diphenylguanidinium rotamers, which on the whole show a better agreement between the two basis sets. The STO-3G basis set predicts the syn-anri (SA) structure to be the most stable in solution for neutral diphenylguanidines, whereas the 4-3 1G basis set favors the A S rotamer. Extended basis set calculations on N-phenylguanidine favor the anti structure as well. The MP2/4-3 1G//4-3 1G correlation corrections applied in solution produce a solvent effect lower than at the SCF level, and only the solvation free energies of AA and AA+ are reduced with respect to that of AS. The cavitation free energy is nearly independent of the basis set. The solvation free energy obtained from the partial charge description of the diphenyl derivatives is analogous to that produced by the STO-3G basis set for the neutral rotamers, while it is about 10 kcal/mol more favorable than the corresponding STO-3G values for the protonated conformers. The free energy a b initio and electrostatic results are in fair qualitative agreement with those obtained with the GB/SA method. The stacking of the aromatic rings is followed by unfavorable electrostatic and favorable van der Waals interactions. The hydration of these rings is slightly favored by the GB/SA method in a nonstacked and rather separated form.
Introduction In the real world most interactions and reactions occur in a condensed phase and, mainly, in aqueous solution. This accounts for the wide effort devoted with a variety of involving either explicit water molecules or a continuum, to the study of water as a solvent. The possibility of computingreliable values for the free energies in solution making use of the reaction field theory on the ab initio description of the solute3 allows the use of limited computer resources in comparison to statistical methods. The results obtained have then been compared with the free energy differences evaluated2c.d by using Monte Carlo free energy perturbation methods,2a and a fine agreement was found: even though the continuum models are intrinsically unable to account for specific solute-solvent H-bond interactions. The solvent reaction field, however, can simulate these specific interactions: the "apparent" charge distribution generated on the cavity surface may well represent either H-bond donor or acceptor water molecules. Thus, the continuum model can follow conformation dependent solutesolvent interactions. The substitution of ab initio charge distributions with partial charge models, resorting to a completely classical approximation in place of a semiclassical one, does not greatly alter the picture of the phenomenon, which remains well-defined even when small changes (such as those depending on conformational degrees of freedom) areinvolved.5 The problems linked to the determination
of the geometries may be solved with the help of molecular mechanical computer codes, whose wide availability encourages their use for virtually every molecular system. On the other hand, their supposed or alleged flexibilityrepresents a source of possible misuse. It is therefore advisable to perform suitable checks between the results obtained with different methods. Physicochemical studies of guanidine derivatives are of particular interest since guanidine is a building block for many compounds of biological and therapeutic interest,68including the amino acid arginine. More recently, diarylguanidines have been reported to be highly active as antagonists at 0 receptor site@ and at the ion-channel associated with the N-methyl-D-aspartic acid (NMDA) receptor.6f As a strong base, guanidine and its derivatives are generally protonated in aqueous solution at physiological pH. Thus, theoretical calculations supporting drug design have studied the ion-pairing of positive ions in water? or conformational analyses have been carried out for several rotamers of the protonated di- and trisubstituted derivatives.8 Due to the large number of atoms in the substituted molecules, the conformational analyses have been performed only at the molecular mechanics level without considering solvent effects.8 Though the large atom number is a limiting factor even in the present study, we report here the results of relative energylfree energy calculations based on ab initio calculationsand compare the results with the combined molecular mechanics/continuum solvent approach implemented in the MacroModel s ~ f t w a r e . ~
Brief Outline of the Methods ?On leave from Chemical Works of Gedeon Richter Ltd., Budapest, Hungary. Abstract published in Advance ACS Absrrucrs, April 15, 1994.
A thorough description of the procedure followed to account for a quantum mechanical solute M immersed in a continuous
0022-3654/94/209~-5422304.50/0 0 1994 American Chemical Society
Relative Stability of N,N'-Diarylguanidines
The Journal of Physical Chemistry, Vol. 98, No. 21, 1994 5423
dielectric solvent (of dielectric constant e) can be found in the sourcepapersof the method.3 Sufficient hereis a short explanation of the computed quantities, which derive from the Hamiltonian of the system consisting of two terms, one related to the unperturbed solute and the other (Vu) to the solute-solvent interaction:
The solute, placed in a cavity inside the dielectric medium, induces an apparent charge distribution u on the cavity surface, depending on the solute total charge distribution, To,and on the apparent charge distribution itself, uo. The trial value for uo (including the mutual polarization of the surface charges and the extracharge compensation) is used to polarize the charge distribution of M, obtaining I" and hence 6'. The next step, leading to I"' and u", is usually sufficient to reach convergence. The main energetic quantity obtainable is the free energy, G:
GSCF= E,,, - l/ZJI'M(r) VJr) dr (with E,,, = E M by
+ Enuf).3f The solvent effect is therefore given sol = &CF
- EO
(3)
where Eo is the energy of the solute in vacuo. When we are dealing with a rigid charge distribution of the solute, in place of eq 3 only
need to be considered in the case of either a continuum or a discrete charge distribution, respectively. The cavitation term is computed separately, with two different methods (see below). The molecular mechanics/continuum solvent method as implemented in the MacroModel program9 is a combined method to calculate energies/free energies in solution. Internal energies, after geometry optimization, are calculated via molecular mechanics. The solvent effect is estimated with the generalized Born/surface area (GB/SA) method.1° In this semianalytical approximation the solvation free energy is the sum of the solutesolvent electrostatic polarization (GP1), dispersion (Gvdw), and cavity formation (Gav) terms. The latter two terms are considered to be proportional to the solvent accessible surface area. The solvent is characterized by its dielectric constant, and the polarization is taken into consideration with a modified Born formulalo using predefined and constant solute charges. The advantage of the method is providing optimized molecular geometries relevant to the solvent environment; thus, no gasphase optimized geometries are needed in solution calculations, as was applied in combined ab initio + Monte Carlo simulations in previous studies.2ad
Computational Details The systems considered are both neutral, HN=C(NH&, and protonated, C+(NH2)3, guanidines with two H's in position N and N', respectively, substituted by two phenyl rings in the smallest compounds. The substituents considered are also made up of either two (naphthyl) or three (anthracenyl) condensed aromatic rings. Depending on the orientation of the substituents with respect to the C=N double bond, the diphenyl neutral compound can be found in the anti-anti (AA), anti-syn (AS), syn-anti (SA), or syn-syn (SS) conformations (see schemes), while the protonated compound shows only three rotamers due to the NH2 group symmetry, namely, SS+, AA+, and AS+ = SA+. The same kinds of rotamers are present in the di-9-anthracenyl case, while
SCHEME 1 N
thedi- 1-naphthyl substitution, analogously to thedi-1-anthracenyl one, produces also the AA' and AA'+ rotamers (03 and 0 4 are either both positive or of opposite sign). The other possible rotamers (when both 0 3 and 0 4 are negative) have not been taken into account, because the overlap of the two aromatic rings is small, leading to very extended structures, presumably unfavorable in aqueous solution. On the smallest compounds (diphenylguanidines) geometry optimizations a.t the SCF level have been carried out using both the STO-3GlI and 4-31GI2 basis sets with either the Gaussian90 or -92 series of programs.13 Ab initio calculations in solution3 at the S C F level have been carried out employing MGPIPC,14 our version of MONSTERGAUSS15that is capable of dealing also with the point charge models of very large molecules, while the MP216 correlation corrections, zero-point energy (ZPE), and thermal corrections have been computed again with Gaussian90 and Gaussian92. The reduction of the rotational entropy upon C2 symmetry and the entropy of mixing for the structures with optical antipodes have been included in the calculation of the free energy. The MP2 calculations in solution have been carried out using the Pisa version17 of HOND08.I8 These codes run on the SGI 4D-420-GTXB and IBM RISC 6000/560-580 workstations at ICQEM. The molecular mechanics calculations have been carried out in Toledo with the MacroModel3.5a software running on a Personal Iris 4D/35 workstation. The AMBER force fieldkg was used with the uall-atomn modelIgb and atomic charges as implemented in MacroModel. The solvation calculation parameters were also taken from the program library. The dielectric constant for the electrostatic interactions of the solute atoms was set to unity, while the continuum water-solute interaction was calculated considering the bulk dielectric constant for water.
Results and Discussion N,N'-Diphenylguanidines. Ab Initio Results. In Scheme 1, the structure of N,N'-diphenylguanidine (AA rigid geometry) is displayed, together with the indication of the four torsional degrees of freedom, all set to 180°: 01 = Crl-Nl-Cc=N and 0 2 = Cr2N2-Cc=N refer to the name arrangement of the conformers, in the order, while 0 3 = Cdngl-Crl-N1-Cc and 0 4 = Cring2-CrZN2-Cc give the deviation of the aromatic rings from the N-C,-N plane. The geometries of the different possible conformers of N,N'-diphenylguanidines have been optimized at the SCF/STO3G and 4-3 1G levels in vacuo. The optimized dihedral angles are compared in Table 1, both for the neutral and protonated diphenylguanidines. The relative energies of the rotamers with respect to the AS conformers, taken as 0 for both the neutral and protonated compounds, are reported in Table 2, under the heading A&F. The AS and SS+ rotamers are the most stable at the 4-31G level, while the SS and AS+ rotamers are the most stable a t the STO-3G level, even though the largest energy difference between the two descriptions of the couples AS/SS and AS+/ SS+ is only about 0.4 and 0.08 kcal/mol, respectively. An outstanding difference between the STO-3G and 4-3 1G descrip tions of the neutral compounds can be found in the syn
5424 The Journal of Physical Chemistry, Vol. 98, No. 21, 1994 TABLE 1: Torsion Angles (des) Optimized in Vacuo at the STO-36 and 4-316 Levels for the Neutral and Protonated N,N'-Diphenylguanidines rotamer a2 a3 a4 PI' mu AS
ss
SA AA AS+
ss+
AA+
AS
ss
SA AA AS+
ss+
AA+ 91 is
120.5 10.1 4.8 121.7 178.4 0.0 159.3
6.9 4.3 -177.5 125.9 -1.0 0.0 159.3
STO-3G 18.7 -49.0 -50.7 21.6 -67.0 -66.7 -46.1
-60.6 -60.4 -53.1 18.0 -64.8 66.7 -46.1
-8.9 143.4 138.6 -8.0 -4.4 177.9 -8.1
140.6 135.6 -34.8 -9.5 175.2 -177.9 -8.1
162.1 10.4 7.0 145.1 180.0 0.0 156.9
-14.3 18.7 -162.9 141.8 0.0 0.0 156.9
4-31G -41.0 -4.2 -4.0 -24.7 -90.0 -90.0 -67.5
-39.6 39.9 50.6 -22.6 -90.0 90.0 -67.5
-9.2 -163.1 -171.1 -14.0 0.0 180.0 -16.0
160.6 -149.0 16.2 -25.5 180.0 180.0 -16.0
the dihedral angle H-NrC,=N
(see Scheme 1).
TABLE 2 Relative Energies (kcaVmo1) in Vacuo and in Solution at the SCF Level (STO-36 and 4-316 Basis Sets) for the Neutral and Protonated N,N'-Diphenylmnidines rotamer
AS'
ss
SA AA ASt
ss+ AA+
ASc
ss
SA AA ASt
d
ss+ AA+
MSCF
AGSCF
0.0 -0.190 2.349 1.849 0.0 0.064 5.525
STO-3G 0.0 -0.562 1.400 2.112 0.0 -1.074 2.240
-6.255 -6.627 -7.204 -5.992 -44.0 19 -45.157 -47.303
0.0 0.424 1.098 4.745 0.0 -0,083 9.942
4-31G 0.0 2.91 1 3.031 5.648 0.0 -0.51 1 6.853
-14.061 -1 1.574 -12.127 -13.158 -46.634 -47.062 -49.724
&F
A G h F 0.0 -0.372 -0.949 0.263
0.0 -1.138 -3.284
0.0 2.488 1.934 0.903 0.0 -0,428 -3.089
E ~ C F= -654.990 110 au; G ~ C F= -655.000 078 au. ESCF = -655.483 238 au; G ~ C = F -655.553 386 au. ESCF= -662.247 558 au; GKF = -662.269 966 au. E ~ C=F-662.682 035 au; G ~ = F-662.756 35 1 au. (I
arrangements of the first ring (that located at the left-hand side in Scheme l), namely, the SS and SA structures: the STO-3G basis set favors gauche forms, whereas the 4-3 1G basis set favors a planar arrangement, with the ring H facing the N lone pair. The different behavior of the two basis sets is reflected by the free energy values in solution ( A G ~ ~ Falso ) , reported in Table 2.
."
AS
SS
SA Conformers
AA
Alagona et al. The solvent effect at the SCF level, is much stronger, as expected, for the N,N'-diphenylguanidinium rotamers, due to the presence of a charged species in solution. Their solvent effect values, however, differ only by about 2 kcal/mol from the STO3G ones, despite a larger deviation from the NC,N plane of the phenyl rings found at the 4-31G level. The basis set effect is more evident for the neutral compounds, where the 4-3 1G solvent effect is about twice as large as the STO-3G one, with, in addition, inversions in the stability order. To assess whether thedifference dependson geometric factors or on theoverestimate of thecharges produced by the split valence shell basis sets (the charge on the H facing the N lone pair is 4 times as large as that found at the STO-3G level for SA), we repeated the calculations making use of the 4-31G basis set at the STO-3G geometries and vice versa. The results are compared in Figure la. The trend clearly depends on the geometry, whereas the range of the values depends on the basis set (and on the level). It is possible in fact to distinguish two families, one related to the 4-3 1G results (below -1 1 kcal/ mol) and the other to the STO-3G results (above -8 kcal/mol). In between there is the solvent effect at the MP2/4-3 1G level on the 4-31G geometry (MP2/4-31G//4-31G), which is located almost exactly halfway between the SCF/4-3 1G//4-3 1G results and the SCF/STO-3G//4-31G ones. The 4-31G geometry produces a lower stabilization upon solvation of the SS and SA structures with respect to the AS one, because the N lone pair is hardly available to interact with the solvent (see Table 1: =4 'for SA at the 4-31G level, whereas = -51O at the STO-3G level). The basis set dependence of the ab initio results for the protonated rotamers is very limited, as can be seen on examining Figure lb. In Tables 3 and 4 the various energetic contributions (correlation corrections at the MP2 level, ZPE and thermal corrections at 298 K and 1 atm) to the internal free energy in the gas phase are reported for the neutral and protonated diphenyl-substituted compounds, respectively. The MP2 correlation corrections (computed on both basis sets for analogy, even though it is not advisable to perform such calculations on a minimal basis set) increase the energy gap between AS and SS,or SA, by about 1k cal/mol with either basis set, while the energy difference between AS and AAdecreases by about 0.8 and 2.8 kcal/mol at the STO3G and 4-3 1G level, respectively. An analogous trend is found for the cations, but the 4-31G basis set effect is almost 2.5 times as much as the STO-3G one for both SS+ and AA+. With the addition of the other terms (Le., ZPE, G = H - TS) the stability order does not change. The inclusion of the solvation free energy slightly favors AA with respect to SS at both the SCF level and MP2 level using the 4-3 1G basis set. SS+and AAt are stabilized with respect to AS+ using the STO-3G basis set, whereas the large solvent effect (especially on the AA+ rotamer) is not
-58
'
1
I
\..
AS+
ss+
AA+
Conformers Figure 1. Solvent effect, @I, on the (a) neutral and (b) protonated diphenylguanidine rotamers, as described at the SCF/STO-3G level (dashed line) and at the SCF/4-31G (solid line) and MP2/4-31G (dash-dotted line) levels. Partial charge descriptions of the solute are indicated by markers (triangles for the STO-3G geometry, squares for the 4-31G one, circles for the MacroModel dihedral angles joined to the STO-3G internal geometry). STO3G//4-31G (long dashed line) and 4-31G//STO-3G calculations (long dash-dotted line) are also reported for the neutral compounds.
Relative Stability of N,N‘-Diarylguanidines
The Journal of Physical Chemistry, Vol. 98, No. 21, 1994 5425
TABLE 3 Energy Contributions (kcal/mol) of the Various Conformers of the Neutral N,N’-Diphenylguanidine with Respect to the AS Conformer’ (Taken as 0) STO-3G SS
SA
4-31G
SS
AA
SA
a Absolute values for AS: MP2 = -655.892 284 au, ZPE = 0.273 554 au, H = 0.285 563 au, TS = 34.0487 kcal/mol at the STO-3G level and MP2 = -663.691 165 au, ZPE = 0.258 557 au, H = 0.270 662 au, TS = 34.0964 kcal/mol, G ~ p = z -663.708 508 au, GELz = -10.8829 kcal/ mol at the 4-31G level.
TABLE 4 Energy Contributions (kcal/mol) of the Conformers of the Protonated N,N’-Diphenylguanidine with Respect to the AS+ Conformer’ (Taken as 0) STO-3G
-TS
+
MP2 ZPE AGhF
+ H - TS
AG:;,
SS+
AA+
SS+
AA+
0.232 0.014 0.012 0.423 0.681 -1.138 -0.457
3.636 -0.192 -0.095 -0.606 2.743 -3.284 -0.541
0.336 0.010 -0,009 0.419 0.756 -0.428 0.328 0.236 0.992
5.322 0.055 0.129 -0.928 4.578 -3.089 1.489 -2.1 17 2.461
AGL2
Absolutevaluesfor AS+: MP2 =-656.378 622au, ZPE = 0.285 336 au, H = 0.298 197 au, TS = 35.55 kcal/mol at the STO-3G level and MP2 = -664.113 815 au, ZPE = 0.273 541 au, H = 0.285 818 au, TS = 34.29 kcal/mol, G M P=~ -664.186 593 au, G& = -46.296 kcal/mol at the 4-31G level.
sufficient to reverse the stability order with the 4-31G basis set. Interestingly enough, the MP2 correlation corrections exert a noticeable stabilizing effect on the AA and AA+ structures both in vacuo and in solution. Electrostatic Results. The AMBER force field makes use of potential-derived charges for the residues, but the systems under scrutiny, although similar to arginine, are nonstandard and required some adjustment from MacroModel. Therefore, a partial charge model, SK,I9C derived from the best fit to the STO-3G electrostatic potentialI3 of the diphenyl systems has been used for comparison. The AGSo1values obtained, however, although
protonated compounds
-0.150 0.150 0.147 -0.5415 0.388 -0.343 0.248
-0,150 0.150 0.147 -0.442 0.799 -0.657 0.362
C,= Clingconnected to N N(H) connected to C, C, (central carbon) =N connected to C, H connected to each N
TABLE 6 Torsion Angles (deg) from MacroModel for the Neutral and Protonated N.N’-DiDhenvleuanides rotamer AS
ss
SA AA AS+
ss+ AA+
*I
173.9 48.9 4.3 153.7 178.9 -13.3 146.5
*2
*3
*4
111’
12’
-1.2 -50.6 172.6 150.4 4x2 13.3 146.5
-55.5 -51.6 -50.6 -35.8 -68.4 -65.9 -47.4
-54.5 -55.3 -54.4 -33.6 -68.3 65.9 -47.4
-17.0 -143.9 171.3 -16.4 -12.0 178.1 -28.6
171.1 121.9 -18.1 -21.1 169.9 -178.1 -28.6
is the dihedral angle H-N,-C,=N
(see Scheme 1).
showing a trend very similar to the STO-3G results, were slightly farther apart from the ab initio ones than those computed using the MacroModel charges, which thus have been used throughout. In Figure 1 are also displayed the results obtained employing the partial charge description (reported in Table 5 ) of the solute on the STO-3G, 4-31G, and MacroModel (the dihedral angles from MacroModel have been combined with the STO-3G bond lengths and angles) geometries (the relevant torsional angles are reported in Table 6 for comparison). The point charge model results turn out to be very close to each other and to the STO-3G ones, despite the considerable geometry differences. The charge distribution therefore seems to produce a sort of balance with its electrostatic effect. This is important in view of the subsequent studies on the other terms of the series (dinaphthyl and dianthracenylderivatives) for which abinitiooptimized geometries are not available. As far as the protonated rotamers are concerned, the values obtained employing the partial charge description of the solute on the STO-3G, 4-31G, and MacroModel geometries, on the contrary, lie well below the a b initio ones: the solvent effect is 7-1 1 kcal/mol more favorable, with a preference for AA+. The larger solvent effect can be due to the fact that the additional proton charge is localized only on the (-NH)*-C,=NHz+ group. The SK model, in fact, which allows some charge to be delivered to C, (the ring carbons connected to the nitrogens), shows a solvent effect close in value to the a b initio one, but the relevant AGSO1 results are slightly worse, as already stated.
29
30 I
28
29
1 27
28
. -Y
neutral compounds
Cring Hring
a vi
4-31G
AGEp2
atom
AA
MP2 0.599 3.349 1.008 1.256 2.044 1.900 ZPE 0.185 0.500 -0.002 -0.077 0.035 0.175 H 0.150 0.366 -0.010 -0.079 -0.002 0.175 -TS 0.065 0.385 0.219 -0.250 0.212 -0.303 MP2 + ZPE + H - TS 0.999 4.600 1.215 0.850 2.289 1.947 -0.372 -0.949 0.263 2.488 1.934 0.903 AGAF 0.687 3.651 1.478 3.338 4.223 2.850 AGFCF 1.872 1.397 0.729 AGPZ 2.722 3.686 2.676 A G ~ ~
MP2 ZPE H
TABLE 5 Partial Charge Models (au) Taken from MacroModel for the Neutral and Protonated N,N’-Diphenylguanidmes
0
E
26
YZ 25 b
24
I
I
AS
SS
1
1
I
1
SA AA AS+ SSC AA+ Conformers Conformers Figure 2. Cavitation free energy, Gav, produced by the Sinanoglu’s (S, plain lines) and Pierotti’s (P, marked lines) formulas on the (a) neutral and (b) protonated diphenylguanidine rotamers, corresponding to the STO-3G (dashed line), 4-31G (solid line), and MacroModel (short dashed line) geometries of the solute.
5426
TABLE 7: Comparison of the Free Energy Results (kcal/ mol) Obtained from the ab Initio Calculahons and from MacroModel for the Neutral and Protonated N,N'-Diphen ylguanidines rotamer
..
The Journal of Physical Chemistry, Vol. 98, No. 21, 1994
AGSCF+ AGav AGSCF+ AGav STO-3G
4-31G
0.0 -0.61 1.21 1.84 0.0 -0.52 1.81
0.0 2.90 3.16 5.26 0.0 -0.16 5.21
Alagona et al.
SCHEME 2
A
N/H
n\,
Eint + (Gav + GvdW
+ Gpod AMBERIIGBISA
~~~~~~
AS
ss
SA AA AS+
ss+ AAt
0.0 1.93 0.99 4.26 0.0 -0.64 5.96
A tentative dissection of the solvent effect into hydrophobic and hydrophilic contributions produces almsot constant hydrophobic Gsol values (on the order of -3 and -14 kcal/mol for the neutral and protonated compounds, respectively, using the STO3G geometry). Thus, the hydrophilic G W 1 values show the same trend as the total values: the ring positions in the various conformers do not affect the ring interaction with solvent, but mainly change the exposed surface of the guanidine moiety. This suggests that the stacking effect (for a more detailed discussion see the last section) is not very important for thediphenyl derivative even though the Cr4 atoms of the two rings are closer than about 6 A, which is required in the present approximation to accommodate solvent in between the aromatic rings. In fact, the +p3 and torsional angles in the STO-3G geometries for AA and AA+ (Table 1) are far from being 90°,a requirement for a faceto-face stacked geometry. Furthermore, the Crl-N1-Cc-N2-Cr2 atoms are not in a common plane, allowing a further departure from the face-to-face arrangement of the rings. As a result, the conformational changes to AS, SS, etc., do not result in a considerable difference in the hydration of the aromatic rings, as reflected by the nearly constant hydrophobic Gm1 values. Cavitation Term. The free energy of cavitation has been considered using both Pierotti'szo and Sinanoglu'sZLformulas, which give analogous results, just shifted by about 3 kcal/mol for both the neutral and protonated compounds (Figure 2). Both formulas depend on a few solvent parameters taken from experimental data and on the cavity surface, the same as that used to evaluate V,. The cavitation free energy is usually almost independent of conformational changes.22 In this case, however, AA (especially in the MacroModel geometry) and AA+ represent an exception to this rule, because the anti-anti arrangement produces a cavity surface smaller than that surrounding the extended conformations. The comparison (reported in Table 7) of the ab initio results with those obtained from the AMBER//GB/SA calculations is performed using the AGSCF AGav term, which is the most consistent with the sum of G,,,, GVdW, G,I, and internal energy
+
-6-
- (a) -8- -__ = - 0
Ea
Y
-P
Q
-10
-
...,.
- - _- - - _ S T W 441G D85 8-310' 6-31 tG'
-14
-
-16 -
-18
-5
- (b) -7 - - - phenyl(3) - - phenyl(4) -vinyl(3) 9- _ _ _ _vinyt(4) diphenylA(3) - - - diphenylA(4) - - dlphenylS(3) - +- diphenylS(4)
-3
.
.-'
/ ~
-12 -
computed in MacroModel. There is a fair agreement between the results for the protonated compounds, especially with those obtained at the 4-31G level. The stability order for the neutral compounds is not univocally determined: the STO-3G basis set favors SS, while both 4-3 1G and AMBER favor AS. SS is second in the stability order at the 4-31G level, whereas S A holds that position for AMBER. Simplified ab Initio Models. In order to understand which geometry, and thus which free energy of solution, was more reliable for the neutral system, we decided to carry out geometry optimizations on a simplified model, obtained upon substituting the right-side phenyl ring bound to N2 (see Scheme 1) by a hydrogen atom. Basis sets besides STO-3G and 4-3 1G were used with polarization (6-31G*23)and diffuse (6-3 l+G*24) functions on the heavy atoms, and a full double-{basis set (D9525) was also considered. The results obtained for anti and syn N-phenylguanidine are displayed in Figure 3a. The trend for the STO3G and 4-3 1G basis sets is analogous to the N,N'-diphenylguanidine one, because both couples AS/SS and AA/SA behave as A/S when described by the same basis set. The slope of the D95 results is steeper than that of the 6-31+G* ones (equal in turn to the 4-3 l G ones), whereas it is smoother for the 6-3 lG* values. In any case, all the extended basis sets have as positive a slope as the 4-31G basis set, in that they favor a similar geometry with the N lone pair facing the H on the ring C. The values in the absence of the second phenyl ring are coincident at the 4-3 1G level with those obtained when the second ring is syn to the C=NH group, while they run parallel a t the STO-3G level (Figure 3b). When the second phenyl ring is anti to the C-NH group, there is a change in the slope with respect to the previous results at both levels, as shown also by the lack of symmetry in the plots of Figure la. If the first phenyl ring is substituted by a vinyl group, the C" values (also displayed in Figure 3b) turn out to be very close to the phenylguanidine ones a t the STO-3G level, and the slope is greatly reduced at the 4-31G one. N,N-Dinaphthyl and Dianthracenyl Guanidines. Scheme 2 shows the structure of N,N'-di- I-naphthylguanidine (AS rigid
c'
-b
-13
/-
,
/'
I
I
-15
-
--
a
7 d I
I
The Journal of Physical Chemistry, Vol. 98, No. 21, 1994 5421
Relative Stability of N,N'-Diarylguanidines
conformations can be obtained from the authors upon request. In all the calculations, however, we used the geometries optimized in solution, which nonetheless show very limited variations ( 5 O a t most in the dihedral angles) with respect to those determined in vacuo. The solvent effect seems thus to be linked to a preferential stabilization of the rotamers more than to a geometrical change in their structures. The CmIvalues, displayed in Figure 4, have a similar trend for all the systems (the diphenyl partial charge model results are plotted for comparison), with just a slightly different behavior for the dinaphthyl derivatives. The increase with the system size of the solvent effect is especially marked for the protonated derivatives, which however show a crossing of dinaphthyl and di-9-anthracenyl curves. Also, the cavitation free energy increases with the system size (Figure 5 ) , as the difference between its Pierotti's and Sinanoglu's estimates. For each system a lower value of the cavitation free energy is obtained again when the rotamer is surrounded by a smaller cavity (AA, AA', AA+, AA'+).
SCHEME 3
SCHEME 4
Comparison between the Free Energy Results
geometry) again with the four torsional degrees of freedom (three set to 180°, one to 0 ' (@2)): 01and @Z are defined as for the diphenyl compound, while @3 = Cbl-Crl-N1-Cc and = CbZC,2-NrCc give the deviation of the condensed aromatic rings from the N-Cc-N plane. Cb indicates the'bridged" atoms, which are shared by the rings. In the point charge models of the neutral and protonated dinaphthyl and dianthracenyl guanidines, the charge on these atoms is equal to 0. Also the C, atoms are now neutral, because in MacroModel their charges were assigned to Ni, which thus bears-0.395 and-0.295 au ofchargein theneutral and protonated systems, respectively. The N,N'-di-1-anthracenyl (SA rigid geometry, @ I = @3 = ; ' 0 @2 = a4 = 180'; Scheme 3) and di-9-anthracenyl (SS rigid geometry, cPI = @2 = @3 = cP4 = 0 ' ; Scheme 4) guanidines can be obtained by condensing one further ring to the naphthyl moiety on the same side of that already added to the phenyl ring to get the 1-naphthyl system or on the opposite side, respectively, producing in the latter case a symmetric substituent. The geometries have been optimized both in vacuo and in solution with MacroModel. The torsional angles defining the various
---6-
blu*menyl
The interest is now focused on the point charge model results in the framework of the continuum3 solvent as applied to the ab initio description of the solute and of the continuum water model AMBER//GB/SAIO deriving from the MacroModel9 calculations. Since both are continuum models, in the following we refer to the former as Cont and to the latter as AMBER//GB/ SA. The comparison of the AMBER//GB/SA results with the a b initio free energies has been carried out in a previous section. As pointed out there, the computed quantities include different contributions. Therefore, in order to make them consistent, we have to add the cavitation free energy to the solvent effect for the Cont results to be compared with the AMBER//GB/SA solvent effect. The total effect then includes the energy of the solute (,Tintra), computed in MacroModel and added to both sets of results. (In the case of the diphenyl compounds this energy is not completely correct, as it is for the other derivatives, because only the dihedral angles from MacroModel have been used for the Cont point charge calculations, not the whole set of internal parameters.) The solvent effect and the total free energy change are plotted for the neutral and protonated diary1 derivatives of guanidine in Figures 6-9. The trend is very similar using either Cont or AMBER//GB/SA, with a slightly larger deviation for the SS rotamersof the dinaphthyl- and dianthracenyl-substituted guanidines and for AA di-9-anthracenyl guanidine, which shows this feature also for the protonated derivatives. In summary, however, the agreement between the results is acceptable, because it is difficult to assign a stability order when the energy differences are only about 1-2 kcal/mol. The approximations involved in the methods in fact do not allow an elevated accuracy.
'
0
--a-d-b9r*ncryl
-68 SA AA AA' AS+ SS+ AA+ AA'+ Conformers Conformers Figure 4. Solvent effect, W', on the MacroModel geometries of the (a) neutral and (b) protonated diarylguanidine rotamers, as described by the partial charges of the solute: diphenyl (solid line), di-1-naphthyl (long dashed line), di-1-anthracenyl (short dashed line), and di-9-anthracenyl (dashed line). -16
AS
ss
Alagona et al.
5428 The Journal of Physical Chemistry, Vol. 98, No. 21, 1994
I
1
AS
I
I
I
SS SA AA Conformers
25
AA'
I AS+ SS+ AA+ AA'+ Conformers
Figure 5. Cavitation free energy, Gav, produced by the Sinanoglu's (plain lines) and Pierotti's (marked lines) formulas on the MacroModel geometries
of the (a) neutral and (b) protonated diarylguanidine rotamers: diphenyl (solid line), di-1-naphthyl (long dashed line), di-1-anthracenyl (short dashed line), and di-9-anthracenyl (dashed line).
-
-
- solvCont
-9-
4 -
- - -& -
0
IotCont
6
sol//QB/SA
- U - tol//QB/SA
A
4
0
Em
-
2 -
-b .$(
0
0 -
0
-2
-2-
-4
- solvCont totCont - - - aoI//QB/SA - - tot//QB/SA
2
Y
Q
---e
-
I
I
AS
SS
I
SA Conformers
- - I* AA
-4
-6
AS+
ssc
AA+
Conformers
Figure 6. Free energy differences with respect to AS (or AS+) diphenylguanidineobtained using the AMBER//GB/SA method (short dashed line) and the continuum solvent (PPIv = A w l + AGaV from MGPIPC, see text, long dashed line) for the (a) neutral and (b) protonated compounds. The total free energy differences (AMBER//GB/SA, dashed line, and MGPIPC, solid line) are also compared. Geometries and internal energies are from MacroModel.
6 4
2
0 -2 -4
-6 -8
AS
SS
-10
SA AA AA' Conformers Conformers Figure 7. Free energy differences with respect to AS (or AS+) di-1-naphthylguanidineobtained using the AMBER//GB/SA method (short dashed line) and the continuum solvent (AGwIv= A w l + AGaV from MGPIPC, see text, long dashed line) for the (a) neutral and (b) protonated compounds.
The total free energy differences (AMBER//GB/SA, dashed line, and MGPIPC, solid line) are also compared. Geometries and internal energies are from MacroMcdel. Stacking - Effect
In systems with sterically close aromatic rings special attention should be paid to the stacking effect. According to a general idea (and many experimental results have been explained by that), stacking of aromatic rings stabilizes the system. Calculations for gas-phase dimers of different aromatic systems, however, revealed that for small rings, like the benzene dimer, a T-shape instead of a face-to-face stacked arrangement is favorable.26The electrostatic model of Hunter and SandersZ7gives a rationale for the different Coulombic interactions in theT-shaped and stacked (face-to-face and offset) arrangements. Partitioning of the energy terms indicate the balance of the repulsive Coulombic and attractive van der Waals interactions for these dimers. Due to
the rapidly increasing stabilization upon the van der Waals interaction for larger ring systems, an offset stacked form becomes the most stable one for the benzene-anthracene dimer.26 Even in single systems with satisfactorily distant aromatic rings and with considerable conformational flexibility the rings can accept a T-shaped arrangement in the solid phase.28 The situation is different for the geometrically rigid 1,s-diphenylnaphthalene systems,Z9where an almost face-to-face stacking is forced by the chemical structure. Our systems are similar to these last ones, given a three-atom chain connecting the two aryl rings in both cases. The diarylguanidines are, however, a conformationally more relaxed system, allowing some offset stacked arrangement in the AA and AA' conformations.
Relative Stability of N,N’-Diarylguanidines
The Journal of Physical Chemistry, Vol. 98, No. 21, 1994 5429
”
6 1
1
4 0
c1 z
.\ \
2
m
Y O
w
-2 Q
-4 -6
4
-
‘.*
-4
-
-6
-
‘A
sol//GB/SA
,
-U-
-8 8. I , , , , , i. n - I
AS
SS
SA AA Conformers
AS+
AA’
solvCon1 tot C o n I
I
lol//GE/SA
\
SS+ AA+ Conformers
AA‘+
Figure 8. Free energy differenceswith respect to AS (or AS+) di-1-anthracenylguanidineobtained using the AMBER//GB/SA method (short dashed line) and the continuum solvent (AG@” = AGW1+ AGwvfrom MGPIPC, see text, long dashed line) for the (a) neutral and (b) protonated compounds. The total free energy differences (AMBER//GB/SA, dashed line, and MGPIPC, solid line) are also compared. Geometries and internal energies are from MacroModel. 6
5
4
=
0
Em 2
2
Y O
solvCon1
-2 Q
-10
a-4
.\
-5
w
c
1
IolConl
\
sol//GB/SA
- U - tot//GB/SA
B
-15 SA AA AA‘ AS+ AA+ Conformers Conformers Figure 9. Free energy differenceswith respect to AS (or AS+) di-9-anthracenylguanidineobtained using the AMBER//GB/SA method (short dashed line) and the continuum solvent (AGmIV= AGM1+ AGQv from MGPIPC, see text, long dashed line) for the (a) neutral and (b) protonated compounds. The total free energy differences (AMBER//GB/SA, dashed line, and MGPIPC, solid line) are also compared. Geometries and internal energies are
-6
AS
ss+
SS
from MacroModel.
Consideration of the factors leading to the stabilization upon ring stacking is possible in a quantum mechanical framework by using the energy-partitioning method proposed by Chalasinski and Szczesniak30for the interaction energies of dimers. Such a partitioning is not available, however, a t the MP2 level for a single molecule. Thus, we can only emphasize the importance of correlation effects in the stacked anti-anti conformations. A more correct description of the electron distribution at the MP2 level as compared to the S C F level leads to a considerable stabilization of the AA and AAf structures for the diphenyl derivative (Table 4). Using, however, a molecular mechanics description of the internal energy within the MacroModel package, specific relative terms can be evalulated. The electrostatic energies of the AA and AA’ conformers relative to the most stable AS forms (or SS for the protonated diphenyl and di-1-naphthyl derivatives) increase from 7 to 10 kcal/mol and from 10 to 14 kcal/mol for the neutral and protonated species, respectively, when going from the diphenyl to the dianthracenyl derivatives. The concomitant change in the van der Waals relative energy is -3 to-10 kcal/mol for both the neutral and the protonated species. Thus, while the two effects compensate for the neutral forms with larger rings, stacking is not favored due to internal terms for the protonated species even in the dianthracenyl derivatives. The effect of hydration in stabilizing the stacked forms is not clear. For the benzene dimer26or aromatic host-guest complexes31 direct contacts of the aromatic rings as compared to solventseparated structures were indicated as favorable in water. Although the overall relative solvation energies for the anti-anti forms are negative for any system studied (sol//GB/SA in Figures 6-9), their contributions from the hydrophobic parts raise questions. Cutting out the guanidine moiety and leaving the arrangements of the aromatic rings in the MacroModel optimized
positions, the broken C-N bonds were replaced by C-H and hydration of the resulting hydrophobic parts was calculated. While the hydration energy difference for the AA and SS “cut-outsn is about 0 for the dibenzene model, the hydration of the SS dianthracene is favorable over the AA arrangement by about 1.4-1.8 kcal/mol. The result for the diphenyl derivative agrees with that obtained by the ab initio continuum model (see previously). More favorable hydration of the syn-syn instead of the anti-anti conformation of the two anthracenyl groups means that the GB/SA method favors the solvent-separated hydration of the aromatics. The cavitation free energies are smaller by a few kcal/molfortheAAthanfortheSSformsofthedianthracenyl systems (Figure 5). Net negative relative hydration free energies of the hydrophobic sites may be due to more negative polarization terms upon the hydration of the SS form. Considering the relatively large ring C and H charges in MacroModel (Table 5 ) , this effect is possible; however, we cannot decide at the present level whether this finding is physically correct. Conclusions The STO-3G solvation free energy is about 5-7 kcal/mol less favorable than the 4-3 1G one for the neutral diphenylguanidines, because of the larger electrostatic effect produced by the 4-31G basis set, but only about 2 kcal/mol less favorable for the protonated compounds, which are overall less sensitive to the basis set employed. The inversion in the stabiliy order in solution found at the STO-3G level for the neutral diphenylguanidines is due to the favorable position with respect to the solvent assumed by the N lone pair when the phenyl ring is syn to the =NH group. At the 4-31G level, on the contrary, the H on the nearby Cringlies almost in the same plane as the N lone pair, preventing
5430 The Journal of Physical Chemistry, Vol. 98, No. 21, 1994
any interaction with water. The trend of the results obtained employing both STO-3G//4-3 1G (STO-3G basis set and 4-31G optimized geometry) and 4-3 lG//STO-3G calculations clearly depends on the geometry. Calculations performed on a simplified model with more extended basis sets are in line with the 4-31G results. The MP2 correlation corrections in solution do not alter the relative positioning of the rotamers, even though they raise GSo1 by about 2.6-3.2 and 0.3-1.3 kcal/mol for the neutral and protonated compounds, respectively. Both in vacuo and in solution, the MP2 corrections fill about half of the energy gap between AA and AS (and AA+/AS+) N,N’-diphenylguanidines. The cavitation free energy is nearly independent of the basis set, but not of the conformation: the AA and AA+ structures in fact present a cavitation free energy slightly lower than that of the other rotamers. The values produced for both neutral and protonated diphenyl systems by Sinanoglu’s formula are about 3 kcal/mol larger than those obtained employing Pierotti’s formula, but their differences (as well as their absolute values) increase with the size of the substituents. The AMBER//GB/SA results are in fair agreement with those from the continuum solvent method3 using the 4-3 1G basis set in the SCF calculations for the neutral and protonated N,N’diphenylguanidine. A comparison with the S C F instead of the MP2 corrected values is appropriate considering that SCF calculations were performed when parametrizing the AMBER force field. The agreement between the a b initio continuum and GB/SAcalculated solvation free energies is also acceptable when fixed atomic point charges are used for further diary1 derivatives. Thus, the AMBER//GB/SA method seems to work reasonably well (mainly for the ionic species) in calculating relative conformational energies/free energies for structures without intramolecular hydrogen bonds32in aqueous solution and may be useful to predict these terms for molecules in which a large number of atoms prevent calculations a t the a b initio level. The stacking of the aromatic rings is disfavored by electrostatic and favored by van der Waals interactions. The stabilization produced by the latter term rapidly increases with the number of aromatic rings involved. The hydration of these rings in a nonstacked and rather separated form is slightly favored by the GB/SA method in contrast to available solution simulations with explicit water models.
Acknowledgment. P. Nagy is grateful to Cambridge Neuroscience for financial support. References and Notes (1) (a) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; O U P New York, 1991. (b) Hansen, J.-P.;McDonald, I. R. TheoryofSimple Liquids;A P London, 1990. (c) Warshel, A. Computer Modelingof Chemical Reactions in Enzymes and Solutions; Wiley: New York, 1991. (d) Tapia, 0. In Quantum Theory of Chemical Reactions; Daudel, R., Pullman, A., Salem, L., Veillard, A,, Eds.; Reidel: Dordrecht, 1980; Vol. 2, p 25. (e) Rivail, J.-L. In New Theoretical Concepts for Understanding Organic Reactions; Bertritn, J., Csizmadia, I. G., Eds.; Kluwer: Dordrecht, 1989; p 219. (0 Cramer, C. J.; Truhlar, D. J. J . Am. Chem. SOC.1991,113, 8305; 1991, 113.9901. (2) (a) Jorgensen, W. L.; Ravimohan, C. J. Chem. Phys. 1985,83,3050. (b) Pearlman, D. A.; Kollman, P. A. J. Chem. Phys. 1989, 90, 2460. (c) Nagy, P. I.; Dunn, W. J., 111; Alagona, G.; Ghio, C. J. Am. Chem. Soc. 1991, 113, 6719; 1992, 114, 4752. (d) Nagy, P. I.; Dunn, W. J., 111; Alagona, G.;
Alagona et al. Ghio, C. J. Phys. Chem. 1993, 97, 4628. (3) (a) Miertus, S.;Scrocco, E.; Tomasi, J. Chem. Phys. 1981,55, 117. (b) Bonaccorsi, R.; Cimiraglia, R.; Tomasi, J. J. Compur. Chem. 1983,4,567. (c) Bonaccorsi, R.; Cimiraglia, R.; Tomasi, J. Chem. Phys. Lett. 1983, 99, 77. (d) Pascual-Ahuir, J. L.; Silla, E.; Tomasi, J.; Bonaccorsi, R. J . Compur. Chem. 1987,8, 778. (4) (a) Alagona, G.; Ghio, C. J . Mol. Struct. (THEOCHEM) 1992,254, 287. (b) Alagona, G.; Ghio, C. J . Mol. Liq., in press. ( 5 ) Alagona, G.; Ghio, C. J. Mol. Struct. (THEOCHEM) 1992, 256, 187. (6) (a) Greenhill, J. V.; Ping Lue. In Progress in Medicinal Chemistry; Ellis, G. P., Luscombe, D. K., Eds.; Elsevier: Amsterdam, 1993; Vol. 30, p 203. (b) Scherz, M. W.; Fialeix, M.; Fischer, J. B.; Reddy, N. L.; Server, A. C.; Sonders, M.; Tester, B. C.; Weber, E.; Wong, S.T.; Keana, J. F. W. J. Med. Chem. 1990, 33, 2421. (c) Reddy, N. L.; Hu, L.-Y.; Cotter, R. E.; Fischer, J. B.; Wong, W. J.; McBurney, R. N.; Weber, E.; Holmes, D. L.; Wong, S.T.; Prasad, R.;Keana, J. F. W. J. Med. Chem. 1994, 37, 260. (7) Boudon, S.;Wipff, G.; Maigret, B. J. Phys. Chem. 1990, 94,6056. (8) Muller, G. W.; Walters, D. E.; DuBois, G. E. J . Med. Chem. 1992, 35, 740. (9) MacroModel V3.5a: Mohamadi, F.; Richards, N. G. J.; Guida, W. C.; Liskamp, R.;Caufield, C.; Chang, G.; Hendrickson, T.; Still, W. C. J. Comput. Chem. 1990, 11,440. (10) Still, W. C.; Tempczyk, A.; Hawley, R. C.; Hendrickson, T. J . Am. Chem. SOC.1990, 112,6127. (1 1) Hehre, W. J.; Stewart, R. F.; Pople, J. A. J. Chem. Phys. 1969.51, 2657. (12) Ditchfield, R.; Hehre, W. J.; Pople, J. A. J . Chem. Phys. 1971, 54, 724. (13) Frisch, M. J.; Trucks, G. W.; Head-Gordon, M.; Gill, P. M. W.; Wong, M. W.; Foreman, J. B.; Johnson, B. G.; Schlegel, H. B.; Robb, M. A.; Replogle, E. S.;Gomperts, R.; Andres, J. L.; Raghavachari, K.; Binkley, J.S.;Gonzalez,C.;Martin,R.L.;Fox,D. J.;Defrces,D. J.;Baker, J.;Stewart, J. J. P.; Pople, J. A. Gaussian9Z; Gaussian Inc.: Pittsburgh, PA, 1992. (14) Ghio, C. MGPIPC (ICQEM-CNR, Pisa, Italy), the point charge version of MGPISA, Bonaccorsi, R.; Cammi, R. (ICQEM-CNR, P i a , Italy). (IS) Peterson, M. R.; Poirier, R. A. Dept. of Chemistry, University of Toronto, Toronto, Ontario, Canada. (16) (a) Maller, C.; Plesset, M. S.Phys. Rev. 1934, 46, 618. (b) Pople, J. A,; Binkley, J. S.; Seeger, R. Int. J. Quantum Chem. 1976, 10s, 1. (c) Krishnan, R.; Frisch, M. J.; Pople, J. A. J . Chem. Phys. 1980, 72, 4244. (17) (a) Persico, M.; Cossi, M. Dip. di Chimica e Chimica Industriale, Universita di Pisa, Pisa, Italy. (b) C o d , M. Thesis, Universita di Pisa, 1991. (18) Clementi, E., Ed.; Modern Techniques in Computational Chemistry: MOTECC-90; ESCOM; Leiden, 1990. (19) (a) Weiner, S.J.; Kollman, P. A.; Case, D. A.; Chandra Singh, U.; Ghio, C.; Alagona, G.; Profeta, S.,Jr.; Weiner, P. J. Am. Chem. Soc. 1984, 106, 765. (b) Weiner, S.J.; Kollman, P. A.; Nguyen, D. T.; Case, D. A. J. Comput. Chem. 1986, 7, 230. (c) Chandra Singh, U.; Kollman, P. A. J . Comput. Chem. 1984, 5, 129. (20) Pierotti, R. A. Chem. Rev. 1976, 76, 717. (21) Sinanoglu, 0. Theor. Chim. Acta 1974, 33, 279. (22) Birnstock, F.; Hofmann, H.-J.; Kdhler, H.-J. Theor.Chim.Acta 1976, 42, 311. (23) (a) Hehre, W. J.; Ditchfield, R.; Pople, J. A. J. Chem. Phys. 1972, 56,2257. (b) Hariharan,P.C.;Pople, J.A. Theor. Chim. Acfa1973,28,213. (24) Clark, T.; Chandrasekhar, J.; Spitznagel, G. W.;von RaguCSchleyer, P.J . Compur. Chem. 1983, 4, 294. (25) Dunning, T. H.; Hay, P. J. Modern Theoretical Chemistry;Plenum: New York, 1976; Chapter 1, pp 1-28. (26) Jorgensen, W. L.; Severance, D. L. J. Am. Chem. SOC.1990, 112, 4768. (27) (a) Hunter, C. A.; Sanders, J. K. M. J. Am. Chem. Soc. 1990,112, 5525. (b) Hunter, C. A. Angew. Chem., Inr. Ed. Engl. 1993, 32, 1584. (28) Declerq, D.; Delbeke, P.; De Schryver, F. C.; Van Meervelt, L.; Miller, R. D. J. Am. Chem. Soc. 1993,115, 5702. (29) (a) Cozzi, F.; Cinquini, M.; Annunziata, R.; Dwyer, T.; Siegel, J. S. J . Am. Chem. SOC.1992,114,5729. (b) Cozzi, F.; Cinquini, M.; Annunziata, R.; Siegel, J. S.J. Am. Chem. Soc. 1993, 115, 5330. (30) Chalasinski, G.; Szczesniak, M. M. Mol. Phys. 1988, 63, 205. (31) Smithrud, D. B.; Diederich, F. J. Am. Chem. Soc. 1990, 112, 339. (32) For shortcomingsof the method related tostructures forming internal hydrogen bonds in the gas phase, see: Nagy, P. I.; Bitar, J. E.;Smith, D. A. J . Comput. Chem., in press.