Re: Pressure Discrepancy

From: Sterling Paramore (paramore_at_hec.utah.edu)
Date: Thu Mar 29 2007 - 13:39:15 CDT

The pressure autocorrelation function you would want to calculate is
<P(t)P(0)> and the decorrelation time would be related to how long it
takes the system to approach a constant value. Allen and Tildesley's
"Computer Simulation of Liquids" has an excellent discussion of what
correlation functions are and even provides a bit of fortran code that
could be used to calculate one. Most stat. mech. books (e.g., McQuarrie
or Chandler) describe what correlation functions are, but don't
necessarily explain how you calculate one from a simulation. I also
imagine that there a number of freely-downloadable utilities you could
find on the web if you're not comfortable coding it yourself.

-Sterling

Morad Alawneh wrote:

> *Regarding the equilibration procedure here is a table that explain
> the 25 ns preparation, before running 10 ns production under NPAT
> ensemble:
> |---------- kcal/mol/A^2 --------|
> *
> *
> * *
> * *Protein* *Lipids* *COM* *Steps*
> *
> * *
> * *non backbone* *backbone* *P* *protein* *small* *large* *huge*
> *min1* *NVT* *fixed* *fixed* *100.0* *
> * *3,500* *5,000* *15,000*
> *heat1* *NV* *fixed* *fixed* *100.0* *
> * *20,000* *20,000* *20,000*
> *eq1* *NVT* *fixed* *fixed* *100.0* *
> * *25,000* *25,000* *25,000*
> *annealing* *NVT* *fixed* *fixed* *100.0* *
> * *200,000* *200,000* *200,000*
> *eq2* *NPAT* *fixed* *fixed* *10.0* *
> * *100,000* *100,000* *100,000*
> *eq3* *NpgT* *fixed* *fixed* *5.0* *
> * *10,000,000* *10,000,000* *15,000,000*
> *eq4* *NPAT* *1.0* *50.0* *1.0* *
> * *100,000* *100,000* *100,000*
> *eq5* *NPAT* *0.5* *20.0* *0.5* *
> * *100,000* *100,000* *100,000*
> *eq6* *NPAT* *0.0* *5.0* *0.0* *
> * *100,000* *100,000* *100,000*
> *eq7* *NPAT* *0.0* *1.0* *0.0* *
> * *100,000* *100,000* *100,000*
> *eq8* *NPAT* *0.0* *0.2* *0.0* *
> * *100,000* *100,000* *100,000*
> *eq9* *NPAT* *0.0* *0.05* *0.0* *
> * *401,500* *400,000* *390,000*
> *eq10* *NPAT* *0.0* *0.00* *0.0* *1000.0* *1,250,000*
> *1,250,000* *1,250,000*
> *
> * *
> * *
> * *
> * *
> * *
> * *
> * *
> * *
> *
> *Total ns* *
> * *
> * *
> * *
> * *
> * *25.00* *25.00* *35.00*
>
>
> *I have many systems with sizes ~6,000, 13,000, 33,000, or 86,000 atoms.
>
> Sterling, would you explain to me how to do pressure autocorrelation
> time calculations.
>
> Another possibility from our discussion is that I am using Multiple
> Timesteps, and yet no one confirmed if that could have an effect or
> not. My understanding that by activating that method it will only
> speed up the process of the integration.
>
> By the way, I had a run for 2.5 ns using Peter's suggestion, and the
> results of that run have the same problem as those who use 500 fs decay.
> **LangevinPistonPeriod 200.0 ;# fs
> LangevinPistonDecay 50.0 ;# fs*
> *
> Bellow my signature I am posting part of my namd script, just to make
> sure that every thing is setup correctly for NPAT ensemble.
>
> Thanks
>
> Morad Alawneh
>
> set Numer_of_Steps 5000000
> set Temp 303.15
> ####################################################################################
> # Input Files
> structure "../../psfgen/$InputName.psf-old_TRP"
> coordinates "../../psfgen/$InputName.pdb-old_TRP"
> bincoordinates [concat "../$InputName\_$Process0\.coor"]
> binvelocities [concat "../$InputName\_$Process0\.vel"]
> parameters "../../psfgen/par_all27_prot_lipid.prm-old_TRP"
> paraTypeCharmm on
> ####################################################################################
> # Output Files
> outputname $OutputName
> binaryoutput yes
> DCDfile $OutputName.dcd
> DCDfreq 1000
> ####################################################################################
> # Timestep Parameters
> numsteps $Numer_of_Steps
> timestep 2.0
> firsttimestep $FirstTimeStep
> ####################################################################################
> # Simulation Space Partitioning
> switching on
> switchDist 8.0
> cutoff 12.0
> pairlistdist 14.0
> margin 5.0
> pairlistsPerCycle 4
> ####################################################################################
> # Basic Dynamics
> exclude scaled1-4
> 1-4scaling 1.0
> seed 2006
> rigidBonds all
> zeroMomentum yes
> ####################################################################################
> # PME Parameters
> PME yes
> PMEGridSizeX 64
> PMEGridSizeY 64
> PMEGridSizeZ 256
> PMEInterpOrder 6
> ####################################################################################
> # Multiple Timestep Parameters
> fullElectFrequency 2
> ####################################################################################
> # Temperature Control
>
> # Langevin Dynamics Parameters
> langevin on
> langevinTemp $Temp
> langevinDamping 5.0
> langevinHydrogen on
> ####################################################################################
> # Boundary Conditions
>
> # Periodic Boundary Conditions
> cellBasisVector1 40.0 00.0 00.0
> cellBasisVector2 00.0 40.0 00.0
> cellBasisVector3 00.0 00.0 90.0
> cellOrigin 0.0 0.0 0.0 ;# default
> wrapWater on
> wrapAll on
> wrapNearest off
> extendedSystem $ESF
> ####################################################################################
> # Pressure Control
> useGroupPressure yes
> useFlexibleCell yes
> useConstantRatio no
> useConstantArea yes
>
> # Nose-Hoover Langevin Piston Pressure Control
> LangevinPiston on
> LangevinPistonTarget 1.01325 ;# 1 atm 1.01325 bar
> LangevinPistonPeriod 200.0 ;# fs
> LangevinPistonDecay 50.0 ;# fs
> LangevinPistonTemp $Temp*

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:44:31 CST