I have general question. If this parameters are not giving the proper results then How come I get the proper results with same parameter if I use 1/2 value of Kbond with the same parameters.
I checked the temp when I use half value of Kbond and it goes max. ~320. Nilesh On Wed, May 11, 2011 4:03 pm, Justin A. Lemkul wrote: > > > 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 didnt 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 dont 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 > > > -- 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
