[gmx-users] Domain decomposition error tied to free energy perturbation

[Ryan Muraglia]

Hello,

I have been attempting to carry out some free energy calculations, but to
verify the sanity of my parameters, I decided to test them on a structure I
knew to be stable -- the lysozyme from Lemkul's lysozyme in water tutorial.

I chose the L75A mutation because it is out on the surface to minimize the
"difficulty of the transformation."
Using my regular mdp file (even with my mutatation topology generated with
the pmx package), my minimization runs to completion with no errors.

Once I introduce the following lines to my mdp file:

"
; Free energy calculations
free_energy = yes
delta_lambda = 0 ; no Jarzynski non-eq
calc_lambda_neighbors = 1 ; only calculate energy to immediate neighbors
(suitable for BAR, but MBAR needs all)
sc-alpha                 = 0.5
sc-coul                  = no
sc-power                 = 1.0
sc-sigma                 = 0.3
couple-moltype           = Protein_chain_A  ; name of moleculetype to
decouple
couple-lambda0           = vdw-q      ; all interactions
couple-lambda1           = vdw     ; remove electrostatics, only vdW
couple-intramol          = no
nstdhdl                  = 100

; lambda vectors ; decharging only.
; init_lambda_state   0   1   2   3   4   5   6   7   8   9   10
init_lambda_state = 00
coul_lambdas =        0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
vdw_lambdas =         0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
bonded_lambdas =      0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ; match
vdw
mass_lambdas =        0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ; match
vdw
temperature_lambdas = 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ; not
doing simulated tempering
"

I notice two things:
1) Running grompp to generate the tpr file takes much longer
2) The minimization fails to run due the following error related to domain
decomposition:

"
Fatal error:
There is no domain decomposition for 4 ranks that is compatible with the
given box and a minimum cell size of 5.51109 nm
Change the number of ranks or mdrun option -rdd
Look in the log file for details on the domain decomposition
"

I noted that it lists a two-body bonded interaction with a strangely large
distance:

"
Initial maximum inter charge-group distances:
    two-body bonded interactions: 5.010 nm, LJC Pairs NB, atoms 1074 1937
  multi-body bonded interactions: 0.443 nm, Proper Dih., atoms 1156 1405
Minimum cell size due to bonded interactions: 5.511 nm
"

Atom 1074 corresponds to a hydrogen off the beta-carbon of proline 70, and
atom 1937 refers to a hydrogen on arginine 128. Neither residue is part of
the protein that is being mutated, and they certainly should not be bonded.
The [bonds] directive in the topology confirms that there should be no
interaction between these atoms.
To force the run to begin to get more information on the nature of the
error, I gave mdrun the -nt 1 option, and got the following warning at the
beginning of the minimization (which goes on to end prematurely prior to
reaching the desired Fmax):

"
WARNING: Listed nonbonded interaction between particles 1 and 195
at distance 2.271 which is larger than the table limit 2.200 nm.
"

I'm at a loss in terms of understanding why the addition of my FEP
parameters is causing this error, and appears to be causing the grompp
parser to decide that there is a bond where there shouldn't be, forcing the
minimimum box size to exceed what makes sense for domain decomposition.

Additional information that may be relevant: I am using the amber99sb
forcefield with explicit tip3p waters. I am attempting steepest descent
minimization. rcoulomb and rvdw are both set to 1.2.

Any advice would be greatly appreciated. Thank you!


-- 
Ryan Muraglia
rmuraglia at gmail.com