[gmx-users] unable to equilibrate protein in membrane with NPT

[Justin A. Lemkul]

Peter C. Lai wrote:
> On 2011-04-08 08:02:14AM -0500, Justin A. Lemkul wrote:
>>
>> Peter C. Lai wrote:
>>> Hello again
>>>
>>> In my protein-membrane-water-ion system (inserted via g_membed) I have run 
>>> a 1ns NVT equilibration with the protein restrained and now I am trying to 
>>> equilibrate with NPT and LINCS/mdrun is crashing after about 5000 iterations.
>>>
>>> Here is my NPT mdp file:
>>>
>>> define                  =-DPOSRES
>>> integrator               = md
>>> ; Start time and timestep in ps
>>> tinit                    = 0
>>> dt                       = 0.002
>>> nsteps                   = 500000
>>> init_step                = 0
>>> comm-mode                = Linear
>>> nstcomm                  = 1
>>> comm-grps                = Protein_POPC SOL_CL
>>> nstxout         = 100           ; save coordinates every 0.2 ps
>>> nstvout         = 100           ; save velocities every 0.2 ps
>>> nstenergy       = 100           ; save energies every 0.2 ps
>>> nstlog          = 100           ; update log file every 0.2 ps
>>> continuation    = yes
>>> constraint_algorithm = lincs    ; holonomic constraints
>>> constraints     = all-bonds     ; all bonds (even heavy atom-H bonds)
>>> lincs_iter      = 1             ; accuracy of LINCS
>>> lincs_order     = 4             ; also related to accuracy
>>> ns_type         = grid          ; search neighboring grid cells
>>> nstlist         = 5             ; 10 fs
>>> rlist           = 1.2           ; short-range neighborlist cutoff (in nm)
>>> rlistlong       = 1.4
>>> rcoulomb        = 1.2           ; short-range electrostatic cutoff (in nm)
>>> rvdw            = 1.2           ; short-range van der Waals cutoff (in nm)
>>> vdwtype         = switch
>>> rvdw_switch     = 0.8
>>> coulombtype     = PME
>>> pme_order       = 4             ; cubic interpolation
>>> fourierspacing  = 0.16          ; grid spacing for FFT
>>> tcoupl          = Nose-Hoover
>>> tc-grps         = Protein POPC SOL_CL
>>> tau-t                    = 0.5 0.5 0.5
>>> ref-t                    = 300 300 300
>>> gen-vel                  = no
>>> pcoupl          = Parrinello-Rahman         ; Pressure coupling on in NPT
>>> pcoupltype      = semiisotropic             ;
>>> tau_p           = 5.0                           ; time constant, in ps
>>> ref_p           = 1.01325 1.01325
>>> compressibility = 4.5e-5        4.5e-5
>>>
>>> Any suggestions?
>> Try using the Berendsen barostat.  P-R allows for wider oscillations that can 
>> lead to instability in incompletely equilibrated systems.
>>
>> -Justin
>>
> 
> Switching to V-rescale and tau_t to 0.1 0.1 0.1 also crashes.

My recommendation was for the barostat, not the thermostat.

> LINCS is complaining about C-H interactions in both the mobile lipid chains
> as well as atoms restrained by posre.itp
> 
> Maybe I will run without LINCS all-bond?

No, that will make your life worse.  Just because the constraints fail doesn't 
mean that you need to disable LINCS.

http://www.gromacs.org/Documentation/Errors#LINCS.2fSETTLE.2fSHAKE_warnings
http://www.gromacs.org/Documentation/Terminology/Blowing_Up

-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

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