[gmx-users] refcoord_scaling

[khuchtumur bumerdene]

Hi,
I have a question about refcoord_scaling option and its role in pressure
coupling. What I'm currently doing is just running the lysozyme tutorial on
my own protein.
I've so far had successful runs when equilibrating on my own machine, which
is currently running gmx4.5.3. NPT came to a stable average of ~1.0 bar.
But when I do it on the cluster with gmx4.5.5, grompp gives the following
error:

WARNING 1 [file npt.mdp]:
  You are using pressure coupling with absolute position restraints, this
  will give artifacts. Use the refcoord_scaling option.

My npt.mdp file is as follows (straight from the tutorial at this stage):

define          = -DPOSRES      ; position restrain the protein
; Run parameters
integrator      = md            ; leap-frog integrator
nsteps          = 150000         ; 2 * 150000 = 300 ps
dt              = 0.002         ; 2 fs
; Output control
nstxout         = 100           ; save coordinates every 0.2 ps
nstvout         = 100           ; save velocities every 0.2 ps
nstenergy       = 100           ; save energies every 0.2 ps
nstlog          = 100           ; update log file every 0.2 ps
; Bond parameters
continuation    = yes           ; Restarting after NVT
constraint_algorithm = lincs    ; holonomic constraints
constraints     = all-bonds     ; all bonds (even heavy atom-H bonds)
constrained
lincs_iter      = 1             ; accuracy of LINCS
lincs_order     = 4             ; also related to accuracy
; Neighborsearching
ns_type         = grid          ; search neighboring grid cells
nstlist         = 5             ; 10 fs
rlist           = 1.0           ; short-range neighborlist cutoff (in nm)
rcoulomb        = 1.0           ; short-range electrostatic cutoff (in nm)
rvdw            = 1.0           ; short-range van der Waals cutoff (in nm)
; Electrostatics
coulombtype     = PME           ; Particle Mesh Ewald for long-range
electrostatics
pme_order       = 4             ; cubic interpolation
fourierspacing  = 0.16          ; grid spacing for FFT
; Temperature coupling is on
tcoupl          = V-rescale     ; modified Berendsen thermostat
tc-grps         = Protein Non-Protein   ; two coupling groups - more
accurate
tau_t           = 0.1   0.1     ; time constant, in ps
ref_t           = 300   300     ; reference temperature, one for each
group, in K
; Pressure coupling is on
pcoupl          = Parrinello-Rahman     ; Pressure coupling on in NPT
pcoupltype      = isotropic     ; uniform scaling of box vectors
tau_p           = 2.0           ; time constant, in ps
ref_p           = 1.0           ; reference pressure, in bar
compressibility = 4.5e-5        ; isothermal compressibility of water,
bar^-1
; Periodic boundary conditions
pbc             = xyz           ; 3-D PBC
; Dispersion correction
DispCorr        = EnerPres      ; account for cut-off vdW scheme
; Velocity generation
gen_vel         = no            ; Velocity generation is off


The only change I made is the nsteps as my box is bigger and thus has more
waters in it.

So after reading the warning I tried to find some info on the
refcoord_scaling, which the only information found by my feeble attempts
were from the mdp options and the manual, which states:

*refcoord_scaling:* *no*The reference coordinates for position restraints
are not modified. Note that with this option the virial and pressure will
depend on the absolute positions of the reference coordinates.*all*The
reference coordinates are scaled with the scaling matrix of the pressure
coupling.*com*Scale the center of mass of the reference coordinates with
the scaling matrix of the pressure coupling. The vectors of each reference
coordinate to the center of mass are not scaled. Only one COM is used, even
when there are multiple molecules with position restraints. For calculating
the COM of the reference coordinates in the starting configuration,
periodic boundary conditions are not taken into account. I don't quite
understand this completely, is the option there to change the coordinates
of the solute as the box size changes, where the 'all' option changes the
coordinates of every atom and the com only changes the coordinates of the
solute com? Please correct me on this.
So before posting this question I ran a few npt equilibrations with the all
and com options chosen with parrinello-rahman from 100 ps, to 300 ps, and
to 1000 ps. All with 10 different initial velocities from the nvt. The
reason for increasing time is that the equilibration did not reach the
target 1 bar and ended at averages of 4-7 bars for the 10 equilibrations.
The oscillation was comparable to the graph seen on the lysozyme tutorial
between 0-400.
So I tried running berendsen instead of parrinello rahman as the manual
says it should help me get to the target pressure better and I can follow
up with parrinello-rahman. But the pressure still didn't reach the target
pressure.

So the questions are,
- is my naive understanding of the refcoord_scaling correct? if not will
you help explain it to me or point me to a link that could help me?
- Also, when pressure oscillations are at 100s of bars, is 3-6 bar
difference from 1 bar going to make any difference? (it should. shouldn't
it?)
- Or have I made some other careless mistake?

If you need more information, I'd be happy to provide them.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20111107/80ed3e4c/attachment.html>