[gmx-users] Re: pressure_coupling

[Dr. Vitaly Chaban]

>
> I want to keep the pressure at 1.0 bar during the NPT simulation. But
> it is fluctuating around 130 bar. So can anyone please inform me
> whether I have missed any keyword in my .mdp file OR is it because of
> the tau_p which I set 1s 1.0 ps.
> Thanks
>
> The .mdp file is given below
>
>
>
> ; 7.3.3 Run Control
> integrator              = md-vv                    ; md integrator
> tinit                   = 0                     ; [ps] starting time for run
> dt                      = 0.001                 ; [ps] time step for integration
> nsteps                  = 15000000                ; maximum number of
> steps to integrate, 0.001 * 15,00,000 = 15 ns
> nstcomm                 = 1                     ; [steps] frequency of
> mass motion removal
> comm_grps               = Protein Non-Protein   ; group(s) for center
> of mass motion removal
>
> ; 7.3.8 Output Control
> nstxout                 = 10000        ; [steps] freq to write
> coordinates to trajectory
> nstvout                 = 10000        ; [steps] freq to write
> velocities to trajectory
> nstfout                 = 10000        ; [steps] freq to write forces
> to trajectory
> nstlog                  = 1000           ; [steps] freq to write
> energies to log file
> nstenergy               = 1000           ; [steps] freq to write
> energies to energy file
> nstxtcout               = 1000           ; [steps] freq to write
> coordinates to xtc trajectory
> xtc_precision           = 1000          ; [real] precision to write
> xtc trajectory
> xtc_grps                = System        ; group(s) to write to xtc trajectory
> energygrps              = System        ; group(s) to write to energy file
>
> ; 7.3.9 Neighbor Searching
> nstlist                 = 1             ; [steps] freq to update neighbor list
> ns_type                 = grid          ; method of updating neighbor list
> pbc                     = xyz           ; periodic boundary conditions
> in all directions
> rlist                   = 1.2           ; [nm] cut-off distance for
> the short-range neighbor list
>
> nsttcouple              = 1
> nstpcouple              = 1
>
> ; 7.3.10 Electrostatics
> coulombtype             = PME           ; Particle-Mesh Ewald electrostatics
> rcoulomb                = 1.2           ; [nm] distance for Coulomb cut-off
>
> ; 7.3.11 VdW
> vdwtype                 = cut-off       ; twin-range cut-off with
> rlist where rvdw >= rlist
> rvdw                    = 1.2           ; [nm] distance for LJ cut-off
> DispCorr                = EnerPres      ; apply long range dispersion
> corrections for energy
>
> ; 7.3.13 Ewald
> fourierspacing          = 0.12          ; [nm] grid spacing for FFT
> grid when using PME
> pme_order               = 4             ; interpolation order for PME, 4 = cubic
> ewald_rtol              = 1e-5          ; relative strength of
> Ewald-shifted potential at rcoulomb
>
> ; 7.3.14 Temperature Coupling
> tcoupl                  = Nose-Hoover                   ; Nose-Hoover
> temperature coupling
> tc_grps                 = Protein    Non-Protein        ; groups to
> couple seperately to temperature bath
> tau_t                   = 1.0        1.0                ; [ps] time
> constant for coupling
> ref_t                   = 300        300                ; [K]
> reference temperature for coupling
>
> ; 7.3.15 Pressure Coupling
> pcoupl                  = MTTK                  ; pressure coupling
> where box vectors are variable
> pcoupltype              = isotropic             ; pressure coupling in
> x-y-z directions
> tau_p                   = 1.0                   ; [ps] time constant
> for coupling
> compressibility         = 4.5e-5                ; [bar^-1] compressibility
> ref_p                   = 1.0                   ; [bar] reference
> pressure for coupling
>
> ; 7.3.17 Velocity Generation
> gen_vel                 = no            ; velocity generation turned off
>
> ; 7.3.18 Bonds
> constraints             = h-bonds
> constraint_algorithm    = SHAKE          ; SHAKE Constraint Solver
> shake_tol               = 1.0e-5
>


You have a correct MDP file for the constant pressure constant
temperature ensemble, if this makes you happy.

One can enumerate a zillion of reasons why the pressure you get  from
g_energy (I hope, you list an average one) is different from that in
the MDP file. The system can be not properly equilibrated, the
statistics can be insufficient, etc.


Dr. Vitaly V. Chaban
MEMPHYS - Center for Biomembrane Physics
Department of Physics, Chemistry and Pharmacy
University of Southern Denmark
Campusvej 55, 5230 Odense M, Denmark