[gmx-users] Re : NPT Simulation average pressure lower than Reference set pressure

[Justin Lemkul]

On 11/1/14 11:56 PM, Agnivo Gosai wrote:
> Dear Users
>
> I am sharing a few plots through Google Drive. If you are able to view the
> files then kindly comment.
>
> My workflow is as follows :-
> 1) Topolgy preparation with AMBER99SB-Parmbsc0 ff in GROMACS 4.5.6. TIP3P
> water model.
> 2) Putting the complex in a box , using spc216 to solvate , adding ions.
> 3) Energy minimization with steepest descent.
> 4) 100 ps NVT equilibration with V-rescale thermostat. Ref T = 300 K , Avg.
> T = 299.7 K and Total Energy also remained constant.
> 5) 1st NPT simulation with Berendsen barostat (Ref. Pressure = 1 bar)
> having the following results :
> ( A)​Statistics over 50001 steps [ 0.0000 through 100.0000 ps ], 1 data sets
> All statistics are over 10001 points
>
> Energy                      Average   Err.Est.       RMSD  Tot-Drift
> ------------------------------
> Pressure                    -5.5706        6.2    84.8085    38.1219  (bar)
> (B)All statistics are over 10001 points
>
> Energy                      Average   Err.Est.       RMSD  Tot-Drift
> ------------------------------
> Density                     1012.13       0.55    1.87172    3.17412
> (kg/m^3)
> 6) 2nd NPT simulation with Parrinello-Rahman barostat (Ref. Pressure = 1
> bar) having the following results :
> These are the average values :-
>
> Energy                      Average   Err.Est.       RMSD  Tot-Drift
> ------------------------------
> Pressure                  -0.299523       0.64     93.313   -1.18081  (bar)
>
> Energy                      Average   Err.Est.       RMSD  Tot-Drift
> ------------------------------------------------------------
> -------------------
> Density                     1012.32       0.26    1.68576    0.82485
> (kg/m^3)
> 7)3rd NPT simulation with Parrinello-Rahman barostat (Ref. Pressure = 1
> bar) having the following results :
> Statistics over 50001 steps [ 0.0000 through 100.0000 ps ], 1 data sets
> All statistics are over 10001 points
>
> Energy                      Average   Err.Est.       RMSD  Tot-Drift
> -------------------------------------------------
> Pressure                    2.03959        0.9    95.9974     -6.788  (bar)
>
> All statistics are over 10001 points
>
> Energy                      Average   Err.Est.       RMSD  Tot-Drift
> ------------------------------------------------------------
> -------------------
> Density                      1012.8       0.15    1.60294   0.243262
> (kg/m^3)
>
> All the NPT equilibrations are of 100 ps and are continued from the
> previous equilibration.
>
> I observe that the average pressure has not yet reached the reference
> pressure and in the 3rd NPT run the Err.Est. , RMSD and Drift has increased.
>
> I have gone through the different comments and suggestions to my earlier
> posts. I am thinking of a 500 ps or 1 ns run with the Parinello Rahman
> barostat. Or am I failing to understand the physical significance of these
> average values and should have stopped after the 2nd NPT equilibration ?​
>   em_potential.png
> <https://docs.google.com/file/d/0B-U8uULVZjfRZ0xfZGVrLVRYQVU/edit?usp=drive_web>
> ​​
>   npt1_density.png
> <https://docs.google.com/file/d/0B-U8uULVZjfRYUpKeGZfenNqRnM/edit?usp=drive_web>
> ​​
>   npt1_pressure.png
> <https://docs.google.com/file/d/0B-U8uULVZjfRa1FSWkNCTkxkM28/edit?usp=drive_web>
> ​​
>   npt1_temp.png
> <https://docs.google.com/file/d/0B-U8uULVZjfRcmFmdHBybjZRTnM/edit?usp=drive_web>
> ​​
>   npt2density.png
> <https://docs.google.com/file/d/0B-U8uULVZjfRaUZKa1JJNWJCcU0/edit?usp=drive_web>
> ​​
>   npt2_pressure.png
> <https://docs.google.com/file/d/0B-U8uULVZjfRTng2eW1iQUxPamc/edit?usp=drive_web>
> ​​
>   npt3_density.png
> <https://docs.google.com/file/d/0B-U8uULVZjfRUGJSQU9MRnhlUTQ/edit?usp=drive_web>
> ​​
>   npt3_pressure.png
> <https://docs.google.com/file/d/0B-U8uULVZjfRcXFJZkt3N1YtdzQ/edit?usp=drive_web>
> ​​
>   nvt_temp.png
> <https://docs.google.com/file/d/0B-U8uULVZjfRLTNNcmlEcFNSWDg/edit?usp=drive_web>
> ​​
>   nvt_totalenergy.png
> <https://docs.google.com/file/d/0B-U8uULVZjfRNFBZd05kMXNoYW8/edit?usp=drive_web>
> ​
>
> Could anybody suggest why the reference pressure is not being reached ?
>

Pressure is incredibly noisy, and you're going to get a distribution of values 
that varies wildly.  Think generally about the statistics: if you were given 
three averages and standard deviations of -5.5 ± 85, -0.3 ± 93, and 2.0 ± 96, 
would you be able to say that they are significantly different from 1?  Of 
course, you can do that analysis, but it's pretty clear that they're not.  I 
have only seen an average pressure come out to 1.0 probably once so far in my 
career, though unfortunately (?) it is the case with one of my tutorials, so now 
everyone expects it to always happen :)

Net result: your results are fine.  Maybe over very long time periods you'll get 
to 1.0 or very close to it, but there is nothing wrong with your simulations and 
you're getting what you should be.

-Justin

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

Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow

Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 629
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201

jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul

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