Loop over replicas when computing well-tempered metadynamics scaling factor

[giacomofiorin]

Fixes #819

[giacomofiorin]

Still working on this (including test, documentation update).

When finished, I want to use it as a test for a CI job that cherry-picks bugfixes into release branches (e.g. namd-3.0, gromacs-2025).

[giacomofiorin]

@HanatoK I didn't understand this comment before, and forgot to ask you about it

[HanatoK]

As I understand, without keepHills, the code only accumulates the hills that are off-grid. However, in theory, those Gaussian hills that have been projected onto the histogram are truncated, which means they do not affect the bias energy of an off-grid position.

With keepHills I think it would be better to compute the bias energy from all hills, no matter whether they are in the histogram.

[giacomofiorin]

If you're doing that, then you would effectively disable use_grids without the user being able to check.

To be honest, at this point trajectories have become so long that the usefulness of keepHills is put into question.

[HanatoK]

No. What I suggested should only affect off-grid positions. If you use both use_grids and keepHills, then to calculate the energy of a position out of the grid, you can sum all hills, and to calculate the energy of a position in the grid, you can just use the grid value. This shouldn't disable use_grids.

[giacomofiorin]

Ok, thanks. The "off-grid hills", which are kept regardless, should also include hills whose center is inside the grid boundaries, but they are within sufficient range that they would affect configurations just outside the grid. At least that was the original logic, assuming that it didn't break along the way.

Looping over all hills would also include hills that are well within the boundaries.

[HanatoK]

In my opinion, the logic of computing bias energy should be the same as colvarbias_meta::calc_energy without adding the new hills:

  bool index_ok = false;
  std::vector<int> curr_bin;

  if (use_grids) {

    curr_bin = values ?
      hills_energy->get_colvars_index(*values) :
      hills_energy->get_colvars_index();

    index_ok = hills_energy->index_ok(curr_bin);

  }

  if ( index_ok ) {
    // index is within the grid: get the energy from there
    for (ir = 0; ir < replicas.size(); ir++) {

      bias_energy += replicas[ir]->hills_energy->value(curr_bin);
      if (cvm::debug()) {
        cvm::log("Metadynamics bias \""+this->name+"\""+
                 ((comm != single_replica) ? ", replica \""+replica_id+"\"" : "")+
                 ": current coordinates on the grid: "+
                 cvm::to_str(curr_bin)+".\n");
        cvm::log("Grid energy = "+cvm::to_str(bias_energy)+".\n");
      }
    }
  } else {
    // off the grid: compute analytically only the hills at the grid's edges
    for (ir = 0; ir < replicas.size(); ir++) {
      calc_hills(replicas[ir]->hills_off_grid.begin(),
                 replicas[ir]->hills_off_grid.end(),
                 bias_energy,
                 values);
    }
  }

[HanatoK]

Is it possible to move the checking out of the for loop, just like colvarbias_meta::calc_energy?

[giacomofiorin]

I moved the initialization of index_ok outside of the loop. As far as moving the check, I'm not sure whether it's worth it (meaning, I don't have enough bandwidth to do this, but wouldn't oppose you if you did).

[HanatoK]

Now I feel the old code summing the hills starting from new_hills_begin is incorrect, and this fixes it. Is that right?

[giacomofiorin]

I'm not sure, sorry. This branch corresponds to the grids-off, multiple-walkers, well-tempered use case.

With this change, the branch is more consistent with the other options. As for whether it was previously buggy, I think we should first check if anybody ever ran this case before investigating further.