I am trying to replicate a LJ gas-liquid nucleation study by Wolde and Frenkel (J.Chem. Phy. 109, 9901 (1998)), where umbrella sampling was used to bias the sampling of cluster of specific sizes: a bias potential of harmonic form W = k(n-n0)^2 (n = size of the largest cluster, n0 = target size) was applied in order to sample unstable clusters that are difficult to observe otherwise.
I want to implement this system on LAMMPS. I saw that umbrella sampling is included in the ‘colvar’ package, but the included functions in the ‘colvar’ library do not quite apply to my application.
I would like to get suggestions on incorporating the bias potential into the simulation:
Hi Junli, you’re right that that specific function is not implemented in Colvars directly, but there are some alternatives.
One possibility is using a custom function:
which uses the Lepton library to compute function of collective variables already implemented, or Cartesian coordinates if you need to start from scratch.
Another possibility is use the coming LAMMPS stable release, which will now carry a better maintained interface with Plumed (http://www.plumed.org/), which includes many functions not implemented by Colvars. You will still have to compile LAMMPS and Plumed separately, but the linking procedure is much simpler than it used to be.
Regarding using existing computes such as compute/chunk as variables, neither library will support it (I’m sure of it for Colvars, and I guess this is true for Plumed as well). Computes are very useful to compute global and local quantities, but because they are not involved in force computation they do not provide atomic gradients for their functions. It would be more appropriate to use a fix object, but those that I am familiar with do not compute a global quantity that can be called a “variable”.
If you know of a fix type that is capable of doing that, it would be a legitimate feature request to be able to use it as a collective variable with an enhanced sampling method.
I should also add that all the above applies in the case of molecular dynamics, where forces are needed.
For Monte Carlo, you can of course use a compute as a variable, for as long as you have a way to generate the random moves. Computing the umbrella potential of an existing compute is also easily done using existing LAMMPS commands.