[gmx-users] particle decomp crashes while domain decomp does not

[Mark Abraham]

On 22/09/2011 8:45 AM, Florian Altvater wrote:
> Hi,
> I don't know if my question somehow got lost, or if I forgot to provide
> some necessary information. Anyway it would be of great help to me, if
> someone could tell me her/his opinion on the problem. I need to use
> particle decomposition as I want to introduce distance restraints from
> 2-6nm, which domain decomposition can't handle apparently.
>
> Thanks again for your help
> Florian
>
> Am 19.09.2011 14:30, schrieb Florian Altvater:
>> Hi there,
>> I haven't found a post describing this problem, so here I am.
>>
>> I'm running GROMACS 4.3.5 with the OPLS-AA force-field and the TIP3P
>> water model on a 170kDa protein dimer with sodium ions to neutralize the
>> charges. I added two new residues to the force-field (which seem to
>> behave normally). Now the problem:
>>
>> When I run with domain decomposition everything works fine during 9ns.
>> The structure seems to equilibrate more or less and the overall domain
>> conformation doesn't seem to alter much.
>>
>> When I switch on particle decomposition, however, the simulation crashes
>> after 2-3ns due to a lincs warning. It seems to happen always at the
>> same residue, Arg318 which gets pulled/pushed out of the helix it is in.
>> The closest modified residue is 326 which doesn't seem to be disturb the
>> system too much.
>>
>> Anyway it seems strange to me that the method of decomposition seems to
>> change the outcome of the simulation a lot (and does not only reflect in
>> the computing time). Now I don't understand too much of the algorithms
>> used. And as I couldn't find anything bout adjusting parameters if using
>> the one or the other I'm a bit in the dark about what to do now.
>>
>> Did you have similar experiences? How do you suggest that I tackle the
>> problem? Do you need more information? In case it is of any relevance
>> here the parameters of the mdp file:
>>
>>      dt                       = 0.002
>>      nsteps                   = 1500000
>>      comm-grps                = system
>>      nstxout                  = 20000
>>      nstvout                  = 20000
>>      nstfout                  = 20000
>>      nstlog                   = 20000
>>      nstcalcenergy            = 1
>>      nstxtcout                = 500
>>      xtc-grps                 = protein
>>      energygrps               = protein
>>      ns_type                  = grid
>>      rlist                    = 1.0
>>      coulombtype              = PME
>>      vdw-type                  = Cut-off
>>      rcoulomb                 = 1.0
>>      rvdw                     = 1.0
>>      DispCorr                 = EnerPres
>>      tcoupl                   = V-rescale
>>      tc-grps                  = Protein Non-Protein
>>      tau_t                    = 0.1    0.1
>>      ref_t                    = 325     325
>>      pcoupl                   = Parrinello-Rahman
>>      pcoupltype               = isotropic
>>      tau_p                    = 2.0
>>      compressibility          = 4.5e-5
>>      ref_p                    = 1.0
>>      andersen_seed            = -1
>>      gen_vel                  = no
>>      gen-temp                 = 325
>>      gen-seed                 = -1
>>      constraints              = all-bonds
>>      constraint_algorithm     = lincs
>>      continuation             = yes
>>      lincs_order              = 4
>>      lincs_iter               = 1
>>      disre                    = simple
>>      disre_weighting          = conservative
>>      disre_fc                 = 800
>>      disre_tau                = 0

The different decomposition algorithm can be enough for the chaotic 
simulation to diverge enough that a marginally stable situation can 
either explode or move to a stable point. However Parrinello-Rahman is 
also unsuitable for equilibration, and that could perhaps amplify the 
former effect. Equilibrate with Berendsen, and if you still encounter PD 
vs DD issues, please try again with the latest 4.5.5 release.

Mark