Nilesh Dhumal wrote:
Hello Justin,
I used 0.5 fs time step and still I got dielectric constant ~2.
This is the md.mdp file I used. I checked the temperature it didn't go
more than ~ 320 K.
Still, that's unacceptable if your ref_t is 295. The fact that the temperature
is still rocketing up suggests the same instability I found.
title = cpeptide MD
cpp = /usr/bin/cpp
integrator = md
dt = 0.001 ; ps !
nsteps = 6500000 ; total 5 ps.
nstcomm = 1000
nstxout = 1000
nstvout = 1000
nstfout = 1
nstlist = 1
ns_type = grid
rlist = 0.6
rcoulomb = 0.6
rvdw = 0.7
These cutoffs make no sense. The paper you've cited used 0.9 nm. If you're
trying to reproduce previous work, use the same settings.
coulombtype = PME
vdwtype = cut-off
Similarly here, you haven't applied dispersion correction, which you should, per
the methods used in the JCP article.
-Justin
pbc = xyz
fourierspacing = 0.12
fourier_nx = 0
fourier_ny = 0
fourier_nz = 0
pme_order = 4
ewald_rtol = 1e-5
optimize_fft = yes
; Berendsen temperature coupling is on
Tcoupl = v-rescale
tau_t = 0.1
tc-grps =system
ref_t = 295
; Pressure coupling is on
Pcoupl = parrinello-rahman
pcoupltype = isotropic
tau_p = 0.5
compressibility = 4.5e-5
ref_p = 1.0
; Generate velocites is on at 300 K.
gen_vel = yes
gen_temp = 295.0
gen_seed = 173529
Nilesh
On Tue, May 10, 2011 6:14 pm, Justin A. Lemkul wrote:
The gmx-developers list is not the forum for these types of questions. I
am replying via gmx-users, which is where the discussion should stay.
I took a few minutes to dig into this. My conclusion is that your system
is not stable. I would encourage you to analyze the temperature and
pressure of your systems that are giving odd results. When I used the
topology you provided, the temperature spiked to over 600 K and the box
began to oscillate extensively as if the system were about to explode.
The resulting epsilon value was about 1.6.
If I reduce the timestep to 0.5 fs, the results are much more reasonable,
and using the Ka and Kb values in the paper (not halved!) I get a much
more sensible result, even after just 100 ps. It looks to me as if the
combination of these parameters and a 1-fs timestep is not entirely
stable. I know the original authors used dt = 1 fs, but the point
remains.
-Justin
Nilesh Dhumal wrote:
Hello,
I am using flexible water model for my system. I am referring a paper
J.
Chem. Phys. 124, 024503 (2006). Author have used Amber type force field.
i.e. 1/2 Kbond (r-req)**2 + 1/2 Kangle(theta-thetaeq)**2.
Kbond = 443153.3808 kJ/mol nm**2
Kangle= 317.5656 kJ/mol rad**2.
I am using olss-aa force field parameters in Gromacs VERSION 4.0.7.
I checked some papers in which author have used oplsaa force field in
Gromacs. 1/2 factor is not in opls force field if I compare opls and
amber.
I didn’t get the proper dielectric constant for water when I used the
parameters reported in paper
(Kbond = 443153.3808 kJ/mol nm**2,Kangle= 317.5656 kJ/mol rad**2)
I half the value of Kb and get the proper dielectric constant (~80) for
water reported in paper. If I half Kangle then I don’t get proper
value. Below are the results for the dielectric constant of water.
Bond length is nm.
Here I have done some analysis. The original value reported in J.
Chem.
Phys. 124, 024503 2006, paper are
Kbond = 443153.3808 kJ/mol nm**2
Kangle = 317.5656 kJ/mol rad**2.
bond length Kbond angle Kangle dielectric constant
0.1012
443153.3808 113.24 317.5656 ~1.9 : orginal value
0.1012 221576.6904 113.24 317.5656 ~80 : 1/2 (Kbond)
0.1012 443153.3808 113.24 158.7828 ~1.58 : 1/2 (kangle)
0.1012 221576.6904 113.24 317.5656 ~1.9 : 1/2
(Kbond)&(Kangle)
Here I pasted spc_fw.itp file used for the simulation.
; Derived from parsing of runfiles/alat.top.orig
[ defaults ]
; nbfunc comb-rule gen-pairs fudgeLJ fudgeQQ
;1 3 yes 0.5 0.5
; comb-rule 3 is square-root sigma, the OPLSAA version
[ atomtypes ]
; full atom descriptions are available in ffoplsaa.atp
; name bond_type mass charge ptype sigma epsilon
opls_1001 OW 8 15.99940 -0.820 A 3.1655e-01
6.503e-01
opls_1002 HW 1 1.00800 0.410 A 0.00e+00 0.00e+00
[ bondtypes ]
; i j func b0 kb
OW HW 1 0.1012 443153.3808 ; J. Chem. Phys.
(2006),124,024503
[ angletypes ]
; i j k func th0 cth
HW OW HW 1 113.24 317.5656 ; J. Chem. Phys.
(2006),124,024503
[ moleculetype ]
; Name nrexcl
WAT 3
[ atoms ]
; nr type resnr residue atom cgnr charge mass typeB
chargeB massB 1 opls_1001 1 WAT OW 1 -0.82
15.99940 ;
2 opls_1002 1 WAT HW1 1 0.41 1.008 ;
3 opls_1002 1 WAT HW2 1 0.41 1.008 ;
[ bonds ]
; i j funct
1 2 1
1 3 1
[ angles ]
; i j k funct
2 1 3 1
In water.top file, I included spc_fw.itp file.
; Include water topology
#include "spc_fw.itp"
I run the simulation 6.5 ns for collecting data and I have total 256
water molecules.
Is there anything wrong with my .itp file?
I am using Gromacs VERSION 4.0.7.
Thanks
Nilesh
--
========================================
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
--
========================================
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
