[gmx-users] Flat-bottomed position restraint set

Oliver Stueker ostueker at gmail.com
Sun Nov 16 18:16:44 CET 2014


Dear Liuyoung,


in your equil.mdp
​ you have ​
:

> ​define          = -D
> ​​
> *​​POSRES​​​*

​​

​and in your topol.top:
​


> #ifdef
> *​​POSRES_WATER*
> ; Position restraint for each water oxygen
> ;[ position_restraints ]
> ;  i funct       fcx        fcy        fcz
> ;  1    1       1000       1000       1000
> [ position_restraints ]
> ;  i funct       g         r(nm)       k
>    1    2        1          3.0       30.0
> ;   2    2        1          3.0       30.0
> ;   3    2        1          3.0       30.0
> #endif


​That means your position_restraints directive will not be included​.

Use whatever word you like for both the define and ifdef line, as long as
it's the same.
e.g. define -DPOSRES  and  #ifdef POSRES
​or define -DPOSRES_FLATBOTTOM and ​#ifdef POSRES_FLATBOTTOM


​Best,
Oliver​



On Sat, Nov 15, 2014 at 11:11 PM, liuyong_1007 at dicp.ac.cn <
liuyong_1007 at dicp.ac.cn> wrote:

> The content of the equil.mdp and topol.top files is shown as follow.
>
> ​​
> equil.mdp:
>
> title       = NVT Equilibration for Na in 56 water
> ; Run parameters
> integrator  = md        ; leap-frog integrator
> nsteps      = 250000000     ; 2 * 50,000 = 100 ps ;(100 ns)
> dt          = 0.002     ; 2 fs
> ​​
> define          = -DPOSRES
> ; Output control
> nstxout     = 500       ; save coordinates every 0.1 ps
> nstvout     = 500       ; save velocities every 1 ps
> nstenergy   = 500       ; save energies every 1 ps
> nstlog      = 500       ; update log file every 1 ps
> nstxtcout                = 500
> xtc-precision            = 1000
> ; Bond parameters
> continuation    = yes           ; Restarting after NVT
> constraint_algorithm = ; holonomic constraints
> constraints =h-angles   ;water both bond and angle constrained
> ;constraints    = hbonds            ; all bonds (even heavy atom-H bonds)
> constrained
> lincs_iter  = 1                 ; accuracy of LINCS
> lincs_order = 4                 ; also related to accuracy
> ; Neighborsearching
> ns_type     = grid      ; search neighboring grid cels
> nstlist     = 1         ; 10 fs
> rlist       = 1.0       ; short-range neighborlist cutoff (in nm)
> ; vdw
> vdw-type                 = Cut-off
> rvdw        = 3.0       ; short-range van der Waals cutoff (in nm)
> ; Electrostatics
> coulombtype = cut-off ;Reaction-Field ;Generalized-Reaction-Field ;cut-off
> rcoulomb    = 3.0       ; short-range electrostatic cutoff (in nm)
> epsilon_rf      = 2.0
> epsilon_r       = 0.5
> cutoff-scheme       = group
> ;pme_order  = 4         ; cubic interpolation
> ;fourierspacing = 0.16      ; grid spacing for FFT
> ; Temperature coupling is on
> ; annealing = single
> ; annealing_time = 0     400    1600    2400  4000  5600   7200  10400
> 13600  16800 ; 230  260  270  300  330   360
> ; annealing_temp  =0     20      40     60    80    100    120   140
>  160   180 ;  200   220   240   260   280   300
> ; annealing_npoints = 10
> ; Temperature coupling is on
> tcoupl      = nose-hoover ;v-rescale; berendsen; nose-hoover    ; More
> accurate thermostat
> tc-grps     = system ; three coupling groups - more accurate
> tau_t       = 0.1    ; time constant, in ps
> ref_t       = 547     ; reference temperature, one for each group, in K
> ; Pressure coupling is on
> pcoupl      = no        ; Pressure coupling on in NPT
> pcoupltype  = semiisotropic         ; uniform scaling of x-y box vectors,
> independent z
> tau_p       = 5.0                   ; time constant, in ps
> ref_p       = 1.0   1.0             ; reference pressure, x-y, z (in bar)
> compressibility = 4.5e-5    4.5e-5  ; isothermal compressibility, bar^-1
> ; Periodic boundary conditions
> pbc         = no        ; 3-D PBC
> ; Velocity generation
> gen_vel     = no    ; assign velocities from Maxwell distribution
> gen_temp    = 133       ; temperature for Maxwell distribution
> gen_seed    = -1        ; generate a random seed
> ; COM motion removal
> ; These options remove motion of the protein/bilayer relative to the
> solvent/ions
> nstcomm         = 1
> comm-mode       = Linear ;ANGULAR;Linear
> comm-grps       =
>
>
>
>
>
>
>
> ################################################################
> topol.top :
>
> ;
> ;   File '250_NA_cen.top' was generated
> ;   By user: onbekend (0)
> ;   On host: onbekend
> ;   At date: Tue Aug 12 09:46:55 2014
> ;
> ;   This is a standalone topology file
> ;
> ;   It was generated using program:
> ;   pdb2gmx_d - VERSION 4.5.5
> ;
> ;   Command line was:
> ;   pdb2gmx_d -f water250_NA_cen.pdb -p 250_NA_cen.top -o 250_NA_cen.gro
> ;
> ;   Force field was read from the standard Gromacs share directory.
> ;
>
> ; Include forcefield parameters
> #include "amber94.ff/forcefield.itp"
>
> [ moleculetype ]
> ; Name            nrexcl
> Ion                 3
>
> [ atoms ]
> ;   nr       type  resnr residue  atom   cgnr     charge       mass
> typeB    chargeB      massB
> ; residue   1 NA  rtp NA   q +1.0
>      1         Na      1     NA     NA      1          1      22.99   ;
> qtot 1
>
> ; Include Position restraint file
> #ifdef POSRES
> ;#include "posre.itp"
> ;[ position_restraints ]
> ;  i funct       g         r(nm)       k
> ;   1    2        1          3.0       30.0
> #endif
>
>
> ; Include water topology
> #include "amber94.ff/spce.itp"
>
> #ifdef POSRES_WATER
> ; Position restraint for each water oxygen
> ;[ position_restraints ]
> ;  i funct       fcx        fcy        fcz
> ;  1    1       1000       1000       1000
> [ position_restraints ]
> ;  i funct       g         r(nm)       k
>    1    2        1          3.0       30.0
> ;   2    2        1          3.0       30.0
> ;   3    2        1          3.0       30.0
> #endif
>
> ; Include topology for ions
> #include "amber94.ff/ions.itp"
>
> [ system ]
> ; Name
> Protein
>
> [ molecules ]
> ; Compound        #mols
> Ion                 1
> SOL               250
>
>
>
>
>
> liuyong_1007 at dicp.ac.cn
>
> From: liuyong_1007 at dicp.ac.cn
> Date: 2014-11-16 10:32
> To: Discussion list for GROMACS users
> Subject: Re: [gmx-users] Flat-bottomed position restraint set
> Dear Erik,
>
> Thanks for your reply! I use a larger cut-off of the non-bonded
> interaction. But water molecules sitll escape to the vacuum after a 220 ns
> MD run.
> The restraints are applied to the O atoms.  Could you please help me to
> check whether the sets of the equil.mdp and topol.top files are right or
> wrong?
> The files are shown in the attachment. Thank you very much !
>
> Best regards,
> Yong Liu
>
> liuyong_1007 at dicp.ac.cn
>
> From: Erik Marklund
> Date: 2014-11-14 19:24
> To: <gmx-users at gromacs.org>
> Subject: Re: [gmx-users] Flat-bottomed position restraint set
> Dear Liuyoung,
>
> Remember that the electrostatic screening is much weaker in gas phase
> systems such as yours. I would use a larger cut-off for the non-bonded
> interactions to effectively have all vs all interactions. Your system is
> fairly small so you won't be simulating for that long anyway. Note that
> gromacs can be told to keep the neighbour list from the first frame, which
> would speed up things a bit.
>
> Kind regards,
> Erik
>
>
> Erik Marklund, PhD
> Postdoctoral Research Fellow, Fulford JRF
>
> Department of Chemistry
> Physical & Theoretical Chemistry Laboratory
> University of Oxford
> South Parks Road
> Oxford
> OX1 3QZ
>
> On 14 Nov 2014, at 07:15, liuyong_1007 at dicp.ac.cn<mailto:
> liuyong_1007 at dicp.ac.cn> wrote:
>
> Hi Justin!
> I try to apply restraint just to the O atom. There is another question
> about  the set of the coulomb and the vdw interaction for the ion-water
> cluster Na(H2O)_250. My set is shown as follow:
>
> ; vdw
> vdw-type      = Cut-off
> rvdw            = 1.0           ; short-range van der Waals cutoff (in nm)
> ; Electrostatics
> coulombtype     = cut-off ;Reaction-Field ;Generalized-Reaction-Field
> ;cut-off ^M
> rcoulomb        = 1.2           ; short-range electrostatic cutoff (in
> nm)^M
> epsilon_rf      = 2.0
> epsilon_r       = 0.5
> cutoff-scheme    = group
>
> The radius of the cluster is about 1.12 nm. Is this set appropriate ?
>
> Best regards,
> Yong Liu
>
>
>
> liuyong_1007 at dicp.ac.cn<mailto:liuyong_1007 at dicp.ac.cn>
>
> From: Justin Lemkul
> Date: 2014-11-14 10:40
> To: gmx-users
> Subject: Re: [gmx-users] Flat-bottomed position restraint set
>
>
> On 11/13/14 7:15 PM, liuyong_1007 at dicp.ac.cn wrote:
> Hi Justin!
> I use the geometric center of the sphere as the reference coordinates.
> However, there are still water molecules escaping to the vacuum. Is there
> other way to aovid this ?
>
>
> This shouldn't happen if things are set up right.  Try applying the
> restraint
> just to the O atom, not all 3 atoms of the water individually.  We do this
> routinely and it works quite well.
>
> -Justin
>
> Best regards,
> Yong Liu
>
>
>
> liuyong_1007 at dicp.ac.cn
>
> From: Justin Lemkul
> Date: 2014-11-12 21:54
> To: gmx-users
> Subject: Re: [gmx-users] Flat-bottomed position restraint set
>
>
> On 11/12/14 8:08 AM, liuyong_1007 at dicp.ac.cn wrote:
> Dear Gromacs users!
>
> I use the flat-bottomed position restraints to avoid the molecules to
> escape from the clusters to the vacuum. The parameters of the restraints on
> water molecules are shown as follow:
> #ifdef POSRES_WATER
> ; Position restraint for each water oxygen
> [ position_restraints ]
> ;  i   funct       g         r(nm)       k
>     1    2        1          3         30.0
>     2    2        1          3         30.0
>     3    2        1          3         30.0
> #endif
> However, water molecules still escape away from the cluster to the vacuum
> after about 100 ns MD run at 550K. What should be done with the restraints
> to avoid the water molecules to escape to the vacuum?
>
>
> Are you using suitable reference coordinates passed to grompp -r?  If
> you're
> not, then the reference coordinates are whatever is in grompp -c, which
> means
> your water molecules can deviate up to 3 nm from that position without
> penalty.
>   For a sphere, the reference coordinates for the water should be the
> geometric
> center of the sphere.
>
> -Justin
>
>
> --
> ==================================================
>
> Justin A. Lemkul, Ph.D.
> Ruth L. Kirschstein NRSA Postdoctoral Fellow
>
> Department of Pharmaceutical Sciences
> School of Pharmacy
> Health Sciences Facility II, Room 629
> University of Maryland, Baltimore
> 20 Penn St.
> Baltimore, MD 21201
>
> jalemkul at outerbanks.umaryland.edu | (410) 706-7441
> http://mackerell.umaryland.edu/~jalemkul
>


More information about the gromacs.org_gmx-users mailing list