Hm, apparently there are constraint algorithms which don't change the ensemble (e.g. http://www.informaworld.com/smpp/content~db=all~content=a750934359). The equipartition theorem requires that the energy depends on no other powers of v but v^2 (http://www.whfreeman.com/modphysics/PDF/8-2c.pdf) which should be the case for the md systems including constraints (lj pot and estatic only depend on the particles' positions). So as far as I can see the equipartition theorem should hold in td equilibrum for these systems. So the constraints might modify the ensemble which could lead to all kinds of trouble?
The water molecules at the interface are oriented but taking the above paragraph into account this shouldn't cause any trouble. By the way, did you find a chance to look at the pulling problem we discussed before last week? Thx, Alex Berk Hess schrieb: > 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.h...@mytum.de > > To: gmx-users@gromacs.org; felix.sedlme...@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@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-requ...@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> > ------------------------------------------------------------------------ > > _______________________________________________ > gmx-users mailing list gmx-users@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-requ...@gromacs.org. > Can't post? Read http://www.gromacs.org/mailing_lists/users.php _______________________________________________ gmx-users mailing list gmx-users@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-requ...@gromacs.org. Can't post? Read http://www.gromacs.org/mailing_lists/users.php
