Hi!

Hope you've had a fantastic weekend. I'm running 2D-REM (replica-exchange) with REST (solute tempering). During the production run, CHARMM was repeatedly killed by "ENERGY CHANGE TOLERANCE EXCEEDED" error. As I know, exceeding the energy change tolerance can indicate model building errors or an incorrect dynamics set-up. But I didn't think there was anything wrong with my system, then I suspected something should be wrong with my dynamic settings. I increased "echeck" to "999.0" (as suggested in a thread in the forum) but got the same error. I'm just wondering whether I missed any points. My scripts are:

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
! Set variables
! -------------
...
set exfrq = 500 ! frequency
set ftemp = 300 ! targeted temperature
...

! Nonbond options
! ---------------
set ctonnb = 10.0 ! the radius at which switch function begins to take effect
set ctofnb = 12.0 ! cut-off radius beyond which interactions are ignored
set cutnb = 14.0 ! cut-off radius for non-bonded list
set cutim = 14.0 ! cut-off radius for image list

! GBSW options
! ------------
set sgamma = 0.000 ! nonpolar surface tension coefficients in kcal/(molxA^2)
set nang = 50 ! number of angular integration points
set conc = 0.000 ! salt concentration in M

! Read topology and parameter files
! ---------------------------------
...

! Read PSF and CRD
! ----------------
...

! Estimate Tmass and Pmass from total system mass
! -----------------------------------------------
scalar mass stat
calc tmass = int ( ?stot / 5.0 )
calc pmass = int ( ?stot / 50.0 ) ! CPT ensemble

! Setup PBC
! ---------
stream @struc.pbc
CRYSTAL DEFINE @XTLtype @A @B @C @alpha @beta @gamma
crystal build cutoff @cutim

! Define selections
! -----------------
define water sele resn TIP3 end
define ion sele resn SOD .or. resn CLA end
define solv sele water .or. ion end
define solu sele .not. solv end
define heavy sele solu .and. ( .not. hydrogen ) end

! Turn-on Image Atoms
! -------------------
image byres xcen 0.0 ycen 0.0 zcen 0.0 sele solv end
image byseg xcen 0.0 ycen 0.0 zcen 0.0 sele solu end

! Turn on faster options and shake
! --------------------------------
faster on ! use fast routines
shake bonh param ! all bonds fixed, shake bond distances read from parameter file

! GBSW
! ----
stream @gbradi ! read in GBSW atomic radius
gbsw hybrid sgamma @sgamma nang @nang conc @conc - ! invoke GBSW for HpHMD use
sele solu end ! include all solute
energy ! call energy

! PHMD
! ----
open unit 23 read form name @phmdinp
PHMD PAR 23 wri 25 PH 7.0 NPRI @exfrq PHFRQ 10 mass 10.0 -
BETA 5.0 TEMP @ftemp LAM -
sele resn ASP .or. resn GLU .or. resn HSP .or. resn LYS .or. resn ARG end

! pH/Temperature replica exchange
! -------------------------------
repd nrep 64 exlm phmd expt temp 300.0 temp 300.0 temp 300.0 temp 300.0 temp 300.0 temp 300.0 temp 300.0 temp 300.0 -
temp 327.1 temp 327.1 temp 327.1 temp 327.1 temp 327.1 temp 327.1 temp 327.1 temp 327.1 -
temp 356.7 temp 356.7 temp 356.7 temp 356.7 temp 356.7 temp 356.7 temp 356.7 temp 356.7 -
temp 389.0 temp 389.0 temp 389.0 temp 389.0 temp 389.0 temp 389.0 temp 389.0 temp 389.0 -
temp 424.2 temp 424.2 temp 424.2 temp 424.2 temp 424.2 temp 424.2 temp 424.2 temp 424.2 -
temp 462.5 temp 462.5 temp 462.5 temp 462.5 temp 462.5 temp 462.5 temp 462.5 temp 462.5 -
temp 504.4 temp 504.4 temp 504.4 temp 504.4 temp 504.4 temp 504.4 temp 504.4 temp 504.4 -
temp 550.0 temp 550.0 temp 550.0 temp 550.0 temp 550.0 temp 550.0 temp 550.0 temp 550.0 -
nrpt 8 freq @exfrq unit 17
open write unit 17 card name stage@stage/@in.phmd.log ! output for each
open write unit 1 card name stage@stage/@in.phmd.out ! replica individually
outu 1
stream hphrex_2dp.cmd ! input for each replica (phmd.cmd_0 ... phmd.cmd_64-1)

stop
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%hphrex_2dp.cmd%%%%%%%%%%%%%%%%%%%%%%%%%%
! Set variables
! -------------
set ph = REPPH
set temp = REPTEMP ! effective temperature

! Reset pH conditions for PHMD of each replica
! --------------------------------------------
PHMD RESPH OLDPH 7.0 NEWPH @ph PKATEMP @ftemp

! Invoke BLOCK
! ------------
calc scalepp = @ftemp / @temp ! scaling coefficient of E_{pp}
calc scalepw = @scalepp ** 0.5 ! scaling coefficient of E_{pw}

block 2
call 2 sele solu end ! solute removed from block 1 and assigned to block 2; block 1 now only has solvent
coeff 1 1 1.0 ! scale E_{ww}
coeff 1 2 @scalepw ! scale E_{pw}
coeff 2 2 @scalepp ! scale E_{pp}
end
energy ! call energy

! Run dynamics
! ------------
if @lstage .eq. 0 then
open read unit 11 card name stage@{lstage}_4/@in.ph@ph.temp@temp.rst
else
open read unit 11 card name stage@lstage/@in.ph@ph.temp@temp.rst
endif

open write unit 12 card name stage@stage/@in.ph@ph.temp@temp.rst
open write unit 13 unform name stage@stage/@in.ph@ph.temp@temp.dcd
open write unit 14 form name stage@stage/@in.ph@ph.temp@temp.ene
open write unit 25 form name stage@stage/@in.ph@ph.temp@temp.lamb

dyna cpt leap restart -
nstep @ntotp timestp @lstep ISEED S1 S2 S3 S4 -
iunread 11 iunwri 12 iuncrd 13 kunit 14 iunvel -1 -
nprint @exfrq iprfrq @exfrq ntrfrq @exfrq -
nsavc @exfrq nsavv 0 isvfrq @exfrq -
iasors 0 iasvel 1 iscvel 0 ichecw 0 -
firstt @ftemp finalt @ftemp tstruc @ftemp -
twindl 0.0 twindh 0.0 -
echeck 999.0 -
pconstant pmass @pmass pref 1 pgamma 20 tbath @ftemp -
hoover reft @ftemp tmass @tmass -
elec atom cdie vdw vatom vswitch -
ctonnb @ctonnb ctofnb @ctofnb cutnb @cutnb cutim @cutim -
inbfrq -1 ihbfrq -1 imgfrq -1 -
ewald pmew fftx 64 ffty 64 fftz 64 kappa 0.34 spline order 6

stop
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Note: REPPH & REPTEMP will be replaced with different values for each replica.


Hope someone can help me to figure out any possible errors in the scripts. Thanks in advance.

Best

Last edited by shane yue; 08/04/14 12:06 PM.

Zhi (Shane) Yue
Department of Chemistry
University of Chicago