[gmx-users] non isotropic kinetic energy

[Berk Hess]

I have been looking at very similar things recently.

Your result could actually be correct.
If the water orders at the interface you should have a difference
in Ekin components, since the constraints, which remove kinetic energy
are not ordered isotropically in space.
You should be able to quantify this by looking at the order of the water
normal.

Berk

> Date: Thu, 10 Sep 2009 18:11:16 +0200
> From: alexander.herz at mytum.de
> To: gmx-users at gromacs.org; felix.sedlmeier at ph.tum.de
> Subject: Re: [gmx-users] non isotropic kinetic energy
> CC: 
> 
> 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
> >
> >
> 

_________________________________________________________________
See all the ways you can stay connected to friends and family
http://www.microsoft.com/windows/windowslive/default.aspx
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20090911/642d8921/attachment.html>