[gmx-users] maxh not causing termination

[Ben Reynwar]

On Tue, Dec 20, 2011 at 7:16 PM, Mark Abraham <Mark.Abraham at anu.edu.au> wrote:
> On 12/19/2011 1:51 PM, Ben Reynwar wrote:
>>
>> I'm having a problem with gromacs not terminating as expected when
>> using the maxh option.
>> It is an REMD simulation with 32 replicas.
>> I'm specifying -maxh 24 and as expected see the following in the stderr
>> output.
>>
>> Step 773882: Run time exceeded 23.760 hours, will terminate the run
>>
>> Step 773876: Run time exceeded 23.760 hours, will terminate the run
>>
>> Step 773880: Run time exceeded 23.760 hours, will terminate the run
>>
>> etc
>>
>> However I can see that the output files continued to be written for
>> another hour until at 25 hours the simulation was terminated by the
>> queueing system.
>> No checkpoint files were produced.  The output files show that the
>> simulation continued until about step 797000.
>>
>> I've done similar things previously without running into this problem.
>>  Anyone have any ideas for what stupid mistake I could be making?
>
>
> Perhaps none. What GROMACS version is this? Does the latest version have the
> same behaviour?
>
> Mark
>
> --
> This message has been scanned for viruses and
> dangerous content by MailScanner, and is
> believed to be clean.
>
> --
> gmx-users mailing list    gmx-users at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
> Please don't post (un)subscribe requests to the list. Use the www interface
> or send it to gmx-users-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Sorry for the slow response.  I've been away on holiday.

The version used is 4.5.5 with a bug fix applied from
http://lists.gromacs.org/pipermail/gmx-developers/2011-October/005405.html
It was a replica exchange simulation.

I've narrowed down the problem to be specific to using infinite
cutoff.  With a finite cutoff there are no problems.

The initial run works fine and creates checkpoints.  However a run
started from one of these checkpoints fails to create checkpoints
itself.

The diff between the an mdp that works fine and one that doesn't is:

53,57c53,58
< ; [nm] cut-off distance for the short-range neighbor list.
< ; Set to zero for infinite cut-off.
< rlist = 0
< ; Recalculating pair-list is not necessary with an infinite cut-off
< nstlist = 0
---
> ; [nm] cut-off distance for the short-range neighbor list
> ; For Generalized Born this must be equal to the cut-off length for
> ; the born radius calculation.
> rlist = 1.6
> ; [steps] freq to update neighbor list
> nstlist = 1
67,68c68,69
< ; infinite coloumb cut-off radius
< rcoulomb = 0
---
> ; coloumb cutoff radius
> rcoulomb = 1.6
76,77c77,78
< ; Infinite VsW cutoff radius.
< rvdw = 0
---
> ; Increasing VdW cutoff to same as everything else.
> rvdw = 1.6
144c145
< rgbradii = 0
---
> rgbradii = 1.6

The complete offending mdp file will follow.
Also if there's anything stupid in here unrelated to my current
problem, please let me know too :).

; 7.3.2 Preprocssing
; ------------------
; Apply constaint to alpha carbon atoms in alpha-crystalline domain.
define = -DALPHACRYST_POSRES -DCHIRRES

; 7.3.3 Run Control
; -----------------
; group(s) for center of mass motion removal
comm_grps = System
; Do Langevin dynamics.
integrator = sd
; maximum number of steps to integrate
nsteps = 100000000
; remove center of mass translation and rotation around centre of mass
comm_mode = Angular
; [ps] time step for integration
dt = 0.002
; [steps] frequency of mass motion removal
nstcomm = 10
; [ps] starting time for run
tinit = 0

; 7.3.4 Langevin Dynamics
; -----------------------
; Use PID to seed random number.
ld_seed = -1

; 7.3.8 Output Control
; --------------------
; [steps] freq to write velocities to trajectory
nstvout = 0
; [steps] freq to write energies to log file
nstlog = 1000
; [steps] freq to write energies to energy file
nstenergy = 1000
; group(s) to write to xtc trajectory
xtc_grps = System
; [real] precision to write xtc trajectory
xtc_precision = 1000
; [steps] freq to write coordinates to xtc trajectory
nstxtcout = 1000
; [steps] freq to write coordinates to trajectory
nstxout = 0
; group(s) to write to energy file
energygrps = System
; [steps] freq to write forces to trajectory
nstfout = 0

; 7.3.9 Neighbour Searching
; -------------------------
; [nm] cut-off distance for the short-range neighbor list.
; Set to zero for infinite cut-off.
rlist = 0
; Recalculating pair-list is not necessary with an infinite cut-off
nstlist = 0
; method of updating neighbor list
ns_type = grid
; no periodic boundary conditions
pbc = no

; 7.3.10 Electrostatics
; ---------------------
; apply a cut-off to electostatic
coulombtype = cut-off
; infinite coloumb cut-off radius
rcoulomb = 0

; 7.3.11 VdW
; ----------
; Dispersion correction makes no sense without box size.
DispCorr = no
; twin-range cut-off with rlist where rvdw >rlist
vdwtype = cut-off
; Infinite VsW cutoff radius.
rvdw = 0

; 7.3.13 Ewald
; ------------
; [nm] grid spacing for FFT grid when using PME
fourierspacing = 0.12
; relative strength of Ewald-shifted potential at rcoulomb
ewald_rtol = 1e-05
; interpolation order for PME, 4 cubic
pme_order = 4

; 7.3.14 Temperature Coupling
; ---------------------------
; Temperature.
ref_t = 600.0
; temperature coupling frequency
nsttcouple = 1
; Thermostat is already taken care of by langevin dynamics.
tcoupl = no
; couple everything to same bath
tc_grps = System
; [ps] time constant for coupling
tau_t = 1

; 7.3.15 Pressure Coupling
; ------------------------
; [bar] reference pressure for coupling
ref_p = 1.0
; pressure coupling in x-y-z directions
pcoupltype = isotropic
; [ps] time constant for coupling
tau_p = 2.0
; no pressure coupling if using generalized born
pcoupl = no
; [bar^-1] compressibility
compressibility = 4.5e-05

; 7.3.17 Velocity Generation
; --------------------------
; velocity generation turned off
gen_vel = no

; 7.3.18 Bonds
; ------------
; apply constraints to the start configuration
continuation = yes
; All bonds are constraints
constraints = allbonds
; [degrees] maximum angle that a bond can rotate before LINCS will complain
lincs_warnangle = 30
; highest order in the expansion of the contraint coupling matrix
lincs_order = 4
; number of iterations to correct for rotational lengthening
lincs_iter = 1
; LINear Constraint Solver
constraint_algorithm = LINCS

; 7.3.22 NMR refinement (dihedral restraints not in manual)
; ---------------------------------------------------------
; scaling factor for force constant for dihedral constraints
dihre-fc = 10
; Set chirality and omega angle restraints.
dihre = yes

; 7.3.27 Implicit Solvent
; -----------------------
; The cutoff radius for calculating the Born radii (must be equal to rlist).
rgbradii = 0
; Use a Generalized Born Solvent
implicit_solvent = GBSA
; Use the Onufriev-Bashford-Case method to calculate the Born radii
gb_algorithm = OBC
; Dielectric constant for implicit solvent.
; Default is 80 but we use 78.5 for consistency with amber.
gb_epsilon_solvent = 78.5