[gmx-users] Negative average pressure value after NPT equilibration

[Mark Abraham]

Hi,

Negative pressure means the system wants to contract. But your simulations
are probably too short to get a reliable measurement of pressure, nor to
equilibrate it if you were not already very close to the right density. See
background at e.g. http://www.gromacs.org/Documentation/Terminology/Pressure
 and
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/lysozyme/07_equil2.html


Mark

On Sat, Oct 22, 2016 at 1:34 PM Dario Akaberi <dario.aka17 at gmail.com> wrote:

> Dear users
>
> I am trying to perform the MD simulation of a protein in water using the
> amberFF99sb-ildn force field.
>
> To so so I created the protein topology specifying TIP3P water. then I
> created a dodecahedron box with -d 1.0, I filled it using spc216.gro. and
> added 7 NA ions. This resulted in a system of  8837 water molecules plus
> the 7 ions and the protein. Minimization were carried out using conjugate
> gradient method (gromacs is compiled with double precision) which converged
> to emtol 0.01 in 48000 steps.The following nvt ensemble worked well (I
> think) and the resulting average temperature was 299.702 (ref temperature
> was 300 K) with Err. Est. = 0.32  RMSD 4.39922 and tot-Drift 1.58106.
> Now the problems begin. After NVT equilibration I performed an NPT
> equilibration step of 100 ps which resulted in an average pressure of -50.
> I then tried to equilibrate for a longer time: 150 and 200 ps with the same
> result. I also tried a 500 ps NPT equilibration which resulted in an
> average  pressure of -60 and average density 998.935 with Err. Est = 0.12,
> RMSD = 2.38 and Tot-Drift = 0.83 /kg/m^3).
>
> This is the NPT.mdp file I used:
>
> define                    = -DPOSRES
> integrator               = md        ; leap-frog integrator
> nsteps                   = 50000             ; 2 * 50000 = 100 ps
> dt                          = 0.002        ; 2 fs
>
> nstxout                   = 500        ; save coordinates every 1.0 ps
> nstvout                   = 500        ; save velocities every 1.0 ps
> nstenergy               = 500        ; save energies every 1.0 ps
> nstlog                     = 500        ; update log file every 1.0 ps
>
> continuation                = yes        ; Restarting after NVT
> constraint_algorithm        = lincs            ; holonomic constraints
> constraints                = all-bonds            ; all bonds (even heavy
> atom-H bonds) constrained
> lincs_iter                = 1                ; accuracy of LINCS
> lincs_order                = 4                ; also related to accuracy
>
> cutoff-scheme              = Verlet
> ns_type                   = grid
> nstlist                   = 10
> rlist                       = 1.0
> rcoulomb                = 1.0
> rvdw                    = 1.0
>
> coulombtype                = PME
> pme_order                = 4
> fourierspacing                = 0.16
>
> ; Temperature coupling is on
> tcoupl                = V-rescale                    ; modified Berendsen
> thermostat
> tc-grps                    = Protein Non-Protein
> tau_t                    = 0.1      0.1
> ref_t                    = 300       300            ; reference
> temperature, one for each group, in K
>
> pcoupl                    = Parrinello-Rahman            ; Pressure
> coupling on in NPT
> pcoupltype                = isotropic                    ; uniform scaling
> of box vectors
> tau_p                    = 2.0                ; time constant, in ps
> ref_p                    = 1.0                ; reference pressure, in bar
> compressibility             = 4.5e-5                    ; isothermal
> compressibility of water, bar^-1
> refcoord_scaling            = com                          * I get an error
> without this
>
> ; Periodic boundary conditions
> pbc                    = xyz                ; 3-D PBC
> ; Dispersion correction
> DispCorr                = EnerPres                    ; account for cut-off
> vdW scheme
> ; Velocity generation
> gen_vel                    = no                ; Velocity generation is off
>
> I don't know what I'm doing wrong, and unfortunately, I am new to molecular
> dynamics simulation which doesn't help. I tried to find a solution and
> found some other cases where people had a negative pressure after NPT
> equilibration step but it was like an average pressure of -2 not -60!
> Any help would be very appreciated.
>
> Thanks,
>
> Dario
> --
> Gromacs Users mailing list
>
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
>