Dear Gromacs users,
I am having a very strange problem with "Nose Hoover thermostat with
Parrinello-Rahman barostat" NPT simulations for a system of an ion of
charge +2 in flexible water molecules. Flexible water is taken from J.
Chem. Phys. 124, 024503 (2006); http://dx.doi.org/10.1063/1.2136877.
The ion is [UO2]2+ with charge on U=2.5 and each O has -0.25.
Atoms very soon get frozen after the simulation starts and they remain
frozen, they do not move at all. Only simulation box oscillates, which
causes atoms to oscillate little bit but they do not move at all.
Pressure, temperature etc. all get messed up.
Here is some output from g_energy:
"
Statistics over 41 steps [ 0. through 400. ps ], 3 data sets
All statistics are over 40001 points
Energy Average Err.Est. RMSD Tot-Drift
---
Potential -23648.5 202891.04 -117.711 (kJ/mol)
Temperature 306.6840.1108.192 -0.566776 (K)
Pressure538.6254.489041.8 79.925 (bar)
"
I am using this "Nose Hoover thermostat with Parrinello-Rahman
barostat" after a well equilibration configuration obtained with
Berendsen's thermostat and barostat.
This problem is happening only when Parrinello-Rahman barostat used
with Nose Hoover thermostat.
I also noted the following observations:
If I change either one of the thermostat or barostat with Berendsen,
frozen problem disappears.
If I change flexible water to rigid water, problem disappears.
If my system has no ion, no problem happens.
If I play with cut offs, rlist, then also sometimes problem disappears.
The problem is independent on the size of box and number of water molecules.
As expected, Diffusion constant is very close to zero.
Here is my input MDP file:
integrator = md
dt = 0.001
tinit = 0
nsteps = 40
nstcomm = 10
pbc = xyz
;constraint_algorithm= lincs
;constraints = all-bonds
continuation= yes ; continuing from equilibrated NPT
nstxout = 200
nstvout = 5
nstfout = 0
nstlog = 5000
nstenergy= 2000
nstxtcout= 0
xtc-precision= 1000
nstlist = 10
ns_type = grid
rlist= 1.1
coulombtype = PME_switch
rcoulomb_switch = 0.6
rcoulomb = 0.9
vdw-type = switch
rvdw-switch = 0.6
rvdw = 0.9
DispCorr= EnerPres
fourierspacing = 0.12
fourier_nx = 0
fourier_ny = 0
fourier_nz = 0
pme_order = 6
ewald_rtol = 1e-6
optimize_fft= yes
tcoupl = nose-hoover
tc-grps = System
tau_t = 2
ref_t = 298.15
pcoupl = Parrinello-Rahman
pcoupltype = isotropic
tau_p = 2.0 ; time constant, in ps
ref_p = 1.0 ; reference pressure (in bar)
compressibility = 4.5e-5; isothermal compressibility, bar^-1
gen_vel = no
Any help will be highly appreciated. Please let me know if someone
needs any additional information to solve this issue.
Thanks,
Surya Prakash Tiwari
--
gmx-users mailing listgmx-users@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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists