Greetings once again fellow Gromacs-users, I have seen a substantial and disturbing change in calculated free energy differences between systems run with Gromacs 3.3 and Gromacs 3.3.1. I am trying to determine if these differences are due to the new software, new parameters, insufficient sampling, or something else, and I don't have enough computer time to investigate every possible combination. I turn to the list for insight.
With Gromacs 3.3, I was using thermodynamic integration over 21 equally spaced and independent (delta_lambda=0) windows with the AMBER99 force field ported by the Pande group. I was running two slightly different DNA systems and using thermodynamic cycles to find the quantity of interest by comparing free energy differences. I was using soft cores and incorrectly using sc_alpha=1.51 because I did not realize that the default sc_power was 1 instead of 2 as in the original soft core method. Nonetheless, my calculations seemed to be in reasonable agreement with experiment and with previous calculations I had performed using NWChem. Recent calculations performed via a slightly different method using AMBER99 within AMBER itself are also in reasonable harmony with these results. When Gromacs 3.3.1 was released, I was chastised by the software for not explictly setting sc_power, and that's when I realized that the default was 1 and that I was using an incorrect value of sc_alpha. I undertook to rerun my calculations with 3.3.1, the same sc_power=1, and the corrected sc_alpha=0.7. I am now getting qualitatively wrong results, though. According to experiment and previous calculations, I should be seeing a difference between my two systems of about -10 kj/mol. In my current simulations I am seeing a difference of about +20 kj/mol. This is after 300 ps of equilibration in each window followed by 500 ps of energy statistics collection. The starting structure for each window was taken from the final frame of a 2 ns ordinary MD production run with free_energy=no and init_lambda=0. The obvious possibility seems to be that the difference in soft cores has nudged the trajectories so that the luck of the draw has me collecting less representative statistics. Is such a large change between free energy cycles credibly attributable to slight sampling differences, though? 30 kj/mol seems pretty substantial. On the other hand, my systems do involve a vast number of degrees of freedom (350-700 DNA atoms, ~20 of which have altered B-states, ~40,000 water atoms, 11 or 22 Na+) and the absolute free energy changes within each system are on the order of 150 kj/mol. Is it possible that the change in sc_alpha alone is responsible for the difference seen with my older runs and newer ones? I have looked at the log files and parameters don't seem to be different between the old systems and the new ones. I know that there was some change in how B-state parameters were assigned with 3.3.1, but for a while I ran a CVS version of Gromacs that included this B-state change and I didn't see any alteration in free energy calculations. Has there been any other change to free energy calculations in 3.3.1? If there are no obvious things to check I will probably be trying fewer windows, with nonlinear spacing, and see if that gives me more believable results that converge faster. Matt Ernst Washington State University _______________________________________________ gmx-users mailing list gmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please don't post (un)subscribe requests to the list. Use the www interface or send it to [EMAIL PROTECTED] Can't post? Read http://www.gromacs.org/mailing_lists/users.php
