[gmx-users] Problems with coupling groups in gas phase simulations

[Sascha Rehm]

Dear Gromacs-User (and Developers),

I tried to set up and simulate some systems with 1000 Argon atoms and 
roughly a thousand water molecules that were randomly distributed (like 
gaseous water in argon gas).
The starting systems are set up with OPLS/AA and 100nm box side length, 
so the systems are much too large and are expected to shrink in the 
first few picoseconds.

Two different coupling groups were used (Ar and SOL) with two different 
compressibilities (1.27e-2 for Argon and 5.3e-5 for Water).
The systems are heated up in the first ns to 318K and then continued for 
another 100 ns. During that time the boxes shrink only to a box size of 
99.5 nm, leading to a totally wrong density.

When I set up the systems with only one couling group ("System") the 
systems shrink to roughly 38nm and show a correct density in a few 
hundred picoseconds but then the density starts to oszilate every few 
ns. I guess this happens because the compressibility is set to 1.2e-2 
and the water starts to form droplets (which should have a 
compressibility of 5.3e-5)? Increasing tau_p from 0.6 to 2.0 does not 
show big effects on the oszilation.

1. Is there an explanation why the systems behave strange when using 2 
coupling groups? My search on the mailing list did not lead to an answer 
(also, there are only a few topics on simulating gas or condensing water 
droplets)

2. Are there reasons, why I should not use only one coupling group for 
the whole system? E.g. in protein-solvent simulations, the manual says 
that using only one coulpling group leads to an increase of temperature 
in the protein and a decrease in the solvent.


Hopefully someone can help me out of that problem...

Best Regards,

Sascha Rehm

---------------------------------------------------------
The used MDP file looks like the following:

integrator		= md
dt			 = 0.002
nsteps  		 = 5000000

annealing = single
annealing_npoints = 3
annealing_time = 0 100 5000000
annealing_temp = 5 318.15 318.15

nstxout 		 = 0
nstvout 		 = 0

nstenergy		 = 1000
nstxtcout		 = 1000

;group(s) to write to energy file
energygrps		 = system

nstlist 		 = 10
ns_type 		 = grid

;cut-off distance for the short-range neighbor list
rlist			 = 1.5

;treatment of electrostatic interactions
coulombtype		 = Cut-off
rcoulomb		 = 1.5

;treatment of van der waals interactions
vdwtype			 = cut-off
rvdw			 = 1.5

; Periodic boudary conditions in all the directions
pbc                      = xyz

;Temperature coupling
tcoupl  		 = berendsen
tc-grps 		 = system
tau_t			 = 0.1
ref_t			 = 318.15

;Pressure coupling
Pcoupl  		 = berendsen
Pcoupltype               = isotropic
tau_p			 = 0.6
compressibility 	 = 1.27e-2
ref_p			 = 1.0

;Velocity generation
gen_vel 		 = yes
gen_temp		 = 5
gen_seed		 = 173529





-- 
Dipl. Bioinf. Sascha Rehm
Institut fuer Technische Biochemie
Allmandring 31
D-70569 Stuttgart
http://www.itb.uni-stuttgart.de/