From frj@dou.dk Sun May 1 07:21:41 1994 Received: from danpost.uni-c.dk for frj@dou.dk by www.ccl.net (8.6.4/930601.1506) id GAA00630; Sun, 1 May 1994 06:48:04 -0400 Received: from localhost (ean@localhost) by danpost.uni-c.dk (8.6.4/8.6) id MAA24752 for chemistry@ccl.net; Sun, 1 May 1994 12:48:00 +0200 Date: 1 May 94 12:48 +0200 From: Frank Jensen To: chemistry@ccl.net Message-ID: <1126*frj@dou.dk> Subject: Development of new Force Fields Hi netters I have a couple of questions regarding development of new force fields that I would like to have some oppinions on. Any force field contains diagonal terms like Estreth, Ebend etc. In developing a new force field, one first have to choose how many coupling terms to include, then select a given functional form for each term, and finally parameterize against experimental and/or other computational results. The early "Class I" force fields like MM2 and AMBER have very few (or none at all) coupling terms, and the diagonal terms are typically just a second order Taylor expansion for Estretch and Ebend. The new "Class II" force fields that are being developed, contain many coupling terms, and employ functions with more parameters for each term. That way one can avoid special parameters for small rings etc. and also obtain reasonable vibrational frequencies. Now the essense of my question is: why is it that people insist on just extending the Taylor expansion instead of selecting a functional form which is fundamentally more correct? For example, in the MM3 and the CFF93 force field from Biosym, the stretch energy is written as a fourth order Taylor expansion: Estr = a(r-r0)^2 + b(r-r0)^3 + c(r-r0)^4 This is better than just the second order term, but still become qualitatively wrong at long bond lengths as it does not go towards the dissociation energy. Why not simply select say, a Morse potential instead? The Taylor expression contains four parameter (r0 and a,b,c), the Morse only three. In my oppinion, the Morse potential will be as good as the fourth order expansion near r0, and much improved at long bond length (yes I know it takes longer to evaluate an exp function than a simple multiply, but for any reasonably sized molecule the non-bonded energy dominates the computational time anyway). The same holds for Ebend. MM3 and CFF93 just expand the Taylor to 6th or 4th order. This may be fine for say the CCC bend in propane, but fundamentally wrong for the COC bend in dimethyl ether (or any divalent central atom). The latter has a maximum at 180 degrees which the polynomia only accidentally will have. Why not select a functional form which from the start has a maximum at 180 degees? Another question along the same lines: why do people insist that the non-bonded energy parameters _must_ be broken down to atomic parameters? The non-bonded energy depends on two atoms, just as the strecth energy. Why then assign atomic non-bonded parameters but stretch parameters which depends on both atoms? Maybe one could get much better non-bonded energies if accepting the parameters as two atom dependent, and thereby also avoiding the sometimes rather complicated rules for combining the atomic non-bonded parameters to the actual energy expression. Yes, that would give more parameters to fit, but it would be no worse than the stretch energy, and still much less than the torsional parameters. Frank From ross@cgl.ucsf.EDU Sun May 1 17:23:03 1994 Received: from socrates.ucsf.EDU for ross@cgl.ucsf.EDU by www.ccl.net (8.6.4/930601.1506) id RAA03972; Sun, 1 May 1994 17:06:36 -0400 Received: by socrates.ucsf.EDU (8.6.7/GSC4.24) id OAA16154; Sun, 1 May 1994 14:06:34 -0700 Date: Sun, 1 May 1994 14:06:34 -0700 Message-Id: <199405012106.OAA16154@socrates.ucsf.EDU> From: ross@cgl.ucsf.edu (Bill Ross ) To: chemistry@ccl.net Subject: Re: Force Field Nonbon Param Frank Jensen wrote: > Another question along the same lines: why do people > insist that the non-bonded energy parameters _must_ be broken > down to atomic parameters? The non-bonded energy depends on > two atoms, just as the strecth energy. Why then assign atomic > non-bonded parameters but stretch parameters which depends > on both atoms? Maybe one could get much better non-bonded > energies if accepting the parameters as two atom dependent, > and thereby also avoiding the sometimes rather complicated rules > for combining the atomic non-bonded parameters to the actual > energy expression. Yes, that would give more parameters to fit, but > it would be no worse than the stretch energy, and still much > less than the torsional parameters. The "complicated" combining rules allow for a net simplification, since although all atom type pairs are not necessarily bonded, they all must have nonbonded parameters, i.e. instead of vdw terms for on the order of N atom types, one would need terms for N + N(N-1)/2 type pairs. For example, the current release of Amber has 75 atom types, which can be grouped into 19 unique vdw types. Without the combining rule, those 19 vdw types would require 190 pairwise parameters. In my opinion, the real complication is in working with the different combining rules used by different force fields, not each rule per se. However: One reason for taking at least a partial step in the direction of per-pair vdw potentials is that water force fields such as TIP and SPC use inflated oxygen radii to encompass the hydrogens. Such force fields were developed to reproduce the behavior of pure water; they capture O-O and O-H radial distribution functions reasonably well for the computer time required. Now consider the parameterization of a cation in such a situation, where the radial distribution ion-O first peak is made to fit experimental data (e.g. Aqvist, J Phys Chem '90, 94, 8021-8024). The ion nonbonded parameters will allow for the larger oxygen radius, so generally will result in a smaller ion in relation to atoms parameterized in a more general force field such as Amber. The potential minimum for O-O interactions is at 2 * 1.77 A for TIP3 water, while for an Amber oxygen it is 2 * 1.60 A. For example, Cs+ was parameterized to have an ion-O first peak at 3.10 Angstroms in TIP3 water, but in a complex with Amber oxygens, the first peak was at 2.95 A. Thus, even though one may find the use of per-atom nonbonded parameters with combining rules acceptable for a relatively self-consistent force field such as Amber, the consequences of mixing force fields should be considered in the light of the purpose of each simulation attempted. While I have shown how ion parameterization affects relative free energies in water and DNA complex (paper with C.C. Hardin, accepted by JACS), one may wonder if the disparity of force fields affects solute conformations in subtle ways as well. I believe Frank Jensen's proposal for vdw parameterization would minimize the chances for errors on the part of relatively naive users, and is the only way to achieve consistent interatomic distances. GROMOS is one program / force field that has taken this approach (van Gunsteren & Mark, Euro J Biochem '92, 204, 947-961). However, it may be sufficient and preferable to use different combining rules for inter- and intra-force field interactions, to avoid the profusion of parameters from order N to order N^2. Bill Ross From jxh@ibm12.biosym.com Sun May 1 19:21:49 1994 Received: from ibm12.biosym.com for jxh@ibm12.biosym.com by www.ccl.net (8.6.4/930601.1506) id SAA04823; Sun, 1 May 1994 18:48:54 -0400 Received: by ibm12.biosym.com (AIX 3.2/UCB 5.64/4.03) id AA17687; Sun, 1 May 1994 15:48:26 -0700 Date: Sun, 1 May 1994 15:48:26 -0700 From: jxh@ibm12.biosym.com (Joerg Hill) Message-Id: <9405012248.AA17687@ibm12.biosym.com> To: CHEMISTRY@ccl.net Subject: Re: CCL:Development of new Force Fields Frank Jensen wrote: > why is it that people insist on just extending > the Taylor expansion instead of selecting a functional form > which is fundamentally more correct? For example, in the MM3 and > the CFF93 force field from Biosym, the stretch energy is > written as a fourth order Taylor expansion: > > Estr = a(r-r0)^2 + b(r-r0)^3 + c(r-r0)^4 > > This is better than just the second order term, but still become > qualitatively wrong at long bond lengths as it does not go towards > the dissociation energy. Why not simply select say, a Morse potential > instead? The Taylor expression contains four parameter > (r0 and a,b,c), the Morse only three. The functional form above is mathematically equivalent to a Morse potential if you express c as a function of a and b (sorry, I don't have the full expression at hand), so it is not qualitatively wrong at long bond lengths. But there is a more practical problem with Morse potentials in geometry optimizations. The Morse potential has a very shallow slope at long bond lengths, i.e. the forces in this range are very small. This results in problems with the convergence of geometry optimizations if the start geometry is very poor. > The same holds for > Ebend. MM3 and CFF93 just expand the Taylor to 6th or 4th order. > This may be fine for say the CCC bend in propane, but fundamentally > wrong for the COC bend in dimethyl ether (or any divalent central > atom). The latter has a maximum at 180 degrees which the polynomia > only accidentally will have. Why not select a functional form > which from the start has a maximum at 180 degees? You can easily apply boundary conditions to a 4th order Taylor series expansion which give you the desired shape of the potential (b and c are in this case determined from the boundary conditions and only a is fitted). This is done in the Biosym CFF force field. > Another question along the same lines: why do people > insist that the non-bonded energy parameters _must_ be broken > down to atomic parameters? The non-bonded energy depends on > two atoms, just as the strecth energy. Why then assign atomic > non-bonded parameters but stretch parameters which depends > on both atoms? That is simple a question of how to determine the non-bond para- meters. If you would have parameters for each pair of interaction, you need a large number of data to fit them. Unfortunately, experi- mental information is sparse and the use of ab initio calculations requires methods which include electron correlation to get the dispersion (i.e. at least MP2). For the derivation of force fields this is still to demanding for the computers currently available. Maybe the density functional theory can provide some improvement here. An other point is the transferability of the force field. Force fields break down interactions to bonds, angles, torsions etc. and assume that these ones are transferable. Bond lenghts or angles are (at least in organic chemistry) comparatively constant so that this approach works. The non-bond interactions have a much larger range of possible distances (in solids in principle to infinite). So these interactions can not easily be based on pairwise parameters. Besides, non-bond parameters are a real problem in all force fields which still needs a lot of research. Joerg-R. Hill