Check force on energy
It is common practise to use the potential energy as a collective energy. Some MD codes thus pass the potential energy to PLUMED and PLUMED can then apply forces on this collective variable. We test that any forces that PLUMED applies on the potential energy are correctly passed back to the MD code by doing the following test. We first run a short simulation at $T$ K with a timestep of $\tau$ ps. During the course of this simulation we monitor the potential energy using the following PLUMED input:
e:ENERGYv :Calculate the total potential energy of the simulation box. More detailsVOLUMECalculate the volume of the simulation box. More detailsPrint quantities to a file. More detailsARG=e,vthe input for this action is the scalar output from one or more other actionsFILE=energy1the name of the file on which to output these quantities
We then run a second simulation (starting from identical conditions) at a temperature of $T\alpha$ and with a timestep of $\tau/\sqrt(\alpha)$. The thermostat and barostat relaxation times are similarly divided by $\sqrt(\alpha)$. In the tests that are run on this website we set $\sqrt(\alpha)=1.1$. The PLUMED file above is used when this test is run but a different time series of energy values is recorded as the MD parameters in this second simulation are different.
If PLUMED is working correctly we should be able to recapture the time series of energy values for the first simulation by running an MD simulation with the modified parameters that were used in the second simulation and the following PLUMED input file:
e:ENERGYv :Calculate the total potential energy of the simulation box. More detailsVOLUME# slope is such thatCalculate the volume of the simulation box. More detailsPrint quantities to a file. More detailsARG=ethe input for this action is the scalar output from one or more other actionsFILE=energy2 # slope should be (alpha-1)=0.21the name of the file on which to output these quantitiesRESTRAINTAdds harmonic and/or linear restraints on one or more variables. More detailsAT=0.0the position of the restraintARG=ethe arguments on which the bias is actingSLOPE=0.21specifies that the restraint is linear and what the values of the force constants on each of the variables are
In other words, when forces are passed correctly the time series for the energies and volumes from the first and third of these calculations should be identical.
To determine if PLUMED passes this test we calculate the difference between the time series that were observed in the first and third simulations described above. We then divide this by the difference between the first and second time series.
An NPT version of this calculation is performed as well as an NVT calculation if the virial is passed to PLUMED.
Trajectories
- Input and output files for the unpeturbed calculation are available in this zip archive
-
Input and output files for the peturbed calculation are available in this zip archive
- Input and output files for the peturbed calculation in which a PLUMED restraint is used to undo the effect of the changed MD parameters are available in this zip archive
Results
The table below includes some of the results from the calculation. The columns contain:
- Time series for the energy and volume that were obtained from the simulation at $T$ K, $x_{md}$.
- Time series for the energy and volume that were obtained from the simulation at $\alpha T$ K and in which PLUMED applied a restraint on the energy, $x_{pl}$.
- The absolute value of the difference between the time series of energies and volumes that were obtained from the simulations running at $T$ K and $\alpha T$ K, $\vert x_{md}’-x_{md} \vert$. No PLUMED restraints were applied in either of these simulations.
- The values of $100\frac{\vert x_{md} - x_{pl}\vert }{ \vert x_{md}’-x_{md} \vert}$.
If the PLUMED interface is working correctly the first two sets of numbers should be identical and the final column should be filled with zeros.
Original result | Result with PLUMED | Effect of peturbation | % Difference |
---|---|---|---|
-683.6798 127.9326 | -683.6798 127.9326 | 0.0000 0.0000 | 0.0000 0.0000 |
-682.8973 127.9326 | -682.8973 127.9326 | 0.0003 0.0000 | 0.0000 0.0000 |
-681.5577 127.9326 | -681.5577 127.9326 | 0.0021 0.0000 | 0.0000 0.0000 |
-679.6256 127.9326 | -679.6256 127.9326 | 0.0077 0.0000 | 0.0000 0.0000 |
-677.0148 127.9326 | -677.0148 127.9326 | 0.0212 0.0000 | 0.0000 0.0000 |
-673.6763 127.9326 | -673.6763 127.9326 | 0.0499 0.0000 | 0.0000 0.0000 |
-669.5530 127.9326 | -669.5530 127.9326 | 0.1057 0.0000 | 0.0000 0.0000 |
-664.4800 127.9326 | -664.4800 127.9326 | 0.2081 0.0000 | 0.0000 0.0000 |
-658.4227 127.9326 | -658.4227 127.9326 | 0.3861 0.0000 | 0.0000 0.0000 |
-651.3662 127.9326 | -651.3662 127.9326 | 0.6822 0.0000 | 0.0000 0.0000 |
-643.3903 127.9326 | -643.3903 127.9326 | 1.1500 0.0000 | 0.0000 0.0000 |
-634.6407 127.9326 | -634.6407 127.9326 | 1.8578 0.0000 | 0.0000 0.0000 |
-625.3528 127.9326 | -625.3528 127.9326 | 2.8766 0.0000 | 0.0000 0.0000 |
-616.1227 127.9326 | -616.1227 127.9326 | 4.2402 0.0000 | 0.0000 0.0000 |
-607.3890 127.9326 | -607.3890 127.9326 | 5.9555 0.0000 | 0.0000 0.0000 |
-599.6953 127.9326 | -599.6953 127.9326 | 7.9395 0.0000 | 0.0000 0.0000 |
-593.5389 127.9326 | -593.5389 127.9326 | 10.0368 0.0000 | 0.0000 0.0000 |
-588.9711 127.9326 | -588.9711 127.9326 | 12.0824 0.0000 | 0.0000 0.0000 |
-586.2146 127.9326 | -586.2146 127.9326 | 13.8552 0.0000 | 0.0000 0.0000 |
-585.0657 127.9326 | -585.0657 127.9326 | 15.2282 0.0000 | 0.0000 0.0000 |