AMBER Archive (2003)

Subject: Re: SPC/E water model

From: Ioana Cozmuta (ioana_at_nas.nasa.gov)
Date: Wed May 07 2003 - 03:53:43 CDT

Hi Holger,

Thank you so much for your reply, indeed the numbers make sense now. I
did run additional 50ps of NPT dynamics for the SPC/E water but the van
der Waals energy is still positive. The first O-O peak in the RDF function
is located at 2.9A although literature data indicate that it should be at
3.3A. As the average time for water reorientation is about 5ps, I would
assume that my dynamics time is not at all short. The density looks fine
(it is about 1g/cc) and the potential energy -11.2kcal/mol of water
molecule seems ok compared to the value indicated by Berendsen in his
paper (about -9.9 kcal/mol of water molecule) So is this something
specific to the SPC/E model, i.e. to have positive value for the VDW
energy or some of my input parameters are not chosen well? I am attaching
the input file at the end of this e-mail.

Thank you,
Ioana

*************
Running NPT dynamics at 300K
 &cntrl
  imin = 0, irest = 1, ntx = 5,
  ntt = 1, temp0 = 300.0, tautp = 0.2,
  vlimit = 20.0,
  ntp = 1, taup = 0.2, pres0 = 1.0, comp = 44.6, npscal = 1,
  ntb = 2, ntc = 2, ntf = 2, iwrap = 1,
  cut = 6.5, scnb = 2.0, scee = 1.2,
  nstlim = 50000, nscm = 0, t = 0.0, dt = 0.001,
  ntwe = 500, ntwx = 500, ntpr = 100, ntwr = 500,
 &end

> > Hi amber users,
> >
> > I have a box with 216 rigid water molecules (using SPC/E model) and I am
> > equilibrating it at 300K (using amber 7.0, Irix version).
> > The density of the system is fine, the only thing that puzzles me is the
> > positive VDW energy. I did check with the rdparm program the distances
> > between the oxygen atoms in the model (only for O-O the vdw is defined)
> > and the average is 2.5A. For this distance, using the A (629362.17) and C
> > (625.27) and a 6-12 potential I get that the value should be negative. So
> > in principle the sum over all the atoms within the cutoff should be
> > negative as well. But the value in the output file is positive.
> > What am I missing or doing wrong?
> >
> > A second question: According to Berendsen paper (1987 in which the
> > parameters for SPC/E were developed) the LJ parameters for the O-O
> > interaction are A=0.37 (kJ/mol)^(1/6)*nm and B = 0.3428 (kJ/mol)^(1/12)*nm
> > (V_LJ = (B/r)^12-(A/r)^6).
> > In amber the potential is V_amber = A/r^12-C/r^6 thus
> > A_amber=B_berendsen^12=2.074e-5(kcal/mol)^(1/12)*Angs and C_amber=
> > A_berendsen^6=0.0206 (kcal/mol)^(1/6)*Angs
>
> Using V_ber = (B_ber/r)^12-(A_ber/r)^6 and V_amb = A_amb/r^12-C_amb/r^6
> and your values for A_ber and B_ber from above, I get:
>
> A_amb = B_ber^12 = 0.3428^12 * kJ/mol * nm^12
> = 0.3428^12 * kcal/(mol*4.184) * 10^12 Angstrom^12
> = 629358 kcal/mol * Angstrom^12
>
> C_amb = A_ber^6 = 0.37^6 * kJ/mol * nm^6
> = 0.37^6 * kcal/(mol*4.184) * 10^6 Angstrom^6
> = 613 kcal/mol * Angstrom^6
>
> This is close to what you get from rdparm.
>
> > However the rdparm gives values of 629362.17 and 625.27 respectively.
> > Why are these values different?
>
> Using the A_amb and C_amb values and a distance of 2.5A, I get a
> positive energy value.
>
> >
> > Thank you,
> > Ioana
> >
>
> To your second email:
>
> > The printTypes command gives
> > Type r* eps
> > OW 1.7767 0.1553
> >
> > The vdw energy here corresponding to r = 2.7A is a negative number.
> >
> > I thought that A = eps*(r*)^12 =153.65 and
> > C = 2*eps *(r*)^6 = 9.77
>
> To relate r* and eps to A and C, you have to consider that r*(ij) =
> r*(i) + r*(j) in U = eps(r*(ij)/r)^12 - 2eps(r*(ij)/r)^6 (see the "Van
> der Waals parameters for cations in Amber format" page on the AMBER
> homepage).
>
> In your case, r*(ij) = 2 * r* = 2*1.7767 = 3.5534 Angstrom. This then
> gives
> A = eps * r*(ij)^12 = 629362.17 kcal/mol * Angstrom^12 and
> C = 2eps * r*(ij)^6 = 625.27 kcal/mol * Angstrom^6,
> in agreement to what you get from rdparm.
>
> Using r = 2.7 Angstrom, I also get U > 0 here.
>
> Best regards
>
> Holger
>
> --
> +++++++++++++++++++++++++++++++++++++++++++++
> Dr. Holger Gohlke
> Dept. of Molecular Biology, TPC15
> The Scripps Research Institute
> 10550 N. Torrey Pines Rd.
> La Jolla CA 92037 USA
> phone: +1-858-784-9788
> fax: +1-858-784-8896
> email: gohlke_at_scripps.edu
> +++++++++++++++++++++++++++++++++++++++++++++
>
>