[MMTK] Barostat not working properly (Argon example)

Konrad Hinsen hinsen at cnrs-orleans.fr
Tue Aug 21 11:44:38 UTC 2007

On 21.08.2007, at 11:14, Andreas Kring wrote:

> There are no problem running the argon.py script - it generates a
> argon_npt.nc file. When I extract the pressure (and the  
> temperature) as
> a function of time from the trajectory in this file, I get the  
> following
> result:
> http://www.fys.ku.dk/~akring/argon_NPT_864_PT_10fs.png
> As seen in this plot the pressure it not constant even though a  
> barostat
> object is added to the universe in the argon.py file? Why is this  
> the case?

The quantity called simply "pressure" in the trajectories should more  
precisely be called "instantaneous pressure". It is a microscopic  
quantity, not a thermodynamic one. It fluctuates wildly and even  
becomes negative. Only its time average over times much longer than  
the barostat relaxation time can be interpreted as a thermodynamic  

Note also that neither the Nosé thermostat nor the Andersen barostat  
keeps any quantity of the simulated system strictly constant. The  
fluctuations become negligible only for sufficienly large systems,  
but that limit is in practice never reached in MD simulations. There  
are other techniques (i.e. the Gaussian thermostat) that do keep the  
temperature strictly constant, but they do not reproduce any standard  
thermodynamic ensemble.

> I extract the pressure values (in units of atm) in the following way
> (pressure value number n from the trajectory in this example):
> tr = Trajectory(None,'argon_npt.nc')
> P[n] = tr[n]['pressure']/Units.atm

This is perfectly correct. You can also retrieve the whole pressure  
time series at once using

	p = tr.pressure/Units.atm

Konrad Hinsen
Centre de Biophysique Moléculaire, CNRS Orléans
Synchrotron Soleil - Division Expériences
Saint Aubin - BP 48
91192 Gif sur Yvette Cedex, France
Tel. +33-1 69 35 97 15
E-Mail: hinsen at cnrs-orleans.fr

More information about the mmtk mailing list