[gmx-users] Re: Viscosity calculation using cos_acceleration
James Cannon
jamesresearching at gmail.com
Tue Apr 16 08:29:36 CEST 2013
Alternatively, has anyone else reproduced the viscosity calculation, or
tried to?
If anyone would be so kind as to forward their input files, that might help
me narrow down the problem with my own input files.
Thank you.
James
On 11 April 2013 16:55, James Cannon <jamesresearching at gmail.com> wrote:
> 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
>
On 11 April 2013 16:55, James Cannon <jamesresearching at gmail.com> wrote:
> 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
>
More information about the gromacs.org_gmx-users
mailing list