[gmx-users] pressure coupling

[Shuangxing Dai]

Thanks for the response.
Yes, I know that large fluctuations will be expected. But why the average
was not even correct? Also, how can I reach some state ( some pressre and
temperature) and know/tell that the system is in equilibrium if the average
is not correct with large fluctuation? I mean I do not know the criteria of
Gromacs. Or maybe my procedure is wrong?

Message: 1
Date: Mon, 05 Jul 2010 10:39:30 -0400
From: "Justin A. Lemkul" <jalemkul at vt.edu>
Subject: Re: [gmx-users] pressure coupling
To: Discussion list for GROMACS users <gmx-users at gromacs.org>
Message-ID: <4C31EEA2.7070105 at vt.edu>
Content-Type: text/plain; charset=UTF-8; format=flowed



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
>
Thanks,
Shuangxing Dai
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20100705/70aac7fa/attachment.html>