[gmx-users] Gromacs 2018 beta: g_membed broken?

[Daniel Bauer]

Hello,

Im using g_membed to insert potassium channels into POPC membranes. My
"standard protocol" which worked very well for me on all proteins I used
so far seems to not work with gromacs 2018 (using the exact same input
files).


This is what I am doing:

1) orienting the protein onto the membrane (lambada)

2) Merging the topology files for the membrane and protein

3) grompp: gmx grompp -f g_membed.mdp -n -c input.gro -r input.gro -o
g_membed.tpr --maxwarn 2

4) gmx mdrun -deffnm g_membed -membed embed.dat -mn index.ndx -mp
embedded.top -c embedded.top # selecting protein and membrane group


First thing I noticed is that the number of warnings i have to ignore
during grompp increased by one because of the new PME warning for
charged systems. My protein has a net charge, so I have to choose if i
insert counter ions before or after protein insertion (and I chose to do
it after insertion).

So the warnings I get are:

WARNING 1: You are using Ewald electrostatics in a system with net
charge  // I think I can ignore this since I neutralize right after
insertion.

WARNING 2: Can not exclude the lattice Coulomb energy between energy
groups // g_membed does not support verlet


Finally when I run the insertion simulation on the 2018 beta I get the
following error:

step 0: One or more water molecules can not be settled.
Check for bad contacts and/or reduce the timestep if appropriate.

Back Off! I just backed up step0b.pdb to ./#step0b.pdb.1#

Back Off! I just backed up step0c.pdb to ./#step0c.pdb.1#
Wrote pdb files with previous and current coordinates
step 0Warning: Only triclinic boxes with the first vector parallel to
the x-axis and the second vector in the xy-plane are supported.
         Box (3x3):
            Box[    0]={        -nan,         -nan,         -nan}
            Box[    1]={        -nan,         -nan,         -nan}
            Box[    2]={        -nan,         -nan,         -nan}
         Can not fix pbc.

This seems to happen independent of my system (I tried several of my
membrane patches and proteins that worked previously). All of them
insert into the membrane with 2016 successfully but fail on the 2018
version. So far I tried to increase the number of steps, reduced the
timestep and played arround with xyinit and zinit as well as setting
-DFLEXIBLE in the mdp, but no success.

Any idea whats wrong? Or should I move this issue directly to redmine
since its related to a version update?


This is my embed.dat file:

nxy        = 1000
xyinit        = 0.1
xyend        = 1.0
rad        = 0.22
ndiff        = 0
asymmetry    = no
pieces        = 1
maxwarn        = 1


And my mdp:

; Run parameters
integrator    = md
nsteps        = 1000       
dt            = 0.002      

; I have posres for the protein enabled, although this is irrelevant
because the protein group coords are frozen. Same issue if i disable
this and remove posres!
refcoord_scaling = all

;define = -DFLEXIBLE ; does not help

; membed options
energygrps    =    insertion_group
freezegrps    =    insertion_group
freezedim    =    Y Y Y
energygrp_excl    =    insertion_group insertion_group

; OUTPUT CONTROL OPTIONS
nstxout                  = 0
nstvout                  = 0
nstfout                  = 0
nstlog                   = 1000
nstenergy                = 100  
nstxtcout                = 1000 

; Bond parameters
continuation    = no          
constraint_algorithm = lincs   
lincs_iter    = 2                
lincs_order    = 8

; CHARMM36 params
constraints = h-bonds
cutoff-scheme = Group        ; Verlet is not supported by g_membed
vdwtype = Cut-off
vdw-modifier = force-switch
rlist = 1.2
rvdw = 1.2
rvdw-switch = 0.8
coulombtype = PME
rcoulomb = 1.2
DispCorr = no

; Temperature coupling
tcoupl        = V-rescale
tc-grps        = Other Water_and_ions Protein
tau_t        = 0.1    0.1 0.1
ref-t =  298.0  298.0  298.0


; Pressure coupling
pcoupl                  = Berendsen
pcoupltype              = semiisotropic
tau_p                   = 5.0
compressibility         = 4.5e-5  4.5e-5
ref_p                   = 1.0     1.0


; Periodic boundary conditions
pbc            = xyz


; Velocity generation
gen_vel        = yes 
gen_temp    = 298.0   
gen_seed    = -1  

; COM motion removal
; These options remove motion of the protein/bilayer relative to the
solvent/ions
comm-mode    = Linear
comm-grps    = Other Water_and_ions Protein


Best regards,

Daniel


-- 
Daniel Bauer, M.Sc.

TU Darmstadt
Computational Biology & Simulation
Schnittspahnstr. 2
64287 Darmstadt
bauer at cbs.tu-darmstadt.de

Don't trust atoms, they make up everything.