Pressure Calculation

From: Priyan Amaras (priyan.amaras_at_yahoo.com)
Date: Fri Sep 02 2011 - 12:43:25 CDT

Dear Jim

We are using NAMD to solve an engineering problem, but we are kind of stuck on something
and we am wondering if you could shed some light on this.
 
My questions have to do with the pressure profile calculations.
 
Firstly, I am confused about the sign: It seems like the sign convention used for the system pressure is: “compressive pressure is positive”, whereas the sign convention for the pressure profile is : “tensile pressure is positive”. Am I right?
 
Just to train myself first, I simulated the behavior of a boxful of TIP3P water molecules in nPT with T=300K, p=1 bar, and n=12000 atoms. The calculated system pressure (time-averaged) is very close to the specified (1.5 bar vs. 1 bar), but the average pressure from the pressure profile tensor (Pave=(Pxx+Pyy+Pzz)/3) was about 700 bars tensile. This is with PME and Ewald for the electrostatic part.
 
We repeated the run without using PME and Ewald for the electrostatic part, and using a cutoff of 10A for van der Walls and 22A for the electrostatic part, the system pressure value was about 1.5bar (compressive), but the average pressure from the pressure profile tensor was still about the same (700 bars tensile).
 
We ran another simulation with p=10000 bars, we got about 10000 bars for the system pressure, but about 9250 bars (compressive) for Pave from pressure profile. There is still a difference of about 750 bars.
 
We noticed that the kinetic part of the pressure profile is always positive (i.e., tensile). I would expect this to be always compressive.
 
In another run, we placed some atoms in a box with some arbitrarily assigned properties (sig=1, eps=1, m=1, T=1, etc.) and ran a simulation under NVT.  The kinetic part of the pressure tensor is still tensile!
 
So, we wonder if there is an error in the sign for the kinetic part in NAMD. Can this explain the 700-750 bar difference we see between the system pressure and Pave from the pressure profile calculations? Only a developer of the code can answer this question. we would really appreciate if you could answer this, or forward this email to the person who developed the pressure profile part of the code. It will save us a lot of time.
Thanks a lot.
Priyan

This archive was generated by hypermail 2.1.6 : Mon Dec 31 2012 - 23:20:46 CST