[gmx-users] Constraints not working in pull code (sometimes, sometimes not)

[Justin A. Lemkul]

Krapnik wrote:
> Dear all, 
> 
> we are trying to simulate molecule in specified depth in octanol to 
> calculate logP. 
> for this purpose we have used box with slab of octanol, molecule in 
> specified distance in z-axis and  this mdp: 
> 
> integrator      = md
> dt              = 0.002
> tinit           = 0
> nsteps          = 5125000       ; 10.250 ns
> nstcomm         = 1
> comm_mode       = Linear
> ; Output parameters
> nstxout         = 0     ; every 1 ns
> nstvout         = 0
> nstfout         = 0
> nstxtcout       = 500           ; every 1 ps
> nstenergy       = 500
> ; Bond parameters
> constraints             = all-bonds
> ; Single-range cutoff scheme
> nstlist         = 5
> ns_type         = grid
> rlist           = 1.4
> rcoulomb        = 1.4
> rvdw            = 1.4
> ; PME electrostatics parameters
> coulombtype     = PME
> fourierspacing  = 0.12
> fourier_nx      = 0
> fourier_ny      = 0
> fourier_nz      = 0
> pme_order       = 4
> ewald_rtol      = 1e-5
> optimize_fft    = yes
> ; Berendsen temperature coupling is on in two groups
> Tcoupl          = V-rescale
> tc_grps         = System
> tau_t           = 0.1
> ref_t           = 310
> ; Pressure coupling is on
> Pcoupl          = berendsen
> pcoupltype      = anisotropic   ; anisotropic pressure coupling
> tau_p           = 10
> compressibility = 0 0 4.5e-5 0 0 0  ; since octanol tries to separate 
> from water and shortening xy plane
> ref_p           = 1.0 1.0 1.0 0 0 0 
> ; Pull code
> pull            = constraint    ;
> pull_geometry   = distance      ; Pull along the vector connecting the 
> two groups. Components can be selected with pull_dim.
> pull_dim        = N N Y         ; put potential only to axis z
> pull_start      = yes           ; define initial COM distance > 0
> pull_ngroups    = 1             ; one group to pull
> pull_group0     = OCT           ; reference group (can be empty to set 
> it to 0,0,0)
> pull_group1     = DRG
> pull_rate1      = 0             ; 0.01 nm per ps = 10 nm per ns = 10 m per s
> pull_k1         = 2000          ; kJ mol^-1 nm^-2
> 
> however when I list pull-x.xvg I obtain:
> @ s0 legend "0 Z"
> @ s1 legend "1 dZ"
> 0.0000  6.06871 -2.81855
> 0.0200  6.07793 -2.81861
> 0.0400  6.08787 -2.81856
> 0.0600  6.08462 -2.81855
> 0.0800  6.10028 -2.82007
> 
> The situation is even more strange when almost the same mdp with DOPC 
> membrane works fine as can be seen in pull-x.xvg:
> 0.0000  3.72735 -1.60929
> 0.0200  3.72727 -1.60929
> 0.0400  3.7272  -1.60929
> 0.0600  3.72716 -1.60929
> 0.0800  3.72715 -1.60929
> 0.1000  3.72717 -1.60929
> 
> Differences in the working case in mdp are only this:
> compressibility = 4.5e-5 4.5e-5 4.5e-5 0 0 0                             
>                                           
> ref_p           = 1.0 1.0 1.0 0 0 0            
> pull_group0     = DOPC
> 
> Strangerly enough - the same happens in GMX4.0.7 and in GMX4.5.3
> 
> I wonder where the difference is, as both logs states that 
> Will apply constraint COM pulling in geometry 'distance'
> between a reference group and 1 group
> Pull group 0:  6912 atoms, mass 100624.859
> Pull group 1:    15 atoms, mass   194.194
> 
> Will apply constraint COM pulling in geometry 'distance'
> between a reference group and 1 group
> Pull group 0:  9570 atoms, mass 124631.453
> Pull group 1:    15 atoms, mass   146.146
> 
> Where should I look and repair? 
> 

I don't see any real problem here.  The DOPC membrane probably holds its shape 
better than the octanol slab and thus the reference position (which determines 
the constraint) is less mobile.  As such, the constraint is more stable in DOPC 
than in octanol.  The slight bump in the octanol system amounts to a change of 
0.16%, which I'd say is quite small.  You're also looking at the first few 
frames of the simulation, which may amount to equilibration, but you haven't 
mentioned if you've done anything prior to this run.

As long as there is no systematic change in position, I'd say there's no problem.

-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

========================================