What you should do depends on your question. If you question is represented by 
your original post and by your subject (what is the effect of 
refcoord_scaling?) then you need to repeat both conditions a few times. If you 
just want to do it correctly, then use refcoord_scaling=com (note that this is 
not a guarantee of success if you have setup your system incorrectly or if your 
parameters lead to an inherently unstable association).

If your system is unstable at 300 K, then I do not agree that this makes it 
acceptable to run at 200 K. This may just be hiding some underlying problem 
with your setup. For example, if you simulate at 200 K for a while and then 
move to 300 K, does that lead to a stable system? If so, then the problem was 
with your initial conformation. What has been done in the literature?

You should use a single temperature coupling group. I use the sd integrator 
(velocity Langevin dynamics) but I believe that you could also use velocity 
rescaling as your temperature coupling algorithm (I have not used it myself).

Usually with all bonds constrained by LINCS one uses a 2 fs timestep. I presume 
that you use a 0.5 fs timestep because that makes your conformation more 
stable? Again, this may just be hiding some problem with initial conformations 
and you may be able to move to a 2 fs timestep after some equilibration. As to 
the nstlist setting, it is common to update the neighbour list every 20 fs, but 
that will be dependent on your forcefield and your temperature. What have other 
people done with systems like this?

This will probably be my last comment on this thread. It is not usually 
possible to move from problems in system setup to a stable system with usable 
results within a week. I think that you would be best to aim for tracking down 
the problem and having a reliable simulation system by the time of your 
presentation.

Chris.

-- original message --

Matthias Ernst Matthias.Ernst2 at student.kit.edu 
Thu Sep 6 17:55:43 CEST 2012
Previous message: [gmx-users] Effect of refcoord_scaling
Messages sorted by: [ date ] [ thread ] [ subject ] [ author ]
> The problem could be gen_vel=no if you are not loading velocities that are 
> consistent with your production temperature of 200 K. If, for example, you do 
> not load any velocities at all, then the initial forces will quickly be 
> scaled up to reach 200 K and this can cause large scale 
> deformations/rearrangements.
I have generated velocities the step before (in NpT equilibration) and 
for the same temperature (200K). So, I would not think this causes the 
error as the should be transferred by using the resulting gro-file from 
NpT as input for production MD run.
> There are a number of other strange things (like nstlist=5 with dt=0.0005 -- 
> nstlist seems unnecessarily frequent, the 200 K temperature, and your use of 
> separate temperature coupling groups which does not give the correct 
> ensemble). Still, those are tangents to your main question.
If I repeat it, I should do so with "correct" parameters, I think. Or do 
you suggest to repeat it with the very same parameter set as before?
To address the things you mentioned:
- What should or can nstlist be set to?
- 200K have been chosen because at higher temperatures, the test systems 
used to move even more.
- Concerning the coupling groups: you mean, I should use a group for 
protein and DNA as whole? But isn't restraining the DNA also a 
perturbation of the ensemble? At least, I thought that this way the 
temperature of the protein should not be affected by the restraints on 
the DNA. Tthis was the reason why I put several coupling groups and why 
I asked yesterday if it would be better not to couple the restrained DNA 
to the thermostate or set its temperature to 0.
> Nevertheless, if you don't have time to run any more simulations (as you 
> stated earlier) then none of this is really going to help you. If, on the 
> other hand, you do have time to run other simulations, then I still think you 
> should start with repeating your result.
That's true. I fear, I will not be able to get it finished till next 
week but at least, I can try, start some jobs and see how far I can get.

Thank you very much for your help,
Matthias

> Chris.
>
>   -- original message --
>
> I guess you want to see the production mdp file. Here it is:
>
>
> ; VARIOUS PREPROCESSING OPTIONS
> title                    = Production Simulation
> cpp                      = /lib/cpp
> define                   = -DPOSRES_DNA
>
> ; RUN CONTROL PARAMETERS
> integrator               = md
> tinit                    = 0       ; Starting time
> dt                       = 0.0005 ; 2 femtosecond time step for integration
> nsteps                   = 100000000
>
> ; OUTPUT CONTROL OPTIONS
> nstxout                  = 50000 ; Writing full precision coordinates
> every nanosecond
> nstvout                  = 50000 ; Writing velocities every nanosecond
> nstfout                  = 0     ; Not writing forces
> nstlog                   = 5000  ; Writing to the log file every step
> nstenergy                = 5000  ; Writing out energy information every step
> nstxtcout                = 5000  ; Writing coordinates every 5 ps
> energygrps               = Protein DNA water_and_ions
>
> ; NEIGHBORSEARCHING PARAMETERS
> nstlist                  = 5
> ns-type                  = Grid
> pbc                      = xyz
> rlist                    = 0.9
>
> ; OPTIONS FOR ELECTROSTATICS AND VDW
> coulombtype              = Reaction-Field
> rcoulomb                 = 1.4
> epsilon_rf               = 78
> epsilon_r                = 1
> vdw-type                 = Cut-off
> rvdw                     = 1.4
>
> ; Temperature coupling
> Tcoupl                   = Berendsen
> tc-grps                  = Protein  DNA   water_and_ions
> tau_t                    = 0.1      0.1      0.1
> ref_t                    = 200      200   200
> ; Pressure coupling
> Pcoupl                   = Berendsen
> Pcoupltype               = Isotropic
> tau_p                    = 1.0
> compressibility          = 4.5e-5
> ref_p                    = 1.0
>
> ; GENERATE VELOCITIES FOR STARTUP RUN
> gen_vel                  = no
>
> ; OPTIONS FOR BONDS
> constraints              = all-bonds
> constraint-algorithm     = Lincs
> unconstrained-start      = yes
> lincs-order              = 4
> lincs-iter               = 1
> lincs-warnangle          = 30
>
>
> Thanks for your help,
> Matthias
--
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

Reply via email to