Cenfeng Fu wrote:
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?
The original derivation of OPLS benzene parameters (in the paper you cite)
applied dispersion correction to the energy term, so I would suspect that you
should do the same. No constraints are mentioned in that paper, so it is hard
to say whether or not you should be using them.
You can easily see the effects of constraints by running a few more simulations
(no constraints, then hbonds only, etc).
Another challenge pertains to the use of cutoffs. The paper you cite used a
1.3-nm ring-ring cutoff. I don't know how you would implement that in Gromacs,
but in your original post
(http://lists.gromacs.org/pipermail/gmx-users/2010-October/054499.html) you were
using 1.0-nm cutoffs. Differences in the way nonbonded interactions are
calculated can (and quite often, will) make a difference in your results.
-Justin
Best regards!
Cenfeng Fu
--
========================================
Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
========================================
--
gmx-users mailing list [email protected]
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to [email protected].
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
