Which Functional Form Is Appropriate for ... - ACS Publications

counterpoise BSSE corrections at the MP2/6-31G** level reproduce satisfactorily the corrected dimerization energies of formamide and acetamide around ...
1 downloads 0 Views 67KB Size
Download PDF
J. Phys. Chem. B 2000, 104, 8321-8326

8321

Which Functional Form Is Appropriate for Hydrogen Bond of Amides? Young Kee Kang* Department of Chemistry, Chungbuk National UniVersity, Cheongju, Chungbuk 361-763, Korea ReceiVed: February 29, 2000; In Final Form: May 22, 2000

Four types of the functional form, which are the Morse, Lennard-Jones 10-12, 6-12, and 6-9 types of the potential function used frequently to describe hydrogen bonds, are tested against ab initio MP2/6-31G** calculations for their ability to describe dimerization energies of N-methylacetamide (NMA), formamide, and acetamide, varying the distance between amide hydrogen and carbonyl oxygen. All parameters for four types of the potential fitted to the scaled dimerization energies of NMA around the minimum by including counterpoise BSSE corrections at the MP2/6-31G** level reproduce satisfactorily the corrected dimerization energies of formamide and acetamide around the minima. The Morse function results in the best fit to the scaled MP2/6-31G** energies of three amides both at short and long distances, whereas other functions are satisfactory only at long distances.

1. Introduction The N-H‚‚‚OdC hydrogen bond plays an important part in determining the conformations of proteins.1 Several potential functions for hydrogen bonds of amides, peptides, and proteins have been employed for empirical force fields such as CFF,2 ECEPP,3 CHARMM,4 and AMBER5 in order to account for the exchange repulsion and other second-order components in hydrogen bonding explicitly. These force fields utilized the Morse type2 and Lennard-Jones (L-J) 10-12 type3-5 of potential for hydrogen bonds. On the other hand, some force fields such as the earlier version and class II CFF,6,7 OPLS,8 MM2,9 MM3,10 and recent versions of AMBER11 and CHARMM12 have employed the L-J 6-12,6,8,11,12 6-9,6,7 or 6-exponential9,10 types of nonbonded (or van der Waals) potential for hydrogen bonds. In addition, the angular dependence of the hydrogen bond with the N-H‚‚‚O or H‚‚‚OdC angular terms has been incorporated into the hydrogen bond potentials.2,4,10 Recently, No et al. proposed the hydrogen bond potential in which the angular dependence of the energy of hydrogen bonds is described with only L-J 6-12 type potential functions without additional functions involving bond angles.13 Interestingly, Isaacs et al. reported recently a Compton profile anisotropies of ordinary ice Ih, which strongly supports a partially covalent nature of the hydrogen bond.14 In addition, Dannenberg et al. suggested the extent to which the hydrogen bond can be considered to be covalent by ab initio molecular orbital calculations in graded uniform electric fields of molecules that form chains of hydrogen bonds.15 We report here the parameters of four hydrogen bond potentials fitted to dimerization energies of N-methylaceatmide (NMA) at the MP2/6-31G** level and test their reproducibility for dimerization energies of formamide and acetamide from MP2/6-31G** calculations. 2. Computational Methods 2.A. Quantum Mechanical Calculations. All ab initio and density functional calculations were carried out using the * Phone: +82-43-261-2285. [email protected].

FAX:

+82-43-273-8328.

E-mail:

Gaussian 9416 and Gaussian 9817 packages. Geometry optimizations for monomers and dimers of formamide, acetamide, and NMA were carried out at HF, B3LYP, and MP2 levels of theory with the 6-31G** basis set. The dimer structures of formamide and acetamide optimized at the HF/6-31G** level13 were used as initial structures for optimizations of their monomers and dimers at HF, B3LYP, and MP2 levels. The structure of NMA from electron diffraction measurements18 was used as an initial structure for optimizations of its monomer and dimer at these levels. The dimerization energies of formamide, acetamide, and NMA optimized at HF, B3LYP, and MP2 levels were corrected for the basis set superposition error (BSSE) using the counterpoise method.19 The dimerization energies of each amide were calculated at the MP2/6-31G** level only varying the hydrogen bond distance r(H‚‚‚O) between amide hydrogen and carbonyl oxygen from 1.7 to 2.7 Å by the step of 0.1 Å, with other structural parameters fixed at the optimized values of its dimer at the MP2/6-31G** level. Then, these MP2/6-31G** dimerization energies of each amide were scaled by a factor, which was given by a ratio of dimerization energies with and without BSSE corrections for its optimized dimer at the MP2/6-31G** level, and used to derive the parameters for four types of potential for hydrogen bond. 2.B. Empirical Potential Energy Functions. The intermolecular interaction energy between two monomers A and B is given by the sum of the electrostatic energy (EES), the nonbonded energy (ENB), and the hydrogen bond energy (EHB). All electrostatic and nonbonded interaction energies are calculated using the 1-6-12 potential given by onA onB

EES + ENB )

∑i ∑j (qiqj/rij + Aij/r12ij - Cij/r6ij)

(1)

The partial atomic charges, q, were calculated by the GDAC (geometry-dependent net atomic charges) method of Cho et al.,20 and rij is the distance between the two interacting atoms i on monomer A and j on monomer B. The GDAC method is extended from the MPEOE (modified partial equalization of orbital electronegativity) method of No et al.21-24 In the GDAC

10.1021/jp000772h CCC: $19.00 © 2000 American Chemical Society Published on Web 08/04/2000

8322 J. Phys. Chem. B, Vol. 104, No. 34, 2000

Kang

TABLE 1: Parameters for Atomic Charges and Nonbonded Interactions for Amides charge parametersb

and

EHB,L-J )

nonbonded parametersc

atom typea

ai

bi

ii

rii0

H1 H2 C1 C2 N O

2.648 2.648 2.810 2.382 3.093 5.727

0.429 0.429 0.432 22.999 -0.810 0.646

0.031 0.094 0.042 0.139 0.157 0.226

2.95 2.33 4.15 3.45 3.38 3.05

(m -n n)(r ) m (r ) B )( m - n)

D k ) k

Defined as H1, H2, C1, C2, N, and O for aliphatic hydrogen, amide hydrogen, aliphatic carbon, carbonyl carbon, amide nitrogen, and carbonyl oxygen, respectively. b GDAC parameters taken from ref 20. c Taken from ref 25; units of  and r 0 are kcal/mol and Å, respectively. ii ii

TABLE 2: Dipole Moments and Their Components for Formamide, Acetamide, and NMAa amide

methodb

µ

µx

µy

µz

formamide

GDAC B3LYP exptlc GDAC B3LYP exptld GDAC B3LYP exptld

3.67 3.82 3.71 3.70 3.77 3.75 3.67 3.71 3.71

3.57 3.75 3.61 -0.08 -0.19

0.85 0.74 0.85 -3.70 -3.77

0.00 0.00 0.00 0.00 0.00

0.62 0.66

-3.61 -3.65

0.00 0.00

NMA

a Units in Debye. µ is the total dipole moment, and their components are presented by µx, µy, and µz. b GDAC values were calculated for the optimized structures at the MP2/6-31G** level using parameters of Cho et al. in ref 20. B3LYP values are those optimized at the B3LYP/ 6-31G** level, taken from ref 20. c Taken from Kurland, R. J.; Wilson, E. B. J. Chem. Phys. 1957, 27, 585. d Taken from Meighan, R. M.; Cole, R. H. J. Phys. Chem. 1964, 68, 503.

method, the damping factor fkl of a bond k-l is calculated by the equation

fkl ) 1 - dkl/(rw,k + rw,l)

k

k

0 m k

(5)

0 n k

(6)

Here, k and rk0 are the depth and minimum distance of the potential for the kth hydrogen bond. The Morse type of potential has another parameter, Rk, as shown in eq 3. The earlier CFF force field2 utilized the Morse type and the ECEPP,3 CHARMM,4 and AMBER5 force fields adopted the L-J 10-12 type of potential for hydrogen bond. Some force fields such as the previous version and class II CFF,6,7 OPLS,8 and recent versions of AMBER11 and CHARMM12 have employed the L-J 6-12,6,8,11,12 and 6-96,7 types of nonbonded (or van der Waals) potential for hydrogen bonds. Therefore, four types of hydrogen bond potential including the Morse, L-J 10-12, 6-12, and 6-9 potentials are considered in this work. 2.C. Optimization of Hydrogen Bond Parameters. The scaled dimerization energies of NMA at the MP2/6-31G** level varying the hydrogen bond distance r(H‚‚‚O) from 1.7 to 2.7 Å by the step of 0.1 Å, which are discussed above in section 2.A, were used to derive the parameters for four types of potential for hydrogen bond. The objective function S was the weighted sum of the squared deviations between the empirical and scaled MP2/6-31G** dimerization energies of NMA so that 11

S)

ωi(Eempirical,i - EMP2,i)2 ∑ i)1

(7)

(2)

Here, dkl is the length of a bond connecting two atoms k and l, whose van der Waals radii are rw,k and rw,l, respectively.20 The electronegativity parameters ai and bi were optimized so as to reproduce dipole moments and their components for several amides optimized at the B3LYP/6-31G** level,20 which are listed in Table 1. The dipole moments and their components calculated for amides optimized at the MP2/6-31G** level using GDAC parameters of Cho et al.20 are shown in Table 2 and compared to corresponding B3LYP and experimental values. The nonbonded parameters Aij and Cij may also be expressed in terms of L-J ij and σij as Aij ) 4ijσij12 and Cij ) 4ijσij6. The distance rij0 of the minimum of the nonbonded potential can be expressed as rij0 ) 21/6σij. The nonbonded parameters ii and rii0 used here are taken from ref 25, which were derived to reproduce lattice parameters and lattice energies of a number of molecular crystals, and they are listed in Table 1. The geometric and arithmetic combination rules were applied for the parameters ij and rij0 for a pair of unlike atoms, respectively. For interactions between amide hydrogen and carbonyl oxygen of the amide dimer, the explicit hydrogen bond energy EHB term replaced the nonbonded energy term in eq 1. The functional forms for the hydrogen bond tested here are the Morse and L-J n-m potentials given by

EHB,Morse )

(4)

where

a

acetamide

∑k (Dk/rmk - Bk/rnk)

∑k k{exp(-2Rk(rk - r0k)) 2 exp(-Rk(rk - r0k ))} (3)

where Eempirical,i is the dimerization energy of NMA calculated by the sum of EES, ENB, and EHB in eqs 1 and 3 (or 4) at the ith distance, and EMP2,i is the corresponding scaled MP2/6-31G** dimerization energy. To reproduce better the scaled MP2/631G** energies around the minimum, i.e., r(H‚‚‚O) ) 2.0 Å, the weighting factor ωi of 108 was employed for r(H‚‚‚O) ) 2.0 Å; ωi ) 106 for r(H‚‚‚O) ) 1.9, 2.1, and 2.2 Å, and ωi ) 1 for other distances. The secant-type unconstrained minimization problem solver (SUMSL) algorithm26 was used for the minimization of the objective function S. The energy Eempirical,i of NMA was calculated for its dimer conformation at the ith distance, which was generated from the monomer structure optimized at the MP2/6-31G** level, with the hydrogen bonding geometries of the dimer optimized at the MP2/6-31G** level. The empirical energies for dimers of formamide and acetamide were calculated for the conformations generated by the same method applied to NMA. 3. Results and Discussion 3.A. Optimized Dimer Structures. Colominas et al. reported recently that the symmetrically double hydrogen-bonded dimer of formamide is expected to be the major form in the gas phase from B3LYP/6-31G* calculations on different conformations of formamide.27 However, we focused here only on the single hydrogen-bonded dimers of formamide, acetamide, and NMA because our primary aim of this work is to derive the parameters of hydrogen bond potentials from the dimerization energies of NMA and to apply them to the conformational energy calcula-

Hydrogen Bond Potential for Amide Dimers

J. Phys. Chem. B, Vol. 104, No. 34, 2000 8323 TABLE 3: Structural Parameters for Hydrogen Bonds of Optimized Amide Dimersa method

r(N-H) r(H‚‚‚O) θ(N-H‚‚‚O) θ(H‚‚‚OdC)

HF/6-31G** B3LYP/6-31G** MP2/6-31G**

0.996 1.017 1.012

formamide 2.079 1.966 1.990

175.7 177.9 172.8

144.1 126.8 121.8

HF/6-31G** B3LYP/6-31G** MP2/6-31G**

0.995 1.013 1.008

acetamide 2.086 1.982 1.998

171.9 167.8 167.0

163.6 151.0 153.4

HF/6-31G** B3LYP/6-31G** MP2/6-31G**

N-methylacetamide 0.995 2.083 177.5 1.013 1.983 170.8 1.011 1.977 170.3

174.0 160.3 148.8

a Units are Å and degrees for distances and angles, respectively. See Figure 1 for the definition of structural parameters.

Figure 1. Structures for dimers of (A) formamide, (B) acetamide, and (C) NMA optimized at the MP2/6-31G** level.

tions of peptides and proteins rather than to search the lowest energy conformations of small amides. The dimer structures of formamide, acetamide, and NMA optimized at the MP2/6-31G** level are shown in Figure 1. For the orientation of the (C′)CH3 group in acetamide, the syn conformation with one hydrogen of CH3 eclipsing the C′(O)-N bond was found to be preferred from microwave experiments on its monomer28 and HF/6-31G** optimizations on its dimer.13 In addition, the syn conformation of acetamide monomer was proved to be more stable than the anti conformation with two hydrogens of CH3 staggered with respect to the C′(O)-N bond from optimizations at B3LYP and MP2 levels in this work. However, the optimization of acetamide dimer with syn (C′)CH3 groups at the B3LYP level led to the dimer conformation with anti (C′)CH3 groups, which was also assumed to be the initial conformation for the optimization of the acetamide dimer at the MP2 level. X-ray diffraction experiments on NMA29 indicate that the preferred conformation is anti, which was also

confirmed by previous HF/DZP30 and HF/6-31++G**31 optimizations of the NMA dimer. The structural parameters for hydrogen bonds of the dimers optimized at HF, B3LYP, and MP2 levels with the 6-31G** basis set are summarized in Table 3. The increase in the bond length r(N-H) due to the formation of hydrogen bond is found to be ∼0.005, ∼0.01, and ∼0.01 Å from HF, B3LYP, and MP2 calculations, respectively. The hydrogen-bonded bond lengths r(N-H) of the dimers at the HF level are shorter by about 0.02 Å than those at B3LYP and MP2 levels. The differences in values of r(N-H) from B3LYP and MP2 calculations are within 0.005 Å. In particular, the hydrogen-bond distances r(H‚‚‚O) from HF calculations are longer by about 0.1 Å than those from B3LYP and MP2 calculations. The values of r(H‚‚‚O) for dimers of formamide and acetamide at the B3LYP level are calculated to be shorter than those at the MP2 level, but the reverse trend is obtained for the NMA dimer. The bond angles θ(N-H‚‚‚O) and θ(H‚‚‚OdC) optimized at B3LYP and MP2 levels are similar to each other. The differences in r(H‚‚‚O), θ(N-H‚‚‚ O), and θ(H‚‚‚OdC) from HF, B3LYP, and MP2 calculations may indicate that the electron correlation plays a role in determining the hydrogen-bonded structures of amides. By analyzing X-ray or neutron diffraction data of organic amide crystals, Taylor et al.32 reported the mean values of r(NH), r(H‚‚‚O), and θ(N-H‚‚‚O) to be 1.030 ( 0.016 Å, 1.934 ( 0.123 Å, and 158.3° ( 15.6°, respectively, which are consistent with our calculations at B3LYP and MP2 levels. In addition, they suggested the statistically significant tendency for hydrogen bonds to occur in the directions of the conventionally viewed oxygen sp2 lone pairs,33 which are also confirmed by our B3LYP and MP2 computations on amide dimers. In particular, Baker and Hubbard34 reported the values of r(H‚‚‚O), θ(N-H‚‚‚O), and θ(H‚‚‚OdC) to be 2.00 Å, 156°, and 147° by analyzing 15 proteins with high resolution (1.41.8 Å), which agree with the values obtained from organic amide crystals32,33 and B3LYP and MP2 calculations in this work. Although there are no significant differences in the dimer structures optimized at B3LYP and MP2 levels, the structures optimized at the MP2 level were used in deriving the potential parameters for the hydrogen bond because the bond angles θ(N-H‚‚‚O) and θ(H‚‚‚OdC) optimized at the MP2 level are closer to those of formamide and NMA in crystals.6 3.B. Dimerization Energies of Amides. Table 4 lists the dimerization energies of formamide, acetamide, and NMA. The geometries were optimized at HF, B3LYP, and MP2 levels with the 6-31G** basis set, and the single point energies were

8324 J. Phys. Chem. B, Vol. 104, No. 34, 2000

Kang

TABLE 4: Dimerization Energies of Amides with Single Hydrogen Bondsa formamideb

acetamideb

N-methylacetamideb

method

w/o BSSE

w/BSSE

w/o BSSE

w/BSSE

w/o BSSE

w/BSSE

HF/6-31G**//HF/6-31G** B3LYP/6-31G**//HF/6-31G** MP2/6-31G**//HF/6-31G** B3LYP/6-31G**//B3LYP/6-31G** MP2/6-31G**//MP2/6-31G**

-6.34 -7.27 -7.52 -7.56 -7.80

-5.81 -5.77 -5.73 -6.26 -6.08

-6.73 -7.76 -8.36 -7.94 -8.43

-5.86 -5.63 -5.90 -5.82 -6.04

-6.71 -7.71 -8.76 -7.90 -9.13

-5.78 -5.61 -6.26 -5.72 -6.48

a All energies are in kcal/mol. See Figure 1 for the optimized conformation of each amide at the MP2/6-31G** level. b Left and right values for each amide correspond to the values without and with corrections for BSSE using the counterpoise method.

calculated at B3LYP and MP2 levels for dimer structures optimized at the HF level. The dimerization energies of formamide, acetamide, and NMA optimized at the HF/6-31G** level are -6.34, -6.73, and -6.71 kcal/mol, respectively. The dimerization energy of NMA is quite similar to the value of -6.5 kcal/mol optimized at the HF/DZP level,30 although there are small differences in dimer structures. The corresponding dimerization energies are calculated to be -7.56, -7.94, and -7.90 kcal/mol at the B3LYP level, and -7.80, -8.43, and -9.13 kcal/mol at the MP2 level, respectively. This indicates that there are decreases in dimerization energy by ∼1.2 kcal/mol by including the electron correlation at the B3LYP level and that the MP2 level results in the electron correlation effect larger than 1.5 kcal/ mol on dimerization energies. In particular, it should be noted that the shorter distance r(H‚‚‚O) and the decrease in angles θ(N-H‚‚‚O) and θ(H‚‚‚OdC) may contribute to lowering dimerization energies at B3LYP and MP2 levels (see Table 3). The single point energies at B3LYP and MP2 levels for optimized dimers at the HF level underestimate dimerization energies by ∼0.3 kcal/mol compared to those optimized at B3LYP and MP2 levels, which may be also ascribed to differences in hydrogen bond structures shown in Table 3. The value of -8.76 kcal/mol for NMA at the MP2/6-31G**//HF/ 6-31G** level agrees with that of -8.9 kcal/mol at the MP2/ aug-cc-pVDZ//HF/DZP.30 The dimerization energies with the BSSE correction are also shown in Table 4. The inclusion of BSSE correction at HF, B3LYP, and MP2 levels weakens the dimerization energies by 0.5-0.9, 1.3-2.2, and 1.7-2.7 kcal/mol for formamide, acetamide, and NMA, respectively. The effect of BSSE correction on the dimerization energies of amides increases as the molecule becomes larger, which appears to be remarkable at the MP2 level followed by the B3LYP level. As the level of theory becomes higher, the dimerization energies without the BSSE correction become larger, but the contributions of BSSE correction become more significant. As a result, the differences in the dimerization energies of each amide at HF, B3LYP, and MP2 levels become smaller. For example, the dimerization energies of NMA without the BSSE correction are -6.71, -7.90, and -9.13 kcal/mol at HF, B3LYP, and MP2 levels, respectively, but they lead to -5.78, -5.72, and -6.48 kcal/mol by including the BSSE correction. In addition, the inclusion of BSSE correction results in similar dimerization energies for formamide, acetamide, and NMA; e.g., -6.08, -6.04, and -6.48 kcal/mol at the MP2 level, respectively. This may imply the possibility to derive the hydrogen bond potential that can be commonly used for formamide, acetamide, and NMA. The dimerization energy of -5.78 kcal/ mol for NMA at the HF/6-31G** level with the BSSE correction agrees with the value of -5.8 kcal/mol at the HF/DZP,30 but the value of -6.26 kcal/mol at the MP2/6-31G**//HF/6-31G** level has a difference of 0.9 kcal/mol from the value of -7.2 kcal/mol at the MP2/aug-cc-pVDZ//HF/DZP level.30

TABLE 5: Optimized Hydrogen Bond Parametersa potential

k

rk0

Rk

rms devb

Morse 10-12 6-12 6-9

0.659 0.389 0.450 0.501

2.347 2.115 2.168 2.210

1.631

0.216 4.778 3.157 2.039

a Hydrogen bond parameters were optimized for the scaled MP2/631G** dimerization energies for NMA varying the hydrogen bond distance r(H‚‚‚O) from 1.7 to 2.7 Å by the step of 0.1 Å. Units of k and rk0 are kcal/mol and Å, respectively. Rk is a dimensionless quantity. See eqs 3 and 4 in the text for definition of parameters. b Rms dev is the root-mean-square of the minimum S in eq 7.

For the dimerization energies of each amide calculated at the MP2 level only varying the hydrogen bond distance r(H‚‚‚O) from 1.7 to 2.7 Å by the step of 0.1 Å, the BSSE correction was computed at each distance. However, the distance r(H‚‚‚ O) of the minimum dimerization energy for each amide was shifted to 2.1 Å by including the BSSE correction, which is longer by 0.1 Å than the distance optimized at the MP2 level. Because it is not feasible to locate the dimer structure at the energy minimum with the BSSE correction, the effect of the BSSE correction on the dimerization energy at each distance is introduced by scaling the MP2 dimerization energy by a factor, which is given by a ratio of dimerization energies with and without BSSE corrections for its optimized dimer at the MP2 level, as described above in section of 2.A. These scaling factors are calculated to be 0.779, 0.716, and 0.710 for dimers of formamide, acetamide, and NMA, respectively, using the values of Table 4 optimized at the MP2/6-31G**//MP2/6-31G** level. 3.C. Optimized Hydrogen Bond Parameters. The scaled MP2 dimerization energies of NMA along the distance r(H‚‚‚ O) so as to include the BSSE correction, described in the previous paragraph, were used to derive the parameters of potentials for the N-H‚‚‚OdC hydrogen bond and tested for their ability to describe dimerization energies of formamide and acetamide. NMA is the simplest model for the peptide linkage in peptides and proteins and has been extensively studied in deriving the potential parameters for peptides and proteins.6-8,11,12,35 The optimized parameters of the Morse, L-J 10-12, 6-12, and 6-9 types of potential for the NMA dimer are listed in Table 5. The depth k and minimum rk0 for the Morse potential in eq 3 are calculated to be 0.659 kcal/mol and 2.347 Å, respectively, and the optimized parameter Rk is 1.631. These parameters k and rk0 for the Morse potential are larger than those optimized for the L-J 10-12, 6-12, and 6-9 types of potential. It is found that the values of k and rk0 increase as the power of the repulsive component decreases. The values of k and rk0 for the L-J 1012 potential are computed to be smaller than those of the L-J 6-12 potential. Calculated dimerization energies of NMA along the hydrogen bond distance r(H‚‚‚O) using the optimized hydrogen bond parameters in Table 5 are shown in Figure 2A. The scaled MP2/

Hydrogen Bond Potential for Amide Dimers

Figure 2. Calculated dimerization energies of (A) NMA, (B) formamide, and (C) acetamide along the hydrogen bond distance r(H‚‚‚O) using the optimized parameters for NMA in Table 4. The dimerization energies computed by the Morse, L-J 10-12, 6-12, and 6-9 types of potential are represented by four types of line such as s, ‚‚‚, - - -, and - ‚ -, respectively. The scaled MP2/6-31G** dimerization energies by including the BSSE correction are represented by open circles. See the text for scaling MP2/6-31G** energies.

6-31G** dimerization energies by including the BSSE correction are also presented by open circles in Figure 2A. The rms deviations of the minimum S in eq 7 for NMA are calculated to be 0.2162, 4.7776, 3.1573, and 2.0385 for the Morse, L-J 10-12, 6-12, and 6-9 potentials, respectively (Table 5). Based upon the values of S, the Morse potential is the best among four types of potential and followed by the L-J 6-9, 6-12, and 10-12 potentials. As shown in Figure 2A, the Morse potential (represented by a solid line) results in the best fit to the scaled MP2/6-31G** energies of NMA (represented by open circles) at both short and long distances, whereas other potentials are satisfactory only at long distances. This also indicates that the

J. Phys. Chem. B, Vol. 104, No. 34, 2000 8325 L-J types of potential are too repulsive at the distances shorter than the minimum. To escape the abrupt changes in energy, MM2 and MM3 force fields have adopted the 6-exponential type of potential for the hydrogen bond.9,10 The Morse potential has been known to reasonably express the vibrational levels of diatomic molecules36 and the potential energies of the hydrogen donor and acceptor bonds for hydrogen bond.37 In particular, Scheiner and Duan have tested several functional forms of potential against ab initio calculations for their ability to accurately simulate the energetics of proton transfer and a pair of Morse functions appeared best in that they closely mimic ab initio results.38 The success in use of the Morse potential for the hydrogen bond and proton transfer may imply that the formation and breaking of a hydrogen bond appears to be a relatively gradual process, as seen from our scaled MP2/ 6-31G** calculations for NMA in Figure 2A, and that the hydrogen bond has the covalent character, which is consistent with the recent Compton profile anisotropy experiments on ice Ih14 and ab initio calculations in graded uniform electric fields of molecules that form chains of hydrogen bonds.15 Nevertheless the L-J types of potential have been widely used in conformational energy calculations of amides, peptides, and proteins.3-8 According to our results on NMA in Figure 2A, the validity for use of the L-J potential for hydrogen bond should be restricted at the distances of r(H‚‚‚O) equal to and longer than the minimum, in which all types of potential yield almost the same dimerization energies of NMA. Ferguson and Kollman have suggested that the L-J 6-12 type of potential can be replaced by the 10-12 form in molecular mechanics calculations.39 Calculated dimerization energies of formamide and acetamide along the hydrogen bond distance r(H‚‚‚O) using the optimized parameters for NMA in Table 5 are shown in Figures 2B and 2C, respectively, in which the scaled MP2/6-31G** dimerization energies are represented by open circles. As found for NMA in Figure 2A, the Morse potential (represented by solid lines) gives the best fit to the scaled MP2/6-31G** energies of formamide and acetamide at both short and long distances, whereas other potentials are satisfactory only at long distances. The minimum dimerization energies of formamide and acetamide using the optimized Morse parameters are calculated to be -5.65 and -5.80 kcal/mol at r(H‚‚‚O) ) 1.94 and 1.98 Å, respectively. The corresponding MP2/6-31G** energies with the BSSE correction are -6.08 and -6.04 kcal/mol at r(H‚‚‚O) ) 1.97 and 1.98 Å, respectively (see Tables 3 and 4). This indicates that the Morse type of potential for hydrogen bond optimized for NMA can reproduce the dimerization energies of formamide and acetamide within 0.2-0.4 kcal/mol. 4. Conclusions Four types of the functional form, which are the Morse, Lennard-Jones 10-12, 6-12, and 6-9 types of the potential function used frequently to describe hydrogen bonds, are tested against ab initio MP2/6-31G** calculations for their ability to describe dimerization energies of N-methylacetamide (NMA), formamide, and acetamide, varying the distance between amide hydrogen and carbonyl oxygen. All parameters for four types of the potential fitted to the scaled dimerization energies of NMA around the minimum by including counterpoise BSSE corrections at the MP2/6-31G** level reproduce satisfactorily the corrected dimerization energies of formamide and acetamide around the minima. The Morse function results in the best fit to the scaled MP2/6-31G** energies of three amides both at short and long distances, whereas other functions are satisfactory only at long distances.

8326 J. Phys. Chem. B, Vol. 104, No. 34, 2000 The success in use of the Morse potential for the hydrogen bond may imply that the formation and breaking of a hydrogen bond appears to be a relatively gradual process, as seen from our scaled MP2/6-31G** calculations for amide dimers, and that the hydrogen bond has the covalent character, which is consistent with the recent Compton profile anisotropy experiments on ice Ih. Further tests for the applicability of the Morse potential for hydrogen bonds to polypeptides and proteins are now being carried out. Acknowledgment. This work is supported by grants from KISTEP/MOST, Korea (1999). The author thanks KORDIC for use of the Cray-T3E supercomputer and Mr. Jong Suk Jhon for computational and graphical assistance. References and Notes (1) Pauling, L. The Nature of the Chemical Bond, 3rd ed.; Cornell University Press: Ithaca, NY, 1960; p 498. (2) Hagler, A. T.; Huler, E.; Lifson, S. J. Am. Chem. Soc. 1974, 96, 5319. (3) Momany, F. A.; McGuire, R. F.; Burgess, A. W.; Scheraga, H. A. J. Phys. Chem. 1975, 79, 2361. (4) Brooks, B. R.; Bruccoleri, R. E.; Olafson, B. D.; States, D. J.; Swaminathan, S.; Karplus, M. J. Comput. Chem. 1983, 4, 187. (5) Weiner, S. J.; Kollman, P. A.; Nguyen, D. T.; Case, D. A. J. Comput. Chem. 1986, 7, 230. (6) Hagler, A. T.; Lifson, S.; Dauber, P. J. Am. Chem. Soc. 1979, 101, 5122. (7) Ewig, C. S.; Thacher, T. S.; Hagler, A. T. J. Phys. Chem. B 1999, 103, 6998. (8) Jorgensen, W. L.; Swenson, C. J. J. Am. Chem. Soc. 1985, 107, 569. (9) Allinger, N. L.; Kok, R. A.; Imam, M. R. J. Comput. Chem. 1988, 9, 591. (10) Lii, J.-H.; Allinger, N. L. J. Comput. Chem. 1998, 19, 1001. (11) Cornell, W. D.; Cieplak, P.; Bayly, C. I.; Gould, I. R.; Merz, K. M., Jr.; Ferguson, D. M.; Spellmeyer, D. C.; Fox, T.; Caldwell, J. W.; Kollman, P. A. J. Am. Chem. Soc. 1995, 117, 5179. (12) MacKerell, A. D., Jr.; Bashford, D.; Bellott, M.; Dunbrack, R. L., Jr.; Evanseck, J. D.; Field, M. J.; Fischer, S.; Gao, J.; Guo, H.; Ha, S.; Joseph-McCarthy, D.; Kuchnir, L.; Kuczera, K.; Lau, F. T. K.; Mattos, C.; Michnick, S.; Ngo, T.; Nguyen, D. T.; Prodhom, B.; Reiher, W. E., III.; Roux, B.; Schlenkrich, M.; Smith, J. C.; Stote, R.; Straub, J.; Watanabe, M.; Wio´rkiewicz-Kuczera, J.; Yin, D.; Karplus, M. J. Phys. Chem. B 1998, 102, 3586. (13) No, K. T.; Kwon, O. Y.; Kim, S. Y.; Jhon, M. S.; Scheraga, H. A. J. Phys. Chem. 1995, 99, 3478. (14) Isaacs, E. D.; Shukla, A.; Platzman, P. M.; Hamann, D. R.; Barbiellini, B.; Tulk, C. A. Phys. ReV. Lett. 1999, 82, 600. (15) Dannenberg, J. J.; Haskamp, L.; Masunov, A. J. Phys. Chem. A 1999, 103, 7083. (16) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Gill, P. M. W.; Johnson, B. G.; Robb, M. A.; Cheeseman, J. R.; Keith, T.; Petersson, G. A.; Montgomery, J. A.; Raghavachari, K.; Al-Laham, M. A.; Zakrzewski,

Kang V. G.; Ortiz, J. V.; Foresman, J. B.; Cioslowski, J.; Stefanov, B. B.; Nanayakkara, A.; Challacombe, M.; Peng, C. Y.; Ayala, P. Y.; Chen, W.; Wong, M. W.; Andres, J. L.; Replogle, E. S.; Gomperts, R.; Martin, R. L.; Fox, D. J.; Binkley, J. S.; Defrees, D. J.; Baker, J.; Stewart, J. P.; HeadGordon, M.; Gonzalez, C.; Pople, J. A. Gaussian 94, Revision D.2; Gaussian, Inc.: Pittsburgh, PA, 1995. (17) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Zakrzewski, V. G.; Montgomery, J. A., Jr.; Stratmann, R. E.; Burant, J. C.; Dapprich, S.; Millam, J. M.; Daniels, A. D.; Kudin, K. N.; Strain, M. C.; Farkas, O.; Tomasi, J.; Barone, V.; Cossi, M.; Cammi, R.; Mennucci, B.; Pomelli, C.; Adamo, C.; Clifford, S.; Ochterski, J.; Petersson, G. A.; Ayala, P. Y.; Cui, Q.; Morokuma, K.; Malick, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Cioslowski, J.; Ortiz, J. V.; Baboul, A. G.; Stefanov, B. B.; Liu, G.; Liashenko, A.; Piskorz, P.; Komaromi, I.; Gomperts, R.; Martin, R. L.; Fox, D. J.; Keith, T.; AlLaham, M. A.; Peng, C. Y.; Nanayakkara, A.; Gonzalez, C.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Andres, J. L.; Gonzalez, C.; Head-Gordon, M.; Replogle, E. S.; Pople, J. A.Gaussian 98, Revision A.7; Gaussian, Inc.: Pittsburgh, PA, 1998. (18) Kitano, M.; Fukuyama, T.; Kuchitsu, K. Bull. Chem. Soc. Jpn. 1973, 46, 384. (19) Boys, S. F.; Bernardi, F. Mol. Phys. 1970, 19, 553. (20) Cho, K.-H.; Kang, Y. K.; No, K. T.; Scheraga, H. A. J. Phys. Chem. B, submitted for publication. (21) No, K. T.; Grant, J. A.; Scheraga, H. A. J. Phys. Chem. 1990, 94, 4732. (22) No, K. T.; Grant, J. A.; Jhon, M. S.; Scheraga, H. A. J. Phys. Chem. 1990, 94, 4740. (23) Park, J. M.; No, K. T.; Jhon, M. S.; Scheraga, H. A. J. Comput. Chem. 1993, 14, 1482. (24) Park, J. M.; Kwon, O. Y.; No, K. T.; Jhon, M. S.; Scheraga, H. A. J. Comput. Chem. 1995, 16, 1011. (25) No, K. T.; Kwon, O. Y.; Kim, S. Y.; Cho, K. H.; Yoon, C. N.; Kang, Y. K.; Gibson, K. D.; Jhon, M. S.; Scheraga, H. A. J. Phys. Chem. 1995, 99, 13019. (26) Gay, D. M. ACM Trans. Math. Software 1983, 9, 503. (27) Colominas, C.; Luque, F. J.; Orozco, M. J. Phys. Chem. A 1999, 103, 6200. (28) Kojima, T.; Yano, E.; Nakagawa, K.; Tsunekawa, S. J. Mol. Spectrosc. 1987, 122, 408. (29) Hamzaoui, F.; Baert, F. Acta Crystallogr. 1994, C50, 757. (30) Dixon, D. A.; Dobbs, K. D.; Valentini, J. J. J. Phys. Chem. 1994, 98, 13435. (31) Torii, H.; Tatsumi, T.; Kanazawa, T.; Tasumi, M. J. Phys. Chem. B 1998, 102, 309. (32) Taylor, R.; Kennard, O.; Versichel, W. Acta Crystallogr. 1984, B40, 280. (33) Taylor, R.; Kennard, O.; Versichel, W. J. Am. Chem. Soc. 1983, 105, 5761. (34) Baker, E. N.; Hubbard, R. E. Prog. Biophys. Mol. Biol. 1984, 44, 97. (35) Kang, Y. K.; No, K. T.; Scheraga, H. A. J. Phys. Chem. 1996, 100, 15588. (36) Morse, P. M. Phys. ReV. 1929, 34, 57. (37) Kang, Y. K.; Jhon, M. S. Bull. Korean Chem. Soc. 1981, 2, 8. (38) Scheiner, S.; Duan, X. In Modeling the Hydrogen Bond; Smith, D. A., Ed.; ACS Symposium Series 569; American Chemical Society: Washington, DC, 1994; Chapter 8. (39) Ferguson, D. M.; Kollman, P. A. J. Comput. Chem. 1991, 12, 620.