Dear Gromacs users, This question seems to come up periodically in the mailing list, but none of the previous answers seem helpful in my case.
I'm trying to reproduce the viscosity calculation of SPC water by Berk Hess (JCP 116, 2002) using cos_acceleration in Gromacs. The answer I get is 2 orders of magnitude out. My topology file and parameter file is appended at the bottom of this email. I run g_energy and get Energy Average Err.Est. RMSD Tot-Drift ------------------------------------------------------------------------------- 1/Viscosity 23.4689 0.13 2.39126 -0.255371 (m s/kg) which gives a viscosity of 0.04 kg/(m s), or 40 mPa.s The value quoted in the paper is about 0.4 mPa.s which is around the correct value for water, give or take a bit. So my question is, where is my missing factor of 100? Any advice is much appreciated. Thank you. James -------------------------- Topology file: #include "ffgmx.itp" #include "spc.itp" [ system ] Pure Water [ molecules ] SOL 3456 -------------------------- Parameter file for system (3456 SPC water molecules, 3.75x3.75x7.5 nm3 box): ; Generic mdp file for SPC water equilibration ; Gromacs 4.3.x ; ; T = 300 K ; ; NVT 1.2 ns ; define = ; define here posres etc., e.g. -DPOSRES integrator = md tinit = 0 dt = 0.001 nsteps = 1200000 ; Bond constraints continuation = no ; switch to 'yes' if need to read in velocities etc. constraints = none ; constrain all bond lengths constraint_algorithm = lincs ; default lincs_order = 4 ; default nstxout = 1000 nstvout = 1000 nstfout = 0 nstlog = 1000 nstenergy = 1000 ; Output frequency and precision for xtc file nstxtcout = 1000 xtc-precision = 1000 ; This selects the subset of atoms for the xtc file. You can ; select multiple groups. By default all atoms will be written. xtc-grps = ; Selection of energy groups energygrps = System ; Neighbor list ns_type = grid ; neighlist type nstlist = 5 ; Freq. to update neighbour list rlist = 0.9 ; nm (cutoff for short-range NL) pbc = xyz periodic_molecules = no ; Non-equilibrium MD ;acc_grps = SYSTEM cos_acceleration = 0.025 ; PPM option for viscosity calculation (nm/psĀ²) coulombtype = PME rcoulomb = 0.9 optimize_fft = yes ; affects only PME calculations ; if you use PME, set also rcoulomb = rlist ; van der Waals interactions vdwtype = Cut-off ; Van der Waals interactions rvdw = 0.9 ; nm (LJ cut-off) DispCorr = EnerPres ; long-range dispersion correction to energy and pressure Tcoupl = berendsen tc-grps = System tau_t = 2.5 ref_t = 300.0 ;Pressure coupling Pcoupl = no gen_vel = no -- gmx-users mailing list gmx-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