[gmx-users] Re: the density of liquid benzene is too large

[Cenfeng Fu]

> Several possible reasons:
>
> 1. The force field parameters aren't perfect, so there is some inherent
> disagreement between simulation and reality.  What is the expected value
> for
> this force field?


> 2. You're using the isothermal compressibility for water.  If your system
> is
> pure liquid benzene, I'd think you would want to use the compressibility
> for
> benzene.  I don't know how big the difference would be off-hand, but at
> least
> you'd be convinced that your simulation was set up properly.
>
> 3. Although it won't matter a huge amount, to what temperature does the
> experimental density correspond?  Usually these parameters are given at 25C
> (298
> K).  The temperature you've used is 300 K.  Again, a minor point, but one
> worth
> doing correctly in a simulation even though the such a change would not
> account
> for the magnitude of difference you're currently seeing.
>
> -Justin
>
> Hi Justin,
  Thanks for your suggestions!
  The OPLSAA model for liquid benzene could get a density of 0.873+/-0.001
g/cm^3 (JACS,1990,112,4768) at 298K and 1 atm. I have done some test and
maybe I have found what is the problem. In the previous simulation, I
applied long range dispersion corrections for energy and pressure with
"DispCorr = EnerPres". When I apply the long range dispersion corrections
only for energy with "DispCorr = Ener" (the temperature is set to be 298K
and the pressure is set to be 1 bar), the density of the system is 0.883,
which is closely to the experimental value and the expect value of this
model.
  Now, I have anther question. After the simulation, I want to calculate the
hear capacity of liquid benzene. So I using this command:
g_energy -f *.edr -s *.tpr -o energy.xvg -b 10000 -nmol 600 -nconstr 12
 And these are the results:
Statistics over 5000001 steps [ 10000.0000 through 20000.0000 ps ], 11 data
sets
All statistics are over 500001 points

Energy                      Average   Err.Est.       RMSD  Tot-Drift
-------------------------------------------------------------------------------
Potential                   19.7926      0.036   0.435187 -0.0870648
(kJ/mol)
Kinetic En.                 29.7264    1.7e-05   0.360236 -7.68964e-05
(kJ/mol)
Total Energy                 49.519      0.036   0.575006 -0.0871418
(kJ/mol)
Temperature                 297.999    0.00017    3.61127 -0.000770665  (K)
Pressure                    1.15978     0.0036    171.048 -0.00963121  (bar)
Box-X                       4.45033     0.0015  0.0105071 -0.0038387  (nm)
Box-Y                       4.45033     0.0015  0.0105071 -0.0038387  (nm)
Box-Z                       4.45033     0.0015  0.0105071 -0.0038387  (nm)
Volume                      88.1423      0.089   0.624413  -0.228234  (nm^3)
Density                     883.013       0.89    6.25236    2.28197
(kg/m^3)
Enthalpy                    29712.2         21    345.005   -52.2857
(kJ/mol)

Temperature dependent fluctuation properties at T = 297.999. #constr/mol =
12
Isothermal Compressibility: 0.000107512 /bar
Adiabatic bulk modulus:        9301.25  bar
Heat capacity at constant pressure Cp:    218.791 J/mol K
Thermal expansion coefficient alphaP: 0.000135136 1/K

I got a hear capacity at constant pressure with 218.791 J/(mol K). However,
the experimental value is 135.98 J/(mol K), and the expect value of the
model is 130.54 J/(mol K) (JACS,1990,112,4768). In the new simulation, I
used LINCS for all-bonds. So I think the #nconstr should be 12 in the
g_energy command. Is this value for #nconstr right? If it is wrong, what
value should I use. Or I should not use long range dispersion correction for
energy? Are there other mistakes with my parameters?

Best regards!
Cenfeng Fu
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20101014/23b0c825/attachment.html>