[gmx-users] Can't do parallell free energy calculation with couple_intramol=no -- workaround?
Michael Daily
mdaily.work at gmail.com
Thu Mar 10 21:04:25 CET 2016
Hi Szilard,
Thanks for the advice. I will try out v5.1 accordingly. Re attachments,
I've now copied my mdp and md output below:
### mdp ###
; Run control
integrator = sd ; Langevin dynamics
tinit = 0
dt = 0.002
nsteps = 1000000 ; 2 ns
nstcomm = 100
define=-DPOSRES
; Output control
nstxout = 0
nstvout = 0
nstfout = 0
nstlog = 5000
nstenergy = 5000
nstxtcout = 5000
energygrps = Protein SOL
;energygrp-excl=Protein Protein
; Neighborsearching and short-range nonbonded interactions
cutoff-scheme = verlet
;cutoff-scheme=group
nstlist = 20
ns_type = grid
pbc = xyz
rlist = 1.2
; Electrostatics
coulombtype = PME
rcoulomb = 1.2
; van der Waals
vdwtype = cutoff
vdw-modifier = potential-switch
rvdw-switch = 1.0
rvdw = 1.2
; Apply long range dispersion corrections for Energy and Pressure
DispCorr = EnerPres
; Spacing for the PME/PPPM FFT grid
fourierspacing = 0.12
; EWALD/PME/PPPM parameters
pme_order = 6
ewald_rtol = 1e-06
epsilon_surface = 0
; Temperature coupling
; tcoupl is implicitly handled by the sd integrator
tc_grps = system
tau_t = 1.0
ref_t = 300
; Pressure coupling is on for NPT
Pcoupl = Parrinello-Rahman
Pcoupltype = isotropic
refcoord-scaling = com
tau_p = 1.0
compressibility = 4.5e-05
ref_p = 1.0
; Free energy control stuff
free_energy = yes
init_lambda_state = my_init_lambda
delta_lambda = 0
calc_lambda_neighbors = 1 ; only immediate neighboring windows
; Vectors of lambda specified here
; Each combination is an index that is retrieved from init_lambda_state for
each simulation
coul_lambdas = 0.0000 0.0625 0.1250 0.1875 0.2500 0.3125 0.3750 0.4375
0.5000 0.5625 0.6250 0.6875 0.7500 0.8125 0.8750 0.9375 1.0000
; Options for the decoupling
sc-alpha = 0.5
sc-coul = no ; linear interpolation of Coulomb (none
in this case)
sc-power = 1.0
sc-r-power = 6
sc-sigma = 0.3
couple-moltype = Protein ; name of moleculetype to (de)couple
couple-lambda0 = vdw ; only van der Waals interactions
couple-lambda1 = vdw-q ; everything
couple-intramol = no
;couple-intramol=yes
nstdhdl=5000
separate-dhdl-file=no ;more efficient
; Do not generate velocities
gen_vel = no
; options for bonds
constraints = all-bonds ; we only have C-H bonds here
; Type of constraint algorithm
constraint-algorithm = lincs
; Constrain the starting configuration
; since we are continuing from NPT
continuation = yes
; Highest order in the expansion of the constraint coupling matrix
;lincs-order = 12
table-extension=1.5
### output (part with the error) ###
Initializing Domain Decomposition on 16 ranks
Dynamic load balancing: auto
Will sort the charge groups at every domain (re)decomposition
Initial maximum inter charge-group distances:
two-body bonded interactions: 2.343 nm, LJC Pairs NB, atoms 4 45
multi-body bonded interactions: 0.383 nm, Proper Dih., atoms 1 9
Minimum cell size due to bonded interactions: 2.577 nm
Maximum distance for 5 constraints, at 120 deg. angles, all-trans: 0.728 nm
Estimated maximum distance required for P-LINCS: 0.728 nm
Using 0 separate PME ranks, as there are too few total
ranks for efficient splitting
Scaling the initial minimum size with 1/0.8 (option -dds) = 1.25
Optimizing the DD grid for 16 cells with a minimum initial size of 3.222 nm
The maximum allowed number of cells is: X 0 Y 0 Z 1
-------------------------------------------------------
Program mdrun, VERSION 5.0.7
Source code file:
/home1/01777/mdaily/software/gromacs-5.0.7/src/gromacs/mdlib/domdec.c,
line: 6902
Fatal error:
There is no domain decomposition for 16 ranks that is compatible with the
given box and a minimum cell size of 3.22187 nm
Change the number of ranks or mdrun option -rdd or -dds
Look in the log file for details on the domain decomposition
For more information and tips for troubleshooting, please check the GROMACS
website at http://www.gromacs.org/Documentation/Errors
-------------------------------------------------------
Thanks,
====================================
Michael D. Daily
Research Scientist, Biochemistry and Molecular Biology
University of Texas Medical Branch at Galveston
509-713-0726
On Thu, Mar 10, 2016 at 1:55 PM, Szilárd Páll <pall.szilard at gmail.com>
wrote:
> On Thu, Mar 10, 2016 at 8:14 PM, Michael Daily <mdaily.work at gmail.com>
> wrote:
>
> > Hi,
> >
> > I am trying to do a solvation free energy calculation on a small peptide
> in
> > gromacs 5.0. The peptide is in an extended conformation about 2.5 nm
> long,
> > and I use position restraints on the heavy atoms to keep it in this
> > conformation. Due to the nature of the calculation I set
> > "couple_intramol=no" in the mdp file (attached for lambda=0).
> >
> > However, due to the way couple_intramol is implemented in gromacs
> > (long-range exclusions), the longest bonded interaction is as long as the
> > molecule, leading to a minimum cell size of over 3 nm, which makes domain
> > decomposition all but impossible (see the attached log file). Is there a
> > way to limit long-range exclusions and still get reasonable charging free
> > energies? I have a few ideas but would like some opinions on whether they
> > are physically sound:
> >
> > 1) change "couple_intramol" to yes and generate intramolecular exclusions
> > manually with some cutoff. Since rlist/rcoulomb/rvdw are 1.2 nm for this
> > simulation, I might use a cutoff of 1.5 to provide a buffer.
> > 2) Use "energygrp-excl = Protein Protein" to turn off nonbonded
> > interactions within the peptide. The only problem with this is it forces
> > the use of group instead of Verlet cutoff scheme which seems to slow
> > performance.
> > 3) Revert to Gromacs 4.7 and use particle decomposition
> >
>
> I guess you mean v4.6.7. As far as I know there is the decomposition should
> not interfere with the physics. However, you could instead consider using
> 5.1 with multi-threading. Unless something is unsupported or you happen to
> stumble into code not OpenMP-parallelized, you should be able to run on up
> to 24-32 cores without domain-decomposition (also note that 5.1 support
> running separate PME without DD).
>
>
> > The log file from my (attempted) simulation with these parameters is also
> > attached.
> >
>
> The user's list does not accept attachments, so you'll have to share links
> to files uploaded somewhere else.
>
> --
> Szilárd
>
>
>
> >
> > Thanks,
> >
> > ====================================
> > Michael D. Daily
> > Research Scientist, Biochemistry and Molecular Biology
> > University of Texas Medical Branch at Galveston
> > 509-713-0726
> >
> > --
> > Gromacs Users mailing list
> >
> > * Please search the archive at
> > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> > posting!
> >
> > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> >
> > * For (un)subscribe requests visit
> > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> > send a mail to gmx-users-request at gromacs.org.
> >
> --
> Gromacs Users mailing list
>
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
More information about the gromacs.org_gmx-users
mailing list