The HQBM Module of CHARMM
By Emanuele Paci, 1997/2000
HQBM is an external perturbation designed induce conformational
changes in macromolecules. The time dependent perturbation is designed
to introduce a very small perturbation to the short time dynamics of
the system and does not affect the conservation of the constants of
motion of the system (the conservation of the total energy or of the
suitable conserved quantity when an extended Lagrangian is used can
then be used as a check of the correctness of the forces).
The external perturbation needs:
- a reference (or target) structure
- a reaction coordinate which defines a "distance" from the
reference structure
* Menu:
* Syntax:: Syntax of the HQBM command
* Function:: Purpose of each of the keywords
* Input:: HQBM Input Description
[INPUT HQBM command]
- read the reference structure
OPEN UNIT 1 READ FORMATTED NAME coor0.crd
READ COOR CARD COMP UNIT 1
CLOSE UNIT 1
- call the perturbation choosing a coupling constant [ALPHA], a
reaction coordinate [RC1, RC2 or RC3], and a selection of atoms
which define the reaction coordinate.
HQBM [AWAY] ALPHA real [RC1] atom-selection [READLIST integer] [IUNJ integer]
[EVAL real] [TARGet real]
- energy has to be called immediately after HQBM !!
ENERGY
The following section describes the keywords of the HQBM command.
HQBM introduces a half quadratic perturbation on a given reaction
coordinate (see below)
Meaning of the HQBM parameters
1) AWAY drive the system away from the reference structure.
As an example, if the reaction coordinate measures the deviation from
a reference conformation, the perturbation will increase it.
2) ALPHA is the force constant of the half harmonic potential.
3) RC1, RC2, RC3, RC4 will select other reaction coordinates (which
can be user-redefined)
4) atom-selection: only the atom selected, for both the actual and
reference structure, enter in the above defined reaction coordinate.
Example: if the reference structure has all the coordinates
equal to zero, and only two atoms are selected, the reaction
coordinate is equal to the square of the distance between the two
atoms.
5) In alternative to a atom-selection a file can be read from unit
READLIST which looks like this one
n
i1 j1
.....
in jn
In this case the reaction coordinate is given by the sum of n terms
(r_{ij}(t)-r_{ij}^R)^2 where ij are the pairs given in the list
6) IUNJ: write the output (istep rc(t) max(rc)) on unit IUNJ
7) EVAL: specify the exponent for the RC1 formula. Note that if value
2 is specified with EVAL keyword square root of the original
RC1 formula is used. If you dont specify EVAL keyword the original
is used. When EVAL is specified the general formula is:
\rho(t)=(\frac{1}{N(N-1)}\sum_{i=1}^N\sum_{j\ne i}^N
(r_{ij}(t)-r_{ij}^R)^{eval})^\frac{1}{eval}
8) TARGet: one can specify the target value for the \rho (default = 0)
As an example, RC1 corresponds to
\rho(t)=\frac{1}{N(N-1)}\sum_{i=1}^N\sum_{j\ne i}^N(r_{ij}(t)-r_{ij}^R)^2
where r_{ij}=r_i-r_j and r_i^R are the coordinates of the reference structure.
If the reference structure has all the coordinates zero, the reaction
coordinate will be the radius of gyration. If the reference structure
is the experimental structure the reaction coordinate will be a
measure of the instantaneous "distance" from that, analogous to mean
square deviation.
The method is described in
E. Paci and M. Karplus. Forced unfolding of fibronectin type 3
modules: An analysis by biased molecular dynamics simulations.
J. Mol. Biol., 288: 441-459, 1999.
Input example for RC4:
open read unit 10 card name pairs.dat
coor helix sele type ca end
hqbm alpha 10.0 rc4 readlist 10 eval 2 -
xhelix ?xaxi yhelix ?yaxi zhelix ?zaxi
NIH/DCRT/Laboratory for Structural Biology
FDA/CBER/OVRR Biophysics Laboratory
Modified, updated and generalized by C.L. Brooks, III
The Scripps Research Institute