Gradient evaluation bugs fixed in September 2025

We have discovered errors in the calculation of applied forces of the following collective variables (CVs) in Colvars: eigenvector with rotational fitting enabled, distanceZ with dynamic axis enabled, orientationAngle, orientation, and coordNum or selfCoordNum with pairlist. Restraints along those CVs using Colvars may be ineffective, i.e. produce a trajectory inconsistent with the restraint, or have otherwise unexpected results. Note that most of these issues only affect simulations with non-default and less common settings, which explains why they had not been reported so far.

Below are details:

• eigenvector with rotational fitting: the gradients with respect to the fitting group (fit gradients), which are normally enabled for all CVs, were disabled by default for eigenvector. However, numerical tests show that these fit gradients cannot be ignored (link). This issue does not affect calculations of eigenvector when the keyword rotateToReference is disabled, or of other CVs irrespective of rotational fitting.

• orientation and orientationAngle: because the two quaternions q and -q represent identical orientations, Colvars uses the convention that the first component q0 >= 0, flipping the sign if needed from the calculated quaternion to simplify analysis of trajectories. However, the same correction was not applied to the gradients (link), resulting in the possibility of inconsistently applied forces over a simulation. This issue does not affect the other rotational CVs (orientationProj, tilt, spinAngle, eulerPhi, eulerTheta, eulerPsi). Numerical tests indicate that this only becomes noticeable for angular deviations above 90 degrees, which are unlikely in a restrained simulation.

• coordNum and selfCoordNum with pairlist: when using these two types of coordination number with pairlist enabled, a factor was missing when calculating the gradients of the switching function with respect to the atomic coordinates. Note that the relative errors between the wrong and correct gradients are linearly correlated with the pairlist tolerance. For example, a pairlist tolerance of 0.01 could make the gradients 1% smaller or larger than the correct results. This issue has been fixed (link), and is specific to the pairlist, i.e. does not occur when the parameter “tolerance” is 0 (default).

• distanceZ with dynamic axis: when defining an optional “ref2” group, the reference axis is dynamically defined as the distance vector between “ref” and “ref2”. The gradients with respect to the dynamic axis were incorrectly calculated, and have now been fixed (link). The issue does not affect calculations when the axis is explicitly given using the keyword “axis”, or when it is omitted to use the default axis (0, 0, 1).

The above bugs may affect your simulations only if you are using Colvars to apply biasing forces along the above CVs with the specified options.