[gmx-users] Simulated Annealing Protocol...

rama david ramadavidgroup at gmail.com
Sat Mar 17 15:16:04 CET 2012


Hi Gromacs Specialist,
               Thank you Justin For your reply.
I run one Simulated Annealing) MD after your reply.
My aim is to find protein conformation and study protein
 self assembly study ..

the protocol I follow is as..
1.  Upto energy minimisation I followed the general way
2.  I did nvt at 277  K(position resstrain)
3.  I did npt at 277   k  ((Position Restrain  )
 4. then I did the SA  (I removed Position Restrain ) with following mdp

title                   = simulated run

integrator       = md                    ; md integrator
tinit                = 0                ; [ps] starting time for run
dt                  = 0.002          ; [ps] time step for integration
nsteps           = 5000000      ; maximum number of steps to
integrate, 0.002 * 5000000 10 ns
comm_mode   = Linear         ; remove center of mass translation
comm_grps     = Protein Non-Protein   ; group(s) for center of mass
motion removal

; 7.3.8 Output Control
nstxout                 = 100       ; [steps] freq to write
coordinates to trajectory
nstvout                 = 100       ; [steps] freq to write velocities
to trajectory
nstfout                 = 100       ; [steps] freq to write forces to trajectory
nstlog                  = 100           ; [steps] freq to write
energies to log file
nstenergy               = 100           ; [steps] freq to write
energies to energy file
nstxtcout               = 500           ; [steps] freq to write
coordinates to xtc trajectory
xtc_precision           = 100          ; [real] precision to write xtc
trajectory
xtc_grps                = System        ; group(s) to write to xtc trajectory
energygrps              = System        ; group(s) to write to energy file

; 7.3.9 Neighbor Searching
nstlist                 = 5             ; [steps] freq to update neighbor list
ns_type                 = grid          ; method of updating neighbor list
pbc                     = xyz           ; periodic boundary conditions
in all directions
rlist                   = 1.0           ; [nm] cut-off distance for
the short-range neighbor list

; 7.3.10 Electrostatics
coulombtype             = PME           ; Particle-Mesh Ewald electrostatics
rcoulomb                = 1.0           ; [nm] distance for Coulomb cut-off

; 7.3.11 VdW
vdwtype                 = cut-off       ; twin-range cut-off with
rlist where rvdw >= rlist
rvdw                    = 1.4           ; [nm] distance for LJ cut-off
DispCorr                = EnerPres      ; apply long range dispersion
corrections for energy

; 7.3.13 Ewald
fourierspacing          = 0.16          ; [nm] grid spacing for FFT
grid when using PME
pme_order               = 4             ; interpolation order for PME, 4 = cubic
ewald_rtol              = 1e-5          ; relative strength of
Ewald-shifted potential at rcoulomb

; 7.3.14 Temperature Coupling
tcoupl                  = berendson                   ;
tc_grps                 = Protein    Non-Protein        ; groups to
couple seperately to temperature bath
tau_t                   = 0.1        0.1                ; [ps] time
constant for coupling
ref_t                   = 277        277               ; [K] reference
temperature for coupling
annealing               = single single
annealing_npoint        = 2  2
annealing_time          = 0 50 0 50
annealing_temp          = 277 333 277 300
; 7.3.15 Pressure Coupling
pcoupl                  = Parrinello-Rahman      ;
pcoupltype              = isotropic                    ; pressure
coupling in x-y-z directions
tau_p                   = 2.0                              ;  [ps]
time constant for coupling
compressibility         = 4.5e-5                       ; [bar^-1]
compressibility
ref_p                   = 1.0                                ; [bar]
reference pressure for coupling

; 7.3.17 Velocity Generation
gen_vel                 = no                   ; velocity generation turned off
gen_temp              = 277
; 7.3.18 Bonds
constraints             = all-bonds     ; convert all bonds to constraints
constraint_algorithm    = LINCS         ; LINear Constraint Solver
continuation            = yes           ; apply constraints to the
start configuration
lincs_order             = 4             ; highest order in the
expansion of the contraint coupling matrix
lincs_iter              = 1             ; number of iterations to
correct for rotational lengthening
lincs_warnangle         = 30            ; [degrees] maximum angle that
a bond can rotate before LINCS will complain
continuation              =yes

 so my Queries are
1.    As I not  define =  -DPOSRES in SA mdp  is it sound well???,
should I have to run Production run after my
       SA run(That means  Is SA is substitute to Production run ??)
??? or I have to use   define =  -DPOSRES in my SA mdp ,
       and then I have to do production run
 2.   Or as per Justin recommendation  NVT at 277 ,followed by SA
       (Then I have to use  define = -DPOSRES  in My SA mdp file ,Is it right??)
      then NPT   afterward Production run ..

 So what is your suggestion ???
    All suggestion are appreciated...
Thank you

With Best Regards,

On Fri, Mar 16, 2012 at 7:10 PM, Justin A. Lemkul <jalemkul at vt.edu> wrote:
>
>
> rama david wrote:
>>
>> Dear Gromacs Specialists,
>>                I am very novice to Molecular Simulation study.
>>  I am using GROMACS 4.5.4 version .
>> I completed some GROMACS tutorials , I not found any
>> tutorial on Simulated Annealing..
>> If Any one know the link please give me it..
>>
>> I make my  protocol to work on simulated annealing as follow ...
>> (I am not writing in detail sorry for that )
>> 1.  pdb2gmx ...
>> 2. editconf
>> 3. Solvent Addition
>> 4.  Ion addition
>> 5. Energy minimisation
>> 6. simulated annealing
>>  mdp for simulated annealing is as follow...
>>
>> ; 7.3.3 Run Control
>> title                   = simulated run
>>
>> integrator              = md                    ; md integrator
>> tinit                      = 0                     ; [ps] starting time
>> for run
>> dt                         = 0.002                 ; [ps] time step
>> for integration
>> nsteps                  = 5000000               ; maximum number of
>> steps to integrate, 0.002 * 5000000 10 ns
>> comm_mode         = Linear                ; remove center of mass
>> translation
>> comm_grps           = Protein Non-Protein   ; group(s) for center of
>> mass motion removal
>>
>> ; 7.3.8 Output Control
>> nstxout                 = 100       ; [steps] freq to write
>> coordinates to trajectory
>> nstvout                 = 100       ; [steps] freq to write velocities
>> to trajectory
>> nstfout                 = 100       ; [steps] freq to write forces to
>> trajectory
>> nstlog                  = 100           ; [steps] freq to write
>> energies to log file
>> nstenergy               = 100           ; [steps] freq to write
>> energies to energy file
>> nstxtcout               = 500           ; [steps] freq to write
>> coordinates to xtc trajectory
>> xtc_precision           = 500          ; [real] precision to write xtc
>> trajectory
>> xtc_grps                = System        ; group(s) to write to xtc
>> trajectory
>> energygrps              = System        ; group(s) to write to energy file
>>
>> ; 7.3.9 Neighbor Searching
>> nstlist                 = 5             ; [steps] freq to update neighbor
>> list
>> ns_type                 = grid          ; method of updating neighbor list
>> pbc                     = xyz           ; periodic boundary conditions
>> in all directions
>> rlist                   = 1.0           ; [nm] cut-off distance for
>> the short-range neighbor list
>>
>> ; 7.3.10 Electrostatics
>> coulombtype             = PME           ; Particle-Mesh Ewald
>> electrostatics
>> rcoulomb                = 1.0           ; [nm] distance for Coulomb
>> cut-off
>>
>> ; 7.3.11 VdW
>> vdwtype                 = cut-off       ; twin-range cut-off with
>> rlist where rvdw >= rlist
>> rvdw                    = 1.4           ; [nm] distance for LJ cut-off
>> DispCorr                = EnerPres      ; apply long range dispersion
>> corrections for energy
>>
>> ; 7.3.13 Ewald
>> fourierspacing          = 0.16          ; [nm] grid spacing for FFT
>> grid when using PME
>> pme_order               = 4             ; interpolation order for PME, 4 =
>> cubic
>> ewald_rtol              = 1e-5          ; relative strength of
>> Ewald-shifted potential at rcoulomb
>>
>> ; 7.3.14 Temperature Coupling
>> tcoupl                  = berendson                   ; temperature
>> coupling with Nose-Hoover ensemble
>> tc_grps                 = Protein    Non-Protein        ; groups to
>> couple seperately to temperature bath
>> tau_t                   = 0.1        0.1                ; [ps] time
>> constant for coupling
>> ref_t                   = 5           5               ; [K] reference
>> temperature for coupling
>> annealing               = single single
>> annealing_npoint        = 2  2
>> annealing_time          = 0 20 0 20
>> annealing_temp          = 5 333 5 333
>> ; 7.3.15 Pressure Coupling
>> pcoupl                  =  parrinello-Rahman    ; pressure coupling
>> where box vectors are variable
>> pcoupltype              = isotropic             ; pressure coupling in
>> x-y-z directions
>> tau_p                   = 2.0                   ; [ps] time constant
>> for coupling
>> compressibility         = 4.5e-5                ; [bar^-1] compressibility
>> ref_p                   = 1.0                   ; [bar] reference
>> pressure for coupling
>>
>> ; 7.3.17 Velocity Generation
>> gen_vel                 = yes           ; velocity generation turned off
>> gen_temp                = 5
>> ; 7.3.18 Bonds
>> constraints                  = all-bonds     ; convert all bonds to
>> constraints
>> constraint_algorithm    = LINCS         ; LINear Constraint Solver
>> continuation               = yes           ; apply constraints to the
>> start configuration
>> lincs_order                 = 4             ; highest order in the
>> expansion of the contraint coupling matrix
>> lincs_iter                    = 1             ; number of iterations
>> to correct for rotational lengthening
>> lincs_warnangle         = 30            ; [degrees] maximum angle that
>> a bond can rotate before LINCS will complain
>>
>> My Queries
>>
>> 1. Is My mdp file ok ???..please give me a nice protocol..
>
>
> It is a bad idea to generate velocities and use P-R pressure coupling at the
> same time.  Such an approach is extremely fragile and your simulation is
> likely to crash.
>
> I assume that since you have specified heating over 20 ps of the 10 ns
> trajectory that you wish to have the temperature maintained from then on.  I
> believe Gromacs will handle this elegantly.  Pay attention to everything
> grompp tells you, in any case.
>
>
>> 2.. Should I have to do position restrained MD before SA(simulated
>> annealing)
>>    (If yes then what the temp. should I have to used for NVT and NPT
>> (as in mdp file has lower 5 K and high 300 k  ) )
>
>
> The answer to this question depends on what you're trying to accomplish with
> SA.  Position restraints are used to disfavor major structural changes in
> the solute as the solvent is equilibrated around it.  If SA is your initial
> step in equilibration, I would think position restraints can't hurt.
>
> My general protocol is NVT at the initial temperature, SA to heat to the
> desired temperature, NPT at the final conditions, and then data collection.
>  You should base your decision upon what you are wishing to observe and by
> interpreting how others have done similar things.
>
> -Justin
>
> --
> ========================================
>
> 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    gmx-users at gromacs.org
> 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 gmx-users-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/Support/Mailing_Lists



More information about the gromacs.org_gmx-users mailing list