[gmx-users] non isotropic kinetic energy

[aherz]

Hi,

we took a couple of days to reexamine this and try out a couple of things:

-On small interfacial systems (ca. 2k spc/e water molecules sandwiched
by vacuum)
the difference between ekin_x/y and ekin_z is 20 to 60kJ/mol depending
on the specific system/settings
(e.g. see attached file of averaged ekin (x/y/z) over time).
-On larger interfacial systems the difference grows (up to several
100kJ/mol for ca. 32k spc/e water).

The following settings make no difference (so the error persists even if
we change the setting):
-thermostat (we tryed nose-hoover, berendsen and no thermostat)
-estatics (cut-off instead of pme)
-type of constraint algo for water (so lincs, shake or settle, we
modified the spc/e water molec definition to change this)
-changing system size so that the box length is equal in all 3 dim (as
said, magnifies problem)
-rotating system (error rotates with system)
-use smaller time step of 1fs (no difference)

We ran a simulation with flexible water (so bonds instead of
constraints) , here the error vanishes!
Also, we ran similar sims with amber, where the error does not appear at
all (for constraint water).

We found out about this problem while doing surface tension/contact
angle calculations.
The surf tens calced via virial differs from the surf tens calculated
via pressure tensor (unless one corrects
for the anisotropic kinetic energy). A simple way to see the deficit is
trying to calc the average pressure tensor from the average virial
tensor (using the simple def from the gromacs manual P_ij=2/V(Ekin-Vir_ij):
a) assuming equidistributed kinetic energy (ekin/3 per dim)
b) using the anisotropic ekin which can be extracted from a trr with a
tool (source attached) if velocities have been written to the trr.

All this appears to hint at a problem with the constraints for the
waters and might also hint at an incorrect ensemble as a consequence.

Alex

sample input file (several variants of this file have been used):

title                    = pure water
cpp                      = /usr/bin/cpp
integrator               = md
nsteps                   = 2500000
nstlist                  = 50
nstxout                  = 1000
nstvout                  = 1000
nstlog                   = 1000
dt                       = 0.002
nstcomm                     = 200
comm-mode                 = linear
constraints              = h-bonds
nstenergy                = 100
ns_type                  = grid
coulombtype              = pme
fourierspacing           = 0.12
pme_order                = 4
rlist                    = 0.9
rvdw                     = 0.9
rcoulomb                 = 0.9
tcoupl                   = Berendsen
tc_grps                  = SOL
energygrps               = SOL
tau_t                    = 0.4
ref_t                    = 300



David van der Spoel schrieb:
> aherz wrote:
>> Hi,
>>
>> we're running simulations of a water slab sandwiched by vacuum (ca. 2k
>> waters, 3x3x12 nm^3 box) in NVT
>> and look at the pressure tensor. We are surprised to find that the
>> kinetic energy is equivalent (within statistical uncertainty) for x and
>> y dim(parallel to the interface) but not in z dim. This appears to
>> violate the equipartition theorem?
>>
>> This is gromacs 3.0 and 4.0 using berendsen thermostat.
>
> First, try it with a real thermostat, second please report some numbers.
>
>>
>> Alex
>> _______________________________________________
>> gmx-users mailing list    gmx-users at gromacs.org
>> http://lists.gromacs.org/mailman/listinfo/gmx-users
>> Please search the archive at http://www.gromacs.org/search before
>> posting!
>> Please don't post (un)subscribe requests to the list. Use the www
>> interface or send it to gmx-users-request at gromacs.org.
>> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>
>

-------------- next part --------------
A non-text attachment was scrubbed...
Name: ekin.zip
Type: application/zip
Size: 16877 bytes
Desc: not available
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20090910/9da870f5/attachment.zip>