Dear All, <br><br>I would like to get little feedback of my findings/experience with recent simulations I did with the CHARMM force field in GROMACS and the -chargegroup options with pdb2gmx. I have performed two simulations with a helical peptide (25 AA) in a cubic box filled with explicit TIP3P water. I used the GMX 4.5.3. The first simulation was performed with -chargegroup option "yes" and the second with option "no" (default). <br>
<br>For the first simulation, I used the old version of CHARMM -> gromacs forcefield conversion given with GMX 4.0.5. I have done this MD because i would like to compare some "suspect" results obtained with the same peptide performed 6 months ago with the GMX 4.0.5 version and the old CHARMM->gromacs files available in this distribution. In case of the second simulation, I used the most recent charmm27.ff given in GMX 4.5.3 with a single atom charge groups. <br>
<br>The two simulations were performed during 24 ns with the same protocol (i.e. with same mdp files, below and the same starting files). <br><br>Briefly the protocol i used to equilibrate my system : <br><br>EM (1000 steps with steep) -> NVT at 298 K (with berendsen thermostat and restraints) (100 ps) -> (with nose hoover thermostat and Parinello-Rahman barostat) 400 ps.<br>
<br>Below md.mdp for the two runs ( i am aware that the Coulomb and electrostatic parameters in the mdp are not otimals for the CHARMM ff, but for sake of comparison with others simulations, i didn't change them). <br>
<br>title = KTM17-TIP3 MD<br>constraints = all-bonds<br>integrator = md<br>nsteps = 12000000 ; 24000ps ou 24ns<br>dt = 0.002<br><br>nstlist = 10<br>nstcalcenergy = 10<br>
nstcomm = 10<br><br>continuation = no ; Restarting after NPT<br>vdw-type = cut-off<br>rvdw = 1.0<br>rlist = 0.9<br>coulombtype = PME<br>rcoulomb = 0.9<br>
fourierspacing = 0.12<br>fourier_nx = 0<br>fourier_ny = 0<br>fourier_nz = 0<br>pme_order = 4<br>ewald_rtol = 1e-05<br>optimize_fft = yes<br>
<br>nstvout = 50000<br>nstxout = 50000<br>nstenergy = 20000<br>nstlog = 5000 ; update log file every 2 ps<br>nstxtcout = 1000 ; frequency to write coordinates to xtc trajectory every 2 ps<br>
<br>Tcoupl = nose-hoover:<br>tc-grps = Protein Non-Protein<br>tau-t = 0.4 0.4<br>ref-t = 298 298<br>; Pressure coupling is on<br>Pcoupl = Parrinello-Rahman<br>pcoupltype = isotropic<br>
tau_p = 3.0<br>compressibility = 4.5e-5<br>ref_p = 1.0135<br>gen_vel = no<br><br>I found that the use of the two -chargegrp options influence drastically the MD results especially in my case for the structural properties of the peptides. Indeed, for the first simulation, the peptide keep in its initial helical with a RMSD around 2.0 A along all 24 ns, whereas in the second simulation, the peptide adopt quickly (~10 ns) a unfold conformation with a RMSD value around 5-6.0 A after only 10 ns. Circular dichroism experiments have shown that this peptide adopt an unfold structure in water. So I am sure that the first simulation is not correct (even if the mdp parameters are not optimal) and that my results confirm that chargegrp "yes" option should not be used in MD with CHARMM force field in GROMACS, as it was discussed recently in this mailing list (for example <a href="http://lists.gromacs.org/pipermail/gmx-users/2010-September/054106.html">http://lists.gromacs.org/pipermail/gmx-users/2010-September/054106.html</a>). <br>
<br>SA<br>