Justin Lemkul wrote > On 9/21/12 8:29 AM, Bastien Loubet wrote: >> Dear gmx users, >> >> We recently got a problem with the rerun feature of mdrun, and we request >> your help in order to help to solve it. >> >> We have run a simulation of a large POPC membrane using the coarse >> grained >> Martini force fields. From these simulations we obtained a trr trajectory >> file which contains both the position and the velocity of each martini >> particle every 20 ps. From this trajectory we rerun the simulation using >> the >> -rerun option of the mdrun program but with a different topology for >> which >> we have set all the bonded and non bonded interaction to zero. >> Specifically >> we have set all the force constant to zero in the martini_v2.P.itp and >> the >> martini_v2.0_lipids.itp file, the aim is to calculate the virial due to >> the >> electrostatic forces alone. >> The problem is the following: from the edr file we get from the rerun we >> extract the mean value of the energies using g_energy. For the edr file >> obtained during the simulation (and not the rerun) we obtained for the >> temperature and the kinetic energy, using also a 20ps time interval: >> >> Energy Average Err.Est. RMSD Tot-Drift >> ------------------------------------------------------------------------------- >> Temperature 324.998 0.0049 0.413669 -0.00743801 >> (K) >> Kinetic En. 1.9525e+06 29 2485.21 -44.6815 >> (kJ/mol) >> >> This is to be expected because we have a Nose-Hover temperature coupling >> of >> 325 K. The kinetic energy is also equals to the number of degree of >> freedom >> times half the Boltzmann constant times the temperature. >> However for the rerun, with no bonded and non-bonded interactions, we >> get: >> >> Energy Average Err.Est. RMSD Tot-Drift >> ------------------------------------------------------------------------------- >> Temperature 327.35 0.0053 0.415978 -0.00514402 >> (K) >> Kinetic En. 1.96663e+06 32 2499.08 -30.8973 >> (kJ/mol) >> >> The kinetic energy is again equals to the number of degree of freedom >> times >> half the Boltzmann constant times the temperature, however the >> temperature >> is 2 K off ! This is surprising as only the velocities and the masses of >> the >> particle should enter the kinetic energy. The velocities are taken from >> the >> trr trajectory file and we have checked using the Linux diff utility that >> our topology have the same masses than the original one. >> Another cause can be that we are doing the rerun on a local machine while >> the original simulation was run on a cluster, however we believe that >> this >> kind of numerical error cannot be responsible for a 2 degree Kelvin >> mistake. >> Furthermore we have rerun the simulation using the original topology and >> we >> obtained: >> >> Energy Average Err.Est. RMSD Tot-Drift >> ------------------------------------------------------------------------------- >> Temperature 324.998 0.005 0.411827 -0.00761929 >> (K) >> Kinetic En. 1.9525e+06 30 2474.15 -45.7788 >> (kJ/mol) >> >> The slightly different value are certainly due to the previously >> mentioned >> fact that the reruns were not run on the same machine as the simulation >> and >> numerical rounding errors. This point out that it is really the topology >> that changes the value of the temperature. However we really have no >> masses >> differences between the two topologies and the velocities are draw from >> the >> same trr trajectory file. >> >> Any thought on this problem would be welcome. >> > > Another thing to consider is the fact that an .edr file is accurate over > all > simulation steps, while reruns only do calculations present in the frames > supplied. If you are using frames every 20 ps, that may only represent > less > than 1% of the total data collected in the simulation, depending on the > value of dt. > > -Justin > > -- > ======================================== > > Justin A. Lemkul, Ph.D. > Research Scientist > Department of Biochemistry > Virginia Tech > Blacksburg, VA > jalemkul[at]vt.edu | (540) 231-9080 > http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin > > ======================================== > -- > gmx-users mailing list
> gmx-users@ > http://lists.gromacs.org/mailman/listinfo/gmx-users > * Please search the archive at > http://www.gromacs.org/Support/Mailing_Lists/Search before posting! > * Please don't post (un)subscribe requests to the list. Use the > www interface or send it to > gmx-users-request@ > . > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists Alright, but we did make a rerun with the same trajectory and the original topology and it gives the third results in my original post: the 2 Kelvin difference is not there. So it is really the different topology file that give the temperature difference I believe. -- View this message in context: http://gromacs.5086.n6.nabble.com/Possible-bug-in-the-temperature-calculation-from-rerun-tp5001194p5001196.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- 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/Support/Mailing_Lists/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/Support/Mailing_Lists