[gmx-users] Pot'l energy difference between runs

[p buscemi]

Dear users,

I've run a water membrane model through minimization, nvt, and npt with pressure coupling = surface tension with compression 4.5e-5 4 0 and with 4.5e-5 4.5e-5. I wanted to increase the box size in the z direction to created an air ( really vacuum ) interface. The membrane is in the middle of the box so there are two interfaces.
After the three runs, gmx energy tells me the energy of the system is
Coulomb-14 -1.15883e+06 160 495.912 318.314 (kJ/mol)
Potential -1.9882e+06 7000 14892.4 -45168.4 (kJ/mol)
Kinetic En. 592281 45 2027.77 305.09 (kJ/mol)
Total Energy -1.39592e+06 7000 15068.9 -44863.3 (kJ/mol)

But then with or without increasing the box size and running the output of the above model with the the same npt used above or with nvt and with velocity continuation, an error is produced:
Step 0: The total potential energy is 8.43672e+15, which is extremely high.
The LJ and electrostatic contributions to the energy are 8.43672e+15 and
-192065, respectively. A very high potential energy can be caused by
overlapping interactions in bonded interactions or very large coordinate
values. Usually this is caused by a badly- or non-equilibrated initial
configuration, incorrect interactions or parameters in the topology.

So, the npt.mdp file which works with or without z compression the first time but not the second is given below
I'll try longer eq time, but why would the pot'l energy output change between runs ?
Thanks
Paul Buscemi
BICB, UMN

==============================
define = -DPOSRES ; position restrain the protein
; Run parameters
integrator = sd
nsteps = 500000 ; 2 * 50000 = 100 ps
dt = 0.002;

; Output contro
nstxout = 5000
nstvout = 5000
nstenergy = 5000
nstlog = 5000

; Bond parameters
continuation = yes
constraint_algorithm = lincs
constraints = h-bonds
lincs_iter = 1
lincs_order = 4

; Neighborsearching
cutoff-scheme = Verlet
ns_type = grid
nstlist = 20.0
rvdw = 1.0
vdw-type = cut-off
verlet-buffer-tolerance = 1.4e-4

; Electrostatics
coulombtype = PME
pme_order = 4
fourierspacing = 0.12

; Temperature coupling is on
tcoupl = V-rescale
tc-grps = H2O MEM DrG
tau_t = 0.1 0.1 .1
ref_t = 300 300 300

; Pressure coupling is on
pcoupl = Berendsen ; Parrinello-Rahman
pcoupltype = surface-tension
tau_p = 2.0
ref_p = 1.0 1.0
compressibility = 4.5e-5 4.5e-5 ; or 4.5e-5 0
refcoord_scaling = com

; Periodic boundary conditions
pbc = xyz ; 3-D PBC
; Dispersion correction
DispCorr = EnerPres ; account for cut-off vdW scheme
; Velocity generation
gen_vel = no ; no = Velocity generation is off *******************************
gen-temp = 200