[gmx-users] pressure coupling

[Justin A. Lemkul]

Shuangxing Dai wrote:
> Hi, all,
>      I am trying to do anisotropic coupling using Parrinello-Rahman. I 
> want to reach the equilibrium state at 1000K, 0.1MPa( 1  bar). However, 
> I find that the fluctuation is so large ( p=80.22 bar, T=1007.23 K, 
> fluctuation is 626 bar for pressure and 11.21 K for temperature.)  I run 
> energy minimization first, then Berenderson pressure coupling for 50 ps, 
> finally Parrinello-Rahman pressure coupling for 50 ps. Anyone has 
> experience on anisotropic coupling and know why? 
>   My system has just 4000 atoms. This is my .mdp file:
> 

For a small system, large fluctations in the pressure are to be expected:

http://www.gromacs.org/Documentation/Terminology/Pressure

-Justin

> ;define                   = -DPOSRES
> ; RUN CONTROL PARAMETERS =
> integrator               = sd
> ; start time and timestep in ps =
> tinit                    = 0
> dt                       = 0.001
> nsteps                   = 50000
> ; number of steps for center of mass motion removal =
> nstcomm                  = 100
> ; OUTPUT CONTROL OPTIONS =
> ; Output frequency for coords (x), velocities (v) and forces (f) =
> nstxout                  = 0
> nstvout                  = 0
> nstfout                  = 0
> ; Output frequency for energies to log file and energy file =
> nstlog                   = 10
> nstenergy                = 10
> ; Output frequency and precision for xtc file =
> nstxtcout                = 10
> xtc-precision            = 1000
> ; NEIGHBORSEARCHING PARAMETERS =
> ; nblist update frequency =
> nstlist                  = 20
> ; ns algorithm (simple or grid) =
> ns_type                  = grid
> 
> ;OPTIONS FOR PRESSURE COUPLING
> Pcoupl                   = Parrinello-Rahman
> pcoupltype               = anisotropic
> tau_p                    = 5
> compressibility          = 2.1645e-07  2.1645e-07 2.7322e-07 0 0 0
> ref_p                    = 1 1 1 0 0 0
> ;OPTIONS FOR TEMPERATURE COUPLING
> tc_grps                  = system
> tau_t                    = 1
> ref_t                    = 1000
> ; OPTIONS FOR BONDS     =
> constraints              = hbonds
> ; Type of constraint algorithm =
> constraint-algorithm     = Lincs
> ; Do not constrain the start configuration =
> unconstrained-start      = no
> ; Relative tolerance of shake =
> shake-tol                = 0.0001
> ; Highest order in the expansion of the constraint coupling matrix =
> lincs-order              = 12
> ; Lincs will write a warning to the stderr if in one step a bond =
> ; rotates over more degrees than =
> lincs-warnangle          = 30
> ; Periodic boundary conditions: xyz, no, xy
> pbc                      = xyz
> periodic_molecules       = no
> ; nblist cut-off        
> rlist                    = 1
> 
> ; OPTIONS FOR ELECTROSTATICS AND VDW
> ; Method for doing electrostatics
> coulombtype              = PME
> rcoulomb                 = 1
> ; Method for doing Van der Waals
> vdw-type                 = Cut-off
> ; cut-off lengths       
> rvdw                     = 1
> 
> ; Spacing for the PME/PPPM FFT grid
> fourierspacing           = 0.12
> ; FFT grid size, when a value is 0 fourierspacing will be used
> fourier_nx               = 0
> fourier_ny               = 0
> fourier_nz               = 0
> ; EWALD/PME/PPPM parameters
> pme_order                = 6
> ewald_rtol               = 1e-4
> ewald_geometry           = 3d
> epsilon_surface          = 0
> optimize_fft             = no
> 
> 
> 
> 
> I used the same procedure to reach different state ( like 0.01 K, 100K, 
> 300K, 500K, 800K, 1000K). Only the 0.01 K got the correct average 
> pressure and small fluctuation.
> T 	p(bar) 	delta p 	T (K) 	delta T
> 0.01 	1.03E+00 	2.12E+00 	9.95E-03 	1.21E-04
> 100 	5.48359 	71.5565 	91.198 	1.25218
> 300 	2.62807 	162.926 	290.324 	3.36594
> 500 	16.6356 	1048.76 	492.616 	5.45231
> 800 	59.3245 	594.773 	799.314 	9.18876
> 1000 	80.2257 	626.789 	1007.23 	11.2159
> 
> Thanks in advance,
> Shuangxing Dai
> 

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

Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

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