Next: Using HFB in your research Up: Up Contents Index

Welcome to the HFB Documentation under CHARMM

HFB subcommands (within the TReK module from Stephan Fischer)

Assuming K is the number of beads in the path

FP -: {DENS | INIT [STRT (1) STOP (2)] } -
[GRID (K*4)] [TRNC (MAXP*4)]
WRKP -: [TRNC (K/4*3)] [RTYP (1)] -
[BEAD (K)] [GRID (K*4)] -
[RFRC (0.0)] [STEP (0.0)] -
[IPRF (-1)] [PTRJ] [PRJC] [FEND]

The following keywords apply to ggaHFB method

RCS HRCM -: [IECM (-1)] [IRCM (-1)] [ORCM (-1)] -
[SELE ... END] -
[MASS] [BFIT] -
[PRJC] [FEND] -
[TRNC (-1)] -
[STEP (0.0)] -
[FORC (0.0)]
|[CLEA]
RCS HRAP -: [SELE ... END] -
[MASS] [BFIT] -
[PRJC] [FEND] -
[TRNC (-1)] -
[STEP (0.0)] -
[FORC (0.0)]
|[CLEA]
HHFB -: [TRNC (K/4*3)] -
[BEAD (K)] [GRID (K*4)] -
[RFRC (0.0)] [STEP (0.0)] -
[IPRF (-1)] [PTRJ] [PRJC] [FEND]

Invoking the HFB method under the TReK

The following is how the TREK module is typically invoked before the HFB method can be used.
DEFINE RCS SELE ... END !This defines the RCS
TREK MAXP K SELE RCS END MASS
This sets the MAXP (the maximum number of points (beads) in the path) to K - the desired number of the total beads in the path. The MAXP is important in the context of the HFB method as this number affects path manipulations.
At present the RCS selection is passed into TREK via a single atom selection at the TREK invocation level. The same selection is used as best-fit selection. Only a single selection is allowed that defines reaction coordinate subspace (RCS). In the future, I am considering adding additional flexibility and providing another selection for the best-fit procedure that essentially defines the coordinate system and does not need to coincide with the RCS. This option might be vital for certain applications and definitions of more sound coordinate systems irrespective of the RCS.
The MASS keyword requests mass weighting during the best-fit procedure. The MASS keyword is not related to the biasing potential used. The type of the biasing potential used is determined independently later.

Detailed HFB keywords

FP: Stands for Fourier Path. This keyword calls the original gradient-free implementation of the HFB method (paper #1)
If FP is invoked MAXP must equal the final (desired) number of beads in the path.
Although superceded by gradient-driven method, the FP is still useful in some cases as it has greater flexibility in terms of beads manipulations. Recall that FP simply parameterizes and then redistributes the points that are read in, and therefore depends only on these points. In contrast the enhanced HFB method in addition to the evolved points requires the corresponding reference points to generate a new set of reference points.

Keywords under the original FP method

FP -: {DENS | INIT [STRT (1) STOP (2)] } -
[GRID (K*4)] [TRNC (MAXP*4)]

Anything read into the FP is translated into the center of mass of the specified RCS and then superposed onto the reactant (the first bead) before any other operations. So the reactant orientation is used as the template.

DENS (logical; default FALSE): triggers expansion of the trajectory up to MAXP (assuming a fewer (or equal) number of points is read in initially)
(give an example): assuming you are reading in a trajectory with N structures and want to get K (=MAXP)> N
INIT (logical; default FALSE): triggers interpolation between two immediately adjacent points in the path such that the total number of points becomes MAXP. This is the key option for the activated evolution procedure.
STRT (integer; default 1): designates the first point for the line interpolation
STOP (integer; default 2): designates the second point for the line interpolation
Note: STOP - STRT = 1 (this is a must, could be trivially changed in the future though)
TRNC (integer; default MAXP): the maximum number of points. MAXP allows user to specify the number of basis functions (sines) used in the Fourier series expansion. Usually it is recommended to start with a small number of approximately MAXP/2 and gradually increase the number from there as desired. In optimization it is not recommended to use TRNC = MAXP as this will lead to path beats and possible instabilities

The FP version also does not really care what kind of restraint is used as long as the structures returned are consistent with it.

WRKP: Stands for WoRK Path. This keyword corresponds to the gradient-augmented HFB (gaHFB) method (paper #2) that in addition to the improved path optimization also provides the (free-)energy profile computed via the reversible work integral of either potential or free energy forces in the reactive coordinate space RCS.
Recall that unlike FP the gaHFB method in addition to the raw evolved beads requires the corresponding reference points to generate the new set of reference beads.

Keywords under the enhanced WRKP method

WRKP -: [TRNC (K/4*3)] [RTYP (1)] -
[BEAD (K)] [GRID (K*4)] -
[RFRC (0.0)] [STEP (0.0)] -
[IPRF (-1)] [PTRJ] [PRJC] [FEND]

If WRKP is invoked MAXP must be at least double the number of the actual beads used in the path.

PRJC (logical; default FALSE): invokes gradient projection orthogonal to the path. If not used, the gradient will be used as is, including the part that slides along the path. The latter is advantageous when endpoints are being optimized as well. Note that the endpoints could be constrained (see below).
FEND (logical; default FALSE): constraints the endpoints in the RCS so that they do not move during the optimization. Hence, one does not need to evolve the corresponding end point beads. However, at present the HFB routines still expect the corresponding structures, so one has to trick the HFB by supplying the same end-structures at each optimization step.
PTRJ (logical; default FALSE): requests that the GRID structures corresponding to the points in the energy profile (quadrature points) be outputted as the trajectory. For this task MAXP has to be adjusted to fit the requested number of grid-beads.
RTYP (integer; default 1): specifies the type of restraint used. At present only one type is implemented that is the absolute harmonic restraint with mass weighting. Must make sure that the RCS selection in the TREK call and that of the restraint are identical.
STEP (real; default 0.0): defines the uniform stepping parameter. Only this rather primitive option is currently available at present, but it seems to suffice.
RFRC (real; default 0.0): defines the force constant used in the restraint (might be restraint type dependent)
TRNC (integer; default MAXP): specifes the number of basis functions (sines) used in the Fourier series expansion. During path optimization, it is recommended to start with a small number of approximately MAXP/2 and gradually increase the number from there as desired. Using TRNC = MAXP during optimization will lead to path beats and possible instabilities.
To compute an accurate PMF profile, TRNC should be as big as possible, however it can never exceed the number of beads used to discretize the Fourier path.
IPRF (integer; default -1): requests that the energy profile along the qudrature grid-points be printed to the specified unit (file). Use IPRF 6 for the output file.
GRID (integer; default K*4): specifies the number of the grid point used for the quadrature and also for computing the energy profile and corresponding structures.
Remember to adjust the MAXP parameter to fit the requested number of grid-beads into the output trajectory. At present, failure to do so will result in unpredictable behavior of the executable, and at best a segmentation fault.
BEAD (integer; default K): allows to specify the number of beads in the path explicitly.

That is pretty much it in terms of the available options. Now let's see how this variety can be used to do some cool things like inserting more points along the path that are still equi-spaced. This could be used to restore some lost information if for some reason bead "i" information got lost or if the run failed (in dynamics one can always replace the failed restart file with one of the older ones or of the adjacent neighbor to propagate HFB optimization further). Also in collection this could be used to increase the density of points (number of points per unit length) such that even subtle transitions get resolved.