<html>
<head>
<style>
.hmmessage P
{
margin:0px;
padding:0px
}
body.hmmessage
{
font-size: 10pt;
font-family:Verdana
}
</style>
</head>
<body class='hmmessage'>
Equipartitioning still holds.<br>The issue here is that constraints remove degrees of freedom.<br>If you remove a degree of freedom, there is not potential or kinetic energy<br>anymore in that generalized coordinate.<br>This does not introduce any problems with ensembles or something like that.<br><br>I did not have time yet to implement pull code with pbc.<br>I will do it this week.<br><br>Berk<br><br>> Date: Mon, 14 Sep 2009 16:08:55 +0200<br>> From: alexander.herz@mytum.de<br>> To: gmx-users@gromacs.org<br>> Subject: Re: [gmx-users] non isotropic kinetic energy<br>> <br>> Hm,<br>> <br>> apparently there are constraint algorithms which don't change the ensemble<br>> (e.g. http://www.informaworld.com/smpp/content~db=all~content=a750934359).<br>> The equipartition theorem requires that the energy depends on<br>> no other powers of v but v^2<br>> (http://www.whfreeman.com/modphysics/PDF/8-2c.pdf)<br>> which should be the case for the md systems including constraints<br>> (lj pot and estatic only depend on the particles' positions). So as far<br>> as I can see the<br>> equipartition theorem should hold in td equilibrum for these systems. So<br>> the constraints<br>> might modify the ensemble which could lead to all kinds of trouble?<br>> <br>> The water molecules at the interface are oriented but taking the above<br>> paragraph into<br>> account this shouldn't cause any trouble.<br>> <br>> By the way,<br>> did you find a chance to look at the pulling problem we discussed before<br>> last week?<br>> <br>> Thx,<br>> Alex<br>> <br>> <br>> Berk Hess schrieb:<br>> > I have been looking at very similar things recently.<br>> ><br>> > Your result could actually be correct.<br>> > If the water orders at the interface you should have a difference<br>> > in Ekin components, since the constraints, which remove kinetic energy<br>> > are not ordered isotropically in space.<br>> > You should be able to quantify this by looking at the order of the water<br>> > normal.<br>> ><br>> > Berk<br>> ><br>> > > Date: Thu, 10 Sep 2009 18:11:16 +0200<br>> > > From: alexander.herz@mytum.de<br>> > > To: gmx-users@gromacs.org; felix.sedlmeier@ph.tum.de<br>> > > Subject: Re: [gmx-users] non isotropic kinetic energy<br>> > > CC:<br>> > ><br>> > > Hi,<br>> > ><br>> > > we took a couple of days to reexamine this and try out a couple of<br>> > things:<br>> > ><br>> > > -On small interfacial systems (ca. 2k spc/e water molecules sandwiched<br>> > > by vacuum)<br>> > > the difference between ekin_x/y and ekin_z is 20 to 60kJ/mol depending<br>> > > on the specific system/settings<br>> > > (e.g. see attached file of averaged ekin (x/y/z) over time).<br>> > > -On larger interfacial systems the difference grows (up to several<br>> > > 100kJ/mol for ca. 32k spc/e water).<br>> > ><br>> > > The following settings make no difference (so the error persists even if<br>> > > we change the setting):<br>> > > -thermostat (we tryed nose-hoover, berendsen and no thermostat)<br>> > > -estatics (cut-off instead of pme)<br>> > > -type of constraint algo for water (so lincs, shake or settle, we<br>> > > modified the spc/e water molec definition to change this)<br>> > > -changing system size so that the box length is equal in all 3 dim (as<br>> > > said, magnifies problem)<br>> > > -rotating system (error rotates with system)<br>> > > -use smaller time step of 1fs (no difference)<br>> > ><br>> > > We ran a simulation with flexible water (so bonds instead of<br>> > > constraints) , here the error vanishes!<br>> > > Also, we ran similar sims with amber, where the error does not appear at<br>> > > all (for constraint water).<br>> > ><br>> > > We found out about this problem while doing surface tension/contact<br>> > > angle calculations.<br>> > > The surf tens calced via virial differs from the surf tens calculated<br>> > > via pressure tensor (unless one corrects<br>> > > for the anisotropic kinetic energy). A simple way to see the deficit is<br>> > > trying to calc the average pressure tensor from the average virial<br>> > > tensor (using the simple def from the gromacs manual<br>> > P_ij=2/V(Ekin-Vir_ij):<br>> > > a) assuming equidistributed kinetic energy (ekin/3 per dim)<br>> > > b) using the anisotropic ekin which can be extracted from a trr with a<br>> > > tool (source attached) if velocities have been written to the trr.<br>> > ><br>> > > All this appears to hint at a problem with the constraints for the<br>> > > waters and might also hint at an incorrect ensemble as a consequence.<br>> > ><br>> > > Alex<br>> > ><br>> > > sample input file (several variants of this file have been used):<br>> > ><br>> > > title = pure water<br>> > > cpp = /usr/bin/cpp<br>> > > integrator = md<br>> > > nsteps = 2500000<br>> > > nstlist = 50<br>> > > nstxout = 1000<br>> > > nstvout = 1000<br>> > > nstlog = 1000<br>> > > dt = 0.002<br>> > > nstcomm = 200<br>> > > comm-mode = linear<br>> > > constraints = h-bonds<br>> > > nstenergy = 100<br>> > > ns_type = grid<br>> > > coulombtype = pme<br>> > > fourierspacing = 0.12<br>> > > pme_order = 4<br>> > > rlist = 0.9<br>> > > rvdw = 0.9<br>> > > rcoulomb = 0.9<br>> > > tcoupl = Berendsen<br>> > > tc_grps = SOL<br>> > > energygrps = SOL<br>> > > tau_t = 0.4<br>> > > ref_t = 300<br>> > ><br>> > ><br>> > ><br>> > > David van der Spoel schrieb:<br>> > > > aherz wrote:<br>> > > >> Hi,<br>> > > >><br>> > > >> we're running simulations of a water slab sandwiched by vacuum<br>> > (ca. 2k<br>> > > >> waters, 3x3x12 nm^3 box) in NVT<br>> > > >> and look at the pressure tensor. We are surprised to find that the<br>> > > >> kinetic energy is equivalent (within statistical uncertainty) for<br>> > x and<br>> > > >> y dim(parallel to the interface) but not in z dim. This appears to<br>> > > >> violate the equipartition theorem?<br>> > > >><br>> > > >> This is gromacs 3.0 and 4.0 using berendsen thermostat.<br>> > > ><br>> > > > First, try it with a real thermostat, second please report some<br>> > numbers.<br>> > > ><br>> > > >><br>> > > >> Alex<br>> > > >> _______________________________________________<br>> > > >> gmx-users mailing list gmx-users@gromacs.org<br>> > > >> http://lists.gromacs.org/mailman/listinfo/gmx-users<br>> > > >> Please search the archive at http://www.gromacs.org/search before<br>> > > >> posting!<br>> > > >> Please don't post (un)subscribe requests to the list. Use the www<br>> > > >> interface or send it to gmx-users-request@gromacs.org.<br>> > > >> Can't post? Read http://www.gromacs.org/mailing_lists/users.php<br>> > > ><br>> > > ><br>> > ><br>> ><br>> > ------------------------------------------------------------------------<br>> > See all the ways you can stay connected to friends and family<br>> > <http://www.microsoft.com/windows/windowslive/default.aspx><br>> > ------------------------------------------------------------------------<br>> ><br>> > _______________________________________________<br>> > gmx-users mailing list gmx-users@gromacs.org<br>> > http://lists.gromacs.org/mailman/listinfo/gmx-users<br>> > Please search the archive at http://www.gromacs.org/search before posting!<br>> > Please don't post (un)subscribe requests to the list. Use the <br>> > www interface or send it to gmx-users-request@gromacs.org.<br>> > Can't post? Read http://www.gromacs.org/mailing_lists/users.php<br>> <br>> _______________________________________________<br>> gmx-users mailing list gmx-users@gromacs.org<br>> http://lists.gromacs.org/mailman/listinfo/gmx-users<br>> Please search the archive at http://www.gromacs.org/search before posting!<br>> Please don't post (un)subscribe requests to the list. Use the <br>> www interface or send it to gmx-users-request@gromacs.org.<br>> Can't post? Read http://www.gromacs.org/mailing_lists/users.php<br><br /><hr />Express yourself instantly with MSN Messenger! <a href='http://clk.atdmt.com/AVE/go/onm00200471ave/direct/01/' target='_new'>MSN Messenger</a></body>
</html>