[gmx-users] Viscosity calculation problems (periodic perturbation method)

[Justin A. Lemkul]

There have been several bug fixes related to viscosity calculations since the 
release of version 4.5.3, of which at least one (on 1/10/2011) specifically 
mentions zero viscosity in the .edr and .log files.

I would suggest pulling the latest release-4-5-patches from the git repository 
and try again.

-Justin

Xiang Gu wrote:
> Hi, all,
> 
> I was trying to reproduce the viscosity results of SPC water reported in 
> Berk Hess' paper (JCP, vol.116 No.1, 209-217, 2002) by using the 
> periodic perturbation method. After reading some old tips on the 
> maillist exchanged by Song, David, etc., I tried to set up the NVT 
> simulation just according to Hess has done, and sent it to calculation.
> 
> Probably because not everything is properly understood, I processed the 
> edr file with g_energy_d but found both the columns 2CosZ*Vel-X & 
> 1/Viscosity always zero. I played with several options of Tcoupl and 
> tau_t, still the same happened. Could anybody experienced in this tutor 
> me what has been improperly specified in the mdp file (see below; system 
> is 3456 SPC water molecules contained in a 3.75*3.75*7.5 nm^3 box--just 
> as in Hess' paper Table II, already well equilibrated at 300 K in 
> advance; perhaps it is unnecessary to run this double precisely, but I 
> was using double precision version (4.5.3.d) just because I forgot to 
> switch to the single precision one, however this shouldn't cause the 
> problem I had ...)?
> 
> Thanks and wish to see your suggestions soon!
> 
> Xiang Gu
> 
> ;    Generic mdp file for SPC water equilibration
> ;    Gromacs 4.3.x
> ;
> ;    T = 300 K
> ;
> ;    NVT equilibration run, 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
> 
> ; X/V/F/E outputs ---- change these according to your system ----
> 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         ; xyz(default), no, full(infinite 
> systems)
> 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              ; berendsen, tau_p = 1.0 for 
> faster equilibration
> ;Pcoupltype          =  isotropic       ; semi-isotropic: xy and z 
> separately (CNT)
> ;tau_p               =  1.0             ; ps
> ;compressibility     =  4.5e-5          ; 1/bar (water @ 1 atm, 300 K)
> ;ref_p               =  1.0             ; bar0  4000 1.034798 0.001949
> 
> gen_vel                 = no
> ;gen_temp                = 300.0
> ;gen_seed                = 173529
> 

-- 
========================================

Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

========================================