Re: [lammps-users] Inquiry on LAMMPS Colvar US

Hi Victor, I’m replying to the list, since your questions are general and answers are most useful to everyone.

Dear Dr. Giacomo,

Hi Dr. Giacomo. Victor here and sorry to disturb you again. I have some questions on the PMF output value from the file out.h_pot.ti.pmf from LAMMPS Colvars module Umbrella Sampling method:

  1. So my case is to drag the polymer molecule (polyacrylamide+acrylic acid) away from the solid surface to evaluate the free energy adsorption. At the same time, I want to try another way by placing the molecule at a far distance and bring the molecule close to the surface to see any difference in the PMF curve value. Theoretically, in either case the relative PMF value should be similar? However, as shown in the output file in the attached link:

-the value is quite different for the case of pulling from distance 8 to distance 25 (I forget to change the upper limit boundary to 25) with the case distance 25 to 7.

-For the case of pulling distance 25 to 7 & distance 25 to 3, even though the starting point is the same, the PMF value also differs a lot. … May I know why, and how I should interpret the result?

I removed the link, but what you are describing looks like what is generally described as histeresis: the work from A to B is not the same as the one performed from B to A, because neither is a reversible path. I took a quick look to see that you ran 2 million steps, which is of the order of nanoseconds (the exact time span depends on your force field and setup). Once you factor in that you’re spanning multiple states, that’s mere picoseconds for each window that you try equilibrating overt. This is probably too short of a time for the polymer’s internal degrees of freedom to equilibrate.

You need to go for much longer simulation times, but even then the complexity of the polymer chain’s internal structure (whose changes need to be averaged out at each window) may be prohibitive. Try computing the adsorption/desorption of a small molecule (monomer?) to see when you get reasonable results. Comparing results between WHAM and TI is very useful, but even more so doing the reverse direction as you have tried to.

  1. Last time Dr. you mentioned, I can use WHAM method to post-process the output file of US. I managed to find a slide that also uses the LAMMPS colvars module. The slide, at page 22 shows how to use a WHAM utility tool to analyse the free energy. However, I’m not sure how he got so many separate COLVAR output files with separate window stages, since the module only produces 1 output file- out.colvars.traj. Is that I need to separate the trajectory frame for each window from the out.colvars.traj?
    https://cpb-us-w2.wpmucdn.com/sites.udel.edu/dist/2/334/files/2018/03/hpc_sym_2018_redacted-reducedsize-z6ijvx.pdf
    http://membrane.urmc.rochester.edu/content/wham

Yep, I saw your thread on the Plumed mailing list as well :slight_smile: That’s a useful set of slides.

The gist of it is that each window is a separate trajectory segment, and it depends on one’s preference how these are run. If you run the windows sequentially, for Grossfield’s WHAM implementation you need to extract all the segments where the restraint is piecewise constant, e.g. using the Python script that I pointed you to.

Alternatively, others also run windows concurrently as independent simulations, i.e. they extract frames from a quick simulation like your and extend them as far as they can.

There is no universal rule about this, because it all hinges upon the variable(s) of interest decorrelating quickly enough that each window yields a converged histogram that is also representative of the true system. Again, this is where it helps if you’re not an expert to build a much simplified version of your system and gradually add back complexity while checking that you can still get reliable PMFs along the process.

Giacomo

Hi Dr. Thanks for your reply. Ya I will try again for the case of polymer using the ‘monomer’ segment. I have also tested using a small water molecule, and the PMF value from point A to point B, and point B to point A is not that much different, so the hysteresis problem is not that obvious for a small molecule. I will test again with longer simulation time. However, for 2 different cases starting from same point A (example: initial distance of 15) to a different point B (target distance 5, or target distance 3) pulled real close into the solid surface), the highest PMF value at large separation can be very different with the slight change in target distance value. Is that common?

Yes. I can’t emphasize this enough: if the variables selected do not capture all relevant degrees of freedom, there is no point in even attempting a PMF computation.

I now have some clues how to use the WHAM tool properly,… Just curious how many frames are ideal for each window? For now, I run with total 1ns for 20 window stages, each stage is around 50ps, and with TrajFrequency of 10ps, so each window only has 5 frames. If I know the ideal number frame for each window then I can adjust my total simulation time accordingly…

You need many more frames than 5, obviously, but most of all you need many uncorrelated frames. A histogram may look very smooth and the corresponding PMF have very small “error bars” but still only represent a small subset of the relevant states. Again, one of the most straightforward approaches is doing forward and reverse transitions. Methods other than umbrella sampling (e.g. ABF, metadynamics) incorporate these multiple sweeps into the algorithm itself.

Giacomo

Ok… noted. Thanks again!

Hi Dr. Thanks for your reply. Ya I will try again for the case of polymer using the ‘monomer’ segment. I have also tested using a small water molecule, and the PMF value from point A to point B, and point B to point A is not that much different, so the hysteresis problem is not that obvious for a small molecule. I will test again with longer simulation time. However, for 2 different cases starting from same point A (example: initial distance of 15) to a different point B (target distance 5, or target distance 3) pulled real close into the solid surface), the highest PMF value at large separation can be very different with the slight change in target distance value. Is that common?

I now have some clues how to use the WHAM tool properly,… Just curious how many frames are ideal for each window? For now, I run with total 1ns for 20 window stages, each stage is around 50ps, and with TrajFrequency of 10ps, so each window only has 5 frames. If I know the ideal number frame for each window then I can adjust my total simulation time accordingly…

Thank you!

Regards,
Victor