[gmx-users] segfault with an otherwise stable system when I turn on FEP (complete decoupling)

[Christopher Neale]

Dear Users:

I have a system with water and a drug (54 total atoms; 27 heavy atoms). The system is stable when I simulate it for 1 ns. However, Once I add the following options to the .mdp file, the run dies after a few ps with a segfault. 

free-energy = yes
init-lambda = 1
couple-lambda0 = vdw-q
couple-lambda1 = none
couple-intramol = no
couple-moltype = drug

I do not get any step files or any lincs warnings. If I look at the .xtc and .edr files, there is no indication of something blowing up before the segfault. I have also verified that the drug runs without any problem in vacuum. I get the same behaviour if I remove constraints and use a timestep of 0.5 fs. The segfault is reproducible with v4.6.1 and v4.6.3. I am using the charmm FF, but I converted all UB angles in my drug to type-1 angles and still got the segfault. I also get the segfault with particle decomposition and/or while running a single thread. I am currently using the SD integrator, but I get the same segfault with md and md-vv. Couple-intramol=yes doesn't resolve it, neither does using separate T-coupling groups for the water and drug. Neither does turning off pressure coupling.

Here is the .mdp file that works fine, but gives me a segfault when I add the free energy stuff (above):

constraints = all-bonds
lincs-iter =  1
lincs-order =  6
constraint_algorithm =  lincs
integrator = sd
dt = 0.002
tinit = 0
nsteps = 100000
nstcomm = 1
nstxout = 0
nstvout = 0
nstfout = 0
nstxtcout = 500
nstenergy = 500
nstlist = 10
nstlog=0 ; reduce log file size
ns_type = grid
vdwtype = cut-off
rlist = 0.8
rvdw = 0.8
rcoulomb = 0.8
coulombtype = cut-off
tc_grps             =  System
tau_t               =  1.0
ld_seed             =  -1
ref_t = 310
gen_temp = 310
gen_vel = yes
unconstrained_start = no
gen_seed = -1
Pcoupl = berendsen
pcoupltype = isotropic
tau_p = 4 
compressibility = 4.5e-5 
ref_p = 1.0 

I do realize that some of these settings are not ideal for a production run. I started with the real Charmm cutoffs + PME, etc, (which also gives the segfault) but this is what I am using right now for quick testing.

The only thing keeping me from filing a redmine issue is that if I remove my drug and do the FEP on one of the water molecules (using the FEP code listed above), I have no segfault. Therefore it is clearly related to the drug, whose parameters I built so I may have caused the problem somehow. Nevertheless, the drug runs fine in water and in vacuum without the FEP code, so I can't imagine what could be causing this segfault (also, the fact that it's a segfault means that I don;t get any useful info from mdrun as to what might be going wrong).

Thank you,
Chris.