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:

Click on the labels of the actions for more information on what each action computes
tested on2.10
tested onmaster
e: ENERGYCalculate the total potential energy of the simulation box. More details
v: VOLUMECalculate the volume the simulation box. More details
PRINTPrint quantities to a file. More details ARGthe labels of the values that you would like to print to the file=e,v FILEthe name of the file on which to output these quantities=energy1

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:

Click on the labels of the actions for more information on what each action computes
tested on2.10
tested onmaster
e: ENERGYCalculate the total potential energy of the simulation box. More details
v: VOLUMECalculate the volume the simulation box. More details
# slope is such that 
PRINTPrint quantities to a file. More details ARGthe labels of the values that you would like to print to the file=e FILEthe name of the file on which to output these quantities=energy2
# slope should be (alpha-1)=0.21
RESTRAINTAdds harmonic and/or linear restraints on one or more variables. More details ATthe position of the restraint=0.0 ARGthe values the harmonic restraint acts upon=e SLOPE specifies that the restraint is linear and what the values of the force constants on each of the variables are=0.21

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

  1. Input and output files for the unpeturbed calculation are available in this zip archive

  2. Input and output files for the peturbed calculation are available in this zip archive

  3. 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

Original With PLUMED Effect of peturbation % Difference
-18174.8223 11.6346 -18174.8223 11.6346 0.0000 0.0000 0.0000 0.0000
-18171.8203 11.6346 -18193.2891 11.6346 6.8730 0.0000 312.3615 0.0000
-18193.5547 11.6475 -18233.7578 11.6321 21.9727 0.0054 182.9689 284.9140
-18164.6562 11.6475 -18286.6484 11.6321 17.0625 0.0054 714.9725 284.9140
-18149.8594 11.6475 -18326.3887 11.6321 18.6348 0.0054 947.3116 284.9140
-18130.2344 11.6475 -18349.4688 11.6321 18.6484 0.0054 1175.6180 284.9140
-18111.7656 11.6475 -18356.0703 11.6321 19.2793 0.0054 1267.1867 284.9140
-18100.5645 11.6475 -18352.4355 11.6321 22.2520 0.0054 1131.9056 284.9140
-18099.7383 11.6475 -18346.5566 11.6321 27.4492 0.0054 899.1817 284.9140
-18107.9746 11.6475 -18344.0215 11.6321 32.8652 0.0054 718.2267 284.9140
-18120.1504 11.6475 -18346.2070 11.6321 35.6582 0.0054 633.9541 284.9140
-18130.6543 11.6475 -18351.4434 11.6321 34.4355 0.0054 641.1661 284.9140
-18136.8770 11.6475 -18358.1055 11.6321 30.9453 0.0054 714.9016 284.9140
-18141.5742 11.6475 -18367.0137 11.6321 29.8770 0.0054 754.5597 284.9140
-18150.7441 11.6475 -18381.0078 11.6321 35.5469 0.0054 647.7747 284.9140
-18168.2754 11.6475 -18401.1016 11.6321 47.2656 0.0054 492.5909 284.9140
-18191.3242 11.6475 -18423.2988 11.6321 58.0176 0.0054 399.8350 284.9140
-18211.0723 11.6475 -18439.4434 11.6321 58.6016 0.0054 389.7014 284.9140
-18218.4766 11.6475 -18442.3301 11.6321 44.7539 0.0054 500.1877 284.9140
-18210.0391 11.6475 -18430.7285 11.6321 20.5312 0.0054 1074.8954 284.9140

The table below includes some of the results from the calculation. The columns contain:

  1. Time series for the energy and volume that were obtained from the simulation at $T$ K, $x_{md}$.
  2. 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}$.
  3. 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.
  4. 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.

Graphical representation (beta)

A visualization of the table above:
engvir_v2.10