FW: [gmx-users] charmm and oplsaa simulation results don't match
sarbani chattopadhyay
sarbani_c84 at rediffmail.com
Tue Aug 19 18:43:42 CEST 2008
hi,
Thanks for the suggestions. Actually the microstates are revisited in the 100ns span and one of them is occupied for longer time span than the others.
I had broken the 100ns trajectory(oplsaa) in 3 parts and took starting structures from each of them and run simulations and obtained compatible results. However what's disturbing is that in the charmm run, the configuration of the staring structures changes within the first few picoseconds and then remains in the new configuration for the rest of the simulation time.It means that the starting structure changes into another structure and retains that for the rest of the simulation and that happens to be the most occupied state in the oplsaa run.
I may be wrong, but I think somewhere i am having problem in the way charmm and tip3p water model works ie. it may be that the L-J potentials of the hydrogen atoms are not treated properly.
Surprisingly, I had run charmm with spc water model(by mistake), and although the most occupied state remains the same but the microstates are obtained.
Thanks in advance
Sarbani
On Tue, 19 Aug 2008 Berk Hess wrote :
>
>
>
>
> From: gmx3 at hotmail.com
>To: sarbani_c84 at rediffmail.com
>Subject: RE: [gmx-users] charmm and oplsaa simulation results don't match
>Date: Tue, 19 Aug 2008 10:47:02 +0200
>
>
>
>
>
>
>
>Hi,
>
>There are two main issues here.
>
>One is sampling. I don't know what the correlations times in your systems are.
>You say you obtain a few microstates. But does the system revisit the same
>microstates over the 100 ns? If not, your simulations are certainly too short
>to have full sampling and you can not expect that two simulations with the same
>or different force fields produce the same sampling.
>
>The second issue is the force field. Different force fields can favor different
>configurations. This is especially true for backbone dihedrals. This will have
>stronger effects on small peptides than on larger (less flexible) proteins.
>
>And finally, it seems you are running with a plain Coulomb cut-off.
>This will produce artefacts. You should use PME (or at least a reaction-field).
>
>Berk
>
>Date: Tue, 19 Aug 2008 08:28:25 +0000
> From: sarbani_c84 at rediffmail.com
>To: gmx-users at gromacs.org
>Subject: [gmx-users] charmm and oplsaa simulation results don't match
>
>
>
>
>Hi everybody,
>
> I had run a long 100ns simulation on a tripeptide using opls/aa force field and
>
>spc water model.
>
>I got a few microstates.
>
>Then I tried to run simulation on the same system using charmm force field and tip3p water
>
>model.
>
>However I couldn't find the microstates. Within the first 10 ns the system adopts the final
>
> configuration that was observed in the simulation using oplsaa force field.
>
>
>
>
>
> The process I used to run gromacs with charmm and tip3p water model is
>
>
>
> pdb2gmx -f opls.gro -ff charmm -ter -water tip3p
>
>
>
> charmm27-run-addtop-gromacs-dihe.sh topol.top new.top
>
>
>
>
>
> mv new.top topol.top
>
>
>
>
>
> #### I did modify the top file to include "tip3p-charmm.itp"
>
> for water topology.
>
>
>
>
>
> editconf -f conf.gro -o box.gro -bt cubic -d 0.75
>
>
>
> genbox -cs spc216.gro -cp box.gro -p topol.top -o solvated.gro
>
>
>
> grompp -zero -f em.mdp -c solvated.gro -o em.tpr
>
>
>
>
>
> mdrun -v -s em.tpr
>
>
>
>
>
> Then I generated the "tpr" file for 100ns.
>
>
>
> The em.mdp file includes :
>
>
>
>
>
>
>
>; User spoel (236)
>
>; Wed Nov 3 17:12:44 1993
>
>; Input file
>
>;
>
>cpp = /usr/bin/cpp
>
>constraints = none
>
>integrator = steep
>
>nsteps = 100
>
>;
>
>; Energy minimizing stuff
>
>;
>
>emtol = 2000
>
>emstep = 0.01
>
>
>
>
>
>
>
>nstcomm = 1
>
>ns_type = grid
>
>rlist = 1
>
>rcoulomb = 1.0
>
>rvdw = 1.0
>
>Tcoupl = no
>
>Pcoupl = no
>
>gen_vel = no
>
>
>
>
>
>The 100ns.mdp file includes:
>
>
>
>
>
>;
>
>; User spoel (236)
>
>; Wed Nov 3 17:12:44 1993
>
>; Input file
>
>;
>
>title = Yo
>
>cpp = /usr/bin/cpp
>
>constraints = all-bonds
>
>integrator = md
>
>dt = 0.002 ; ps !
>
>nsteps = 50000000 ; total 10000 ps.
>
>nstcomm = 1
>
>nstxout = 250
>
>nstvout = 1000
>
>nstfout = 0
>
>nstlog = 250
>
>nstenergy = 250
>
>nstlist = 10
>
>ns_type = grid
>
>rlist = 1.0
>
>rcoulomb = 1.0
>
>rvdw = 1.0
>
>; Berendsen temperature coupling is on in two groups
>
>Tcoupl = berendsen
>
>tc-grps = Protein SOL
>
>tau_t = 0.1 0.1
>
>ref_t = 300 300
>
>; Energy monitoring
>
>energygrps = Protein SOL
>
>; Isotropic pressure coupling is now on
>
>Pcoupl = berendsen
>
>Pcoupltype = isotropic
>
>tau_p = 0.5
>
>compressibility = 4.5e-5
>
>ref_p = 1.0
>
>; Generate velocites is off at 300 K.
>
>gen_vel = no
>
>gen_temp = 300.0
>
>gen_seed = 173529
>
>
>
>
>
>Can anyone suggest why am I not getting the microstates with "charmm " force field.
>
>Is there anything more that I needed to do and I missed. I am sorry for such a lengthy query
>
>
>
>
>
>Thanks in advance
>
>
>
>Sarbani
>
>
>
>
>
>
>
>
>Express yourself instantly with MSN Messenger! MSN Messenger
>
>_________________________________________________________________
>Express yourself instantly with MSN Messenger! Download today it's FREE!
>http://messenger.msn.click-url.com/go/onm00200471ave/direct/01/
>_______________________________________________
>gmx-users mailing list gmx-users at gromacs.org
>http://www.gromacs.org/mailman/listinfo/gmx-users
>Please search the archive at http://www.gromacs.org/search before posting!
>Please don't post (un)subscribe requests to the list. Use the
>www interface or send it to gmx-users-request at gromacs.org.
>Can't post? Read http://www.gromacs.org/mailing_lists/users.php
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20080819/00e8dcf8/attachment.html>
More information about the gromacs.org_gmx-users
mailing list