[gmx-users] Re: gmx-users Digest, Vol 74, Issue 84 (anisotropic pressure coupling)

Justin A. Lemkul jalemkul at vt.edu
Tue Jun 15 03:19:59 CEST 2010



Shuangxing Dai wrote:
> Well, the final average pressure is not correct:
> 
> * **<======  ###############  ==>*
> * **<====  A V E R A G E S  ====>*
> * **<==  ###############  ======>*
> *
> *
> *   Energies (kJ/mol)*
> *  Buck.ham (SR)   Coulomb (SR)   Coul. recip.      Potential    Kinetic 
> En.*
> *    1.60633e+06   -5.69243e+06   -3.45003e+06   -7.53612e+06   
>  1.80719e+04*
> *   Total Energy    Temperature Pressure (bar)*
> *   -7.51805e+06    3.62347e+02    1.11758e+04*
> *
> *
> *
> *
> *   Pressure (bar)*
> *    1.21051e+04   -2.94017e-01    1.30456e+00*
> *   -2.94017e-01    9.03597e+03   -8.10850e-01*
> *    1.30456e+00   -8.10850e-01    1.23863e+04*
> *
> *
> *
> *
> * **<======  ###############################  ==>*
> * **<====  R M S - F L U C T U A T I O N S  ====>*
> * **<==  ###############################  ======>*
> *
> *
> *   Energies (kJ/mol)*
> *  Buck.ham (SR)   Coulomb (SR)   Coul. recip.      Potential    Kinetic 
> En.*
> *    3.92376e+05    4.00326e+05    4.71365e+01    3.02978e+04   
>  2.35039e+03*
> *   Total Energy    Temperature Pressure (bar)*
> *    3.01091e+04    4.71259e+01    1.71117e+05*
> *
> *
> *
> *
> *   Pressure (bar)*
> *    1.82323e+05    4.91513e+02    3.07107e+02*
> *    4.91513e+02    2.02402e+05    3.13857e+02*
> *    3.07107e+02    3.13857e+02    1.30179e+05*
> *
> *
> *Also the fluctuation is very large. I think 200ps is long enough. But 
> still I cannot find the equilibrium.*

Then 200 ps is clearly not enough, or your entire procedure is not correct. 
Note too that your temperature (362 K) is not at all close to what you want (300 
K).  You may be setting tau_t incorrectly; check the archives.  This value, when 
using the sd integrator, is not the same as the normal interpretation of 
relaxation time, so instead of 0.1, you may want 1.0 or larger.

> *I do not know whether the parameters are wrong or I ran the anisotropic 
> pressure coupling wrong. Is *Parrinello-Rahman pressure coupling 
> sensitive to the initial pressure? I mean if the initial pressure ( in 
> order of GPa), object pressure is 1 bar ( 1e5 Pa). 

If your system is far from equilibrium, it may be beneficial to relax it with a 
weak coupling method (i.e., Berendsen) first.  Otherwise, a system that is not 
close to equilibrium may oscillate very wildly and not converge to the desired 
value.

-Justin

> *
> *
> Thanks,
> Shuangxing Dai
> 
> 
> On Mon, Jun 14, 2010 at 8:50 PM, <gmx-users-request at gromacs.org 
> <mailto:gmx-users-request at gromacs.org>> wrote:
> 
> 
>     Message: 3
>     Date: Mon, 14 Jun 2010 19:18:14 -0400
>     From: "Justin A. Lemkul" <jalemkul at vt.edu <mailto:jalemkul at vt.edu>>
>     Subject: Re: [gmx-users] anisotropic pressure coupling
>     To: Discussion list for GROMACS users <gmx-users at gromacs.org
>     <mailto:gmx-users at gromacs.org>>
>     Message-ID: <4C16B8B6.4010005 at vt.edu <mailto:4C16B8B6.4010005 at vt.edu>>
>     Content-Type: text/plain; charset=UTF-8; format=flowed
> 
> 
> 
>     Shuangxing Dai wrote:
>      > Hi all:
>      >     I was using Parrinello-Rahman for pressure coupling and here
>     is the
>      > mdp file:
>      >
>      > ; RUN CONTROL PARAMETERS =
>      > integrator               = sd
>      > ; start time and timestep in ps =
>      > tinit                    = 0
>      > dt                       = 0.001
>      > nsteps                   = 200000
>      > ; number of steps for center of mass motion removal =
>      > nstcomm                  = 100
>      > ; NEIGHBORSEARCHING PARAMETERS =
>      > ; nblist update frequency =
>      > nstlist                  = 50
>      > ; ns algorithm (simple or grid) =
>      > ns_type                  = grid
>      >
>      > ;OPTIONS FOR PRESSURE COUPLING
>      > Pcoupl                   = Parrinello-Rahman
>      > pcoupltype               = anisotropic
>      > tau_p                    = 1
>      > compressibility          = 2.1645e-09  2.1645e-09 2.7322e-09 0 0 0
>      > ref_p                    = 1 1 1 0 0 0
>      > ;OPTIONS FOR TEMPERATURE COUPLING
>      > tc_grps                  = system
>      > tau_t                    = 0.1
>      > ref_t                    = 300
>      > ; 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
>      >
>      >  It seems that the temperature coupling is fast and fluctuation
>     is less
>      > that 1%. But the pressure coupling fails. The pressure fluctuates
>     from
>      > +1e3 bar to -1e3 bar. Anyone can help?
>      >
> 
>     I don't call that failure.  Depending on the size of your system,
>     you can
>     probably expect very wide fluctuations in the pressure.  Is the
>     average correct?
> 
>     http://www.gromacs.org/Documentation/Terminology/Pressure
> 
>     -Justin
> 
>      >
>      > Thanks,
>      > 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 <http://vt.edu> | (540) 231-9080
>     http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
> 
>     ========================================
> 
> 
> 

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

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

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



More information about the gromacs.org_gmx-users mailing list