|Gromacs|. However, these equations apply to the angle potential and the
improper dihedral potential as well.
-.. math::
-
- \begin{aligned}
- V_b &=&{\frac{1}{2}}\left[{(1-{\lambda})}k_b^A +
- {\lambda}k_b^B\right] \left[b - {(1-{\lambda})}b_0^A - {\lambda}b_0^B\right]^2 \\
- {\frac{\partial V_b}{\partial {\lambda}}}&=&{\frac{1}{2}}(k_b^B-k_b^A)
- \left[b - {(1-{\lambda})}b_0^A + {\lambda}b_0^B\right]^2 +
- \nonumber\\
- & & \phantom{{\frac{1}{2}}}(b_0^A-b_0^B) \left[b - {(1-{\lambda})}b_0^A -{\lambda}b_0^B\right]
- \left[{(1-{\lambda})}k_b^A + {\lambda}k_b^B \right]\end{aligned}
+.. math:: \begin{aligned}
+ V_b &=&{\frac{1}{2}}\left[{(1-{\lambda})}k_b^A +
+ {\lambda}k_b^B\right] \left[b - {(1-{\lambda})}b_0^A - {\lambda}b_0^B\right]^2 \\
+ {\frac{\partial V_b}{\partial {\lambda}}}&=&{\frac{1}{2}}(k_b^B-k_b^A)
+ \left[b - {(1-{\lambda})}b_0^A + {\lambda}b_0^B\right]^2 +
+ \nonumber\\
+ & & \phantom{{\frac{1}{2}}}(b_0^A-b_0^B) \left[b - {(1-{\lambda})}b_0^A -{\lambda}b_0^B\right]
+ \left[{(1-{\lambda})}k_b^A + {\lambda}k_b^B \right]\end{aligned}
+ :label: eqnfepharmpot
GROMOS-96 bonds and angles
~~~~~~~~~~~~~~~~~~~~~~~~~~
For the proper dihedrals, the equations are somewhat more complicated:
-.. math::
-
- \begin{aligned}
- V_d &=&\left[{(1-{\lambda})}k_d^A + {\lambda}k_d^B \right]
- \left( 1+ \cos\left[n_{\phi} \phi -
- {(1-{\lambda})}\phi_s^A - {\lambda}\phi_s^B
- \right]\right)\\
- {\frac{\partial V_d}{\partial {\lambda}}}&=&(k_d^B-k_d^A)
- \left( 1+ \cos
- \left[
- n_{\phi} \phi- {(1-{\lambda})}\phi_s^A - {\lambda}\phi_s^B
- \right]
- \right) +
- \nonumber\\
- &&(\phi_s^B - \phi_s^A) \left[{(1-{\lambda})}k_d^A - {\lambda}k_d^B\right]
- \sin\left[ n_{\phi}\phi - {(1-{\lambda})}\phi_s^A - {\lambda}\phi_s^B \right]\end{aligned}
+.. math:: \begin{aligned}
+ V_d &=&\left[{(1-{\lambda})}k_d^A + {\lambda}k_d^B \right]
+ \left( 1+ \cos\left[n_{\phi} \phi -
+ {(1-{\lambda})}\phi_s^A - {\lambda}\phi_s^B
+ \right]\right)\\
+ {\frac{\partial V_d}{\partial {\lambda}}}&=&(k_d^B-k_d^A)
+ \left( 1+ \cos
+ \left[
+ n_{\phi} \phi- {(1-{\lambda})}\phi_s^A - {\lambda}\phi_s^B
+ \right]
+ \right) +
+ \nonumber\\
+ &&(\phi_s^B - \phi_s^A) \left[{(1-{\lambda})}k_d^A - {\lambda}k_d^B\right]
+ \sin\left[ n_{\phi}\phi - {(1-{\lambda})}\phi_s^A - {\lambda}\phi_s^B \right]\end{aligned}
+ :label: eqnfeppropdihedral
**Note:** that the multiplicity :math:`n_{\phi}` can not be
parameterized because the function should remain periodic on the
For tabulated bonded interactions only the force constant can
interpolated:
-.. math::
-
- \begin{aligned}
- V &=& ({(1-{\lambda})}k^A + {\lambda}k^B) \, f \\
- {\frac{\partial V}{\partial {\lambda}}} &=& (k^B - k^A) \, f\end{aligned}
+.. math:: \begin{aligned}
+ V &=& ({(1-{\lambda})}k^A + {\lambda}k^B) \, f \\
+ {\frac{\partial V}{\partial {\lambda}}} &=& (k^B - k^A) \, f\end{aligned}
+ :label: eqnfeptabbonded
Coulomb interaction
~~~~~~~~~~~~~~~~~~~
The Coulomb interaction between two particles of which the charge varies
with :math:`{\lambda}` is:
-.. math::
-
- \begin{aligned}
- V_c &=& \frac{f}{{\varepsilon_{rf}}{r_{ij}}}\left[{(1-{\lambda})}q_i^A q_j^A + {\lambda}\, q_i^B q_j^B\right] \\
- {\frac{\partial V_c}{\partial {\lambda}}}&=& \frac{f}{{\varepsilon_{rf}}{r_{ij}}}\left[- q_i^A q_j^A + q_i^B q_j^B\right]\end{aligned}
+.. math:: \begin{aligned}
+ V_c &=& \frac{f}{{\varepsilon_{rf}}{r_{ij}}}\left[{(1-{\lambda})}q_i^A q_j^A + {\lambda}\, q_i^B q_j^B\right] \\
+ {\frac{\partial V_c}{\partial {\lambda}}}&=& \frac{f}{{\varepsilon_{rf}}{r_{ij}}}\left[- q_i^A q_j^A + q_i^B q_j^B\right]\end{aligned}
+ :label: eqnfepcoloumb
where :math:`f = \frac{1}{4\pi \varepsilon_0} = {138.935\,458}` (see
chapter :ref:`defunits`).
sign of :math:`{\frac{\partial E_k}{\partial {\lambda}}}` being
incorrect \ :ref:`99 <refGunsteren98a>`):
-.. math::
-
- \begin{aligned}
- E_k &=& {\frac{1}{2}}\frac{\mathbf{p}^2}{{(1-{\lambda})}m^A + {\lambda}m^B} \\
- {\frac{\partial E_k}{\partial {\lambda}}}&=& -{\frac{1}{2}}\frac{\mathbf{p}^2(m^B-m^A)}{({(1-{\lambda})}m^A + {\lambda}m^B)^2}\end{aligned}
+.. math:: \begin{aligned}
+ E_k &=& {\frac{1}{2}}\frac{\mathbf{p}^2}{{(1-{\lambda})}m^A + {\lambda}m^B} \\
+ {\frac{\partial E_k}{\partial {\lambda}}}&=& -{\frac{1}{2}}\frac{\mathbf{p}^2(m^B-m^A)}{({(1-{\lambda})}m^A + {\lambda}m^B)^2}\end{aligned}
+ :label: eqnfepekin
after taking the derivative, we *can* insert
-:math:`\mathbf{p}` =
-m :math:`\mathbf{v}`, such that:
+:math:`\mathbf{p}` = m :math:`\mathbf{v}`, such that:
.. math:: {\frac{\partial E_k}{\partial {\lambda}}}~=~ -{\frac{1}{2}}\mathbf{v}^2(m^B-m^A)
+ :label: eqnfepekinderivative
Constraints
~~~~~~~~~~~
:math:`k = 1 \ldots K` constraint equations :math:`g_k` for LINCS, then
.. math:: g_k = | \mathbf{r}_{k} | - d_{k}
+ :label: eqnfepconstr
where :math:`\mathbf{r}_k` is the displacement vector
between two particles and :math:`d_k` is the constraint distance between
has a :math:`{\lambda}` dependency by
.. math:: d_k = {(1-{\lambda})}d_{k}^A + {\lambda}d_k^B
+ :label: eqnfepconstrdistdep
Thus the :math:`{\lambda}`-dependent constraint equation is
.. math:: g_k = | \mathbf{r}_{k} | - \left({(1-{\lambda})}d_{k}^A + {\lambda}d_k^B\right).
+ :label: eqnfepconstrlambda
The (zero) contribution :math:`G` to the Hamiltonian from the
constraints (using Lagrange multipliers :math:`\lambda_k`, which are
logically distinct from the free-energy :math:`{\lambda}`) is
-.. math::
-
- \begin{aligned}
- G &=& \sum^K_k \lambda_k g_k \\
- {\frac{\partial G}{\partial {\lambda}}} &=& \frac{\partial G}{\partial d_k} {\frac{\partial d_k}{\partial {\lambda}}} \\
- &=& - \sum^K_k \lambda_k \left(d_k^B-d_k^A\right)\end{aligned}
+.. math:: \begin{aligned}
+ G &=& \sum^K_k \lambda_k g_k \\
+ {\frac{\partial G}{\partial {\lambda}}} &=& \frac{\partial G}{\partial d_k} {\frac{\partial d_k}{\partial {\lambda}}} \\
+ &=& - \sum^K_k \lambda_k \left(d_k^B-d_k^A\right)\end{aligned}
+ :label: eqnconstrfreeenergy
For SHAKE, the constraint equations are
.. math:: g_k = \mathbf{r}_{k}^2 - d_{k}^2
+ :label: eqnfepshakeconstr
with :math:`d_k` as before, so
-.. math::
-
- \begin{aligned}
- {\frac{\partial G}{\partial {\lambda}}} &=& -2 \sum^K_k \lambda_k \left(d_k^B-d_k^A\right)\end{aligned}
+.. math:: \begin{aligned}
+ {\frac{\partial G}{\partial {\lambda}}} &=& -2 \sum^K_k \lambda_k \left(d_k^B-d_k^A\right)\end{aligned}
+ :label: eqnfepshakeconstr2
Soft-core interactions
~~~~~~~~~~~~~~~~~~~~~~
of the regular potentials, so that the singularity in the potential and
its derivatives at :math:`r=0` is never reached:
-.. math::
-
- \begin{aligned}
- V_{sc}(r) &=& {(1-{\lambda})}V^A(r_A) + {\lambda}V^B(r_B)
- \\
- r_A &=& \left(\alpha \sigma_A^6 {\lambda}^p + r^6 \right)^\frac{1}{6}
- \\
- r_B &=& \left(\alpha \sigma_B^6 {(1-{\lambda})}^p + r^6 \right)^\frac{1}{6}\end{aligned}
+.. math:: \begin{aligned}
+ V_{sc}(r) &=& {(1-{\lambda})}V^A(r_A) + {\lambda}V^B(r_B)
+ \\
+ r_A &=& \left(\alpha \sigma_A^6 {\lambda}^p + r^6 \right)^\frac{1}{6}
+ \\
+ r_B &=& \left(\alpha \sigma_B^6 {(1-{\lambda})}^p + r^6 \right)^\frac{1}{6}\end{aligned}
+ :label: eqnfepsoftcore
where :math:`V^A` and :math:`V^B` are the normal “hard core” Van der
Waals or electrostatic potentials in state A (:math:`{\lambda}=0`) and
quickly switch the soft-core interaction to an almost constant value for
smaller :math:`r` (:numref:`Fig. %s <fig-softcore>`). The force is:
-.. math::
-
- F_{sc}(r) = -\frac{\partial V_{sc}(r)}{\partial r} =
- {(1-{\lambda})}F^A(r_A) \left(\frac{r}{r_A}\right)^5 +
- {\lambda}F^B(r_B) \left(\frac{r}{r_B}\right)^5
+.. math:: F_{sc}(r) = -\frac{\partial V_{sc}(r)}{\partial r} =
+ {(1-{\lambda})}F^A(r_A) \left(\frac{r}{r_A}\right)^5 +
+ {\lambda}F^B(r_B) \left(\frac{r}{r_B}\right)^5
+ :label: eqnfepsoftcoreforce
where :math:`F^A` and :math:`F^B` are the “hard core” forces. The
contribution to the derivative of the free energy is:
-.. math::
-
- \begin{aligned}
- {\frac{\partial V_{sc}(r)}{\partial {\lambda}}} & = &
- V^B(r_B) -V^A(r_A) +
- {(1-{\lambda})}\frac{\partial V^A(r_A)}{\partial r_A}
- \frac{\partial r_A}{\partial {\lambda}} +
- {\lambda}\frac{\partial V^B(r_B)}{\partial r_B}
- \frac{\partial r_B}{\partial {\lambda}}
- \nonumber\\
- &=&
- V^B(r_B) -V^A(r_A) + \nonumber \\
- & &
- \frac{p \alpha}{6}
- \left[ {\lambda}F^B(r_B) r^{-5}_B \sigma_B^6 {(1-{\lambda})}^{p-1} -
- {(1-{\lambda})}F^A(r_A) r^{-5}_A \sigma_A^6 {\lambda}^{p-1} \right]\end{aligned}
+.. math:: \begin{aligned}
+ {\frac{\partial V_{sc}(r)}{\partial {\lambda}}} & = &
+ V^B(r_B) -V^A(r_A) +
+ {(1-{\lambda})}\frac{\partial V^A(r_A)}{\partial r_A}
+ \frac{\partial r_A}{\partial {\lambda}} +
+ {\lambda}\frac{\partial V^B(r_B)}{\partial r_B}
+ \frac{\partial r_B}{\partial {\lambda}}
+ \nonumber\\
+ &=&
+ V^B(r_B) -V^A(r_A) + \nonumber \\
+ & &
+ \frac{p \alpha}{6}
+ \left[ {\lambda}F^B(r_B) r^{-5}_B \sigma_B^6 {(1-{\lambda})}^{p-1} -
+ {(1-{\lambda})}F^A(r_A) r^{-5}_A \sigma_A^6 {\lambda}^{p-1} \right]\end{aligned}
+ :label: eqnfepsoftcorederivative
The original GROMOS Lennard-Jones soft-core
function\ :ref:`100 <refBeutler94>` uses :math:`p=2`, but :math:`p=1` gives a smoother
the standard soft-core path described above \ :ref:`101 <refPham2011>`,
:ref:`102 <refPham2012>`. Specifically, we have:
-.. math::
-
- \begin{aligned}
- V_{sc}(r) &=& {(1-{\lambda})}V^A(r_A) + {\lambda}V^B(r_B)
- \\
- r_A &=& \left(\alpha \sigma_A^{48} {\lambda}^p + r^{48} \right)^\frac{1}{48}
- \\
- r_B &=& \left(\alpha \sigma_B^{48} {(1-{\lambda})}^p + r^{48} \right)^\frac{1}{48}\end{aligned}
+.. math:: \begin{aligned}
+ V_{sc}(r) &=& {(1-{\lambda})}V^A(r_A) + {\lambda}V^B(r_B)
+ \\
+ r_A &=& \left(\alpha \sigma_A^{48} {\lambda}^p + r^{48} \right)^\frac{1}{48}
+ \\
+ r_B &=& \left(\alpha \sigma_B^{48} {(1-{\lambda})}^p + r^{48} \right)^\frac{1}{48}\end{aligned}
+ :label: eqnnewsoftcore
This “1-1-48” path is also implemented in |Gromacs|. Note that for this
path the soft core :math:`\alpha` should satisfy
biod.pnpi.spb.ru / alexxy / gromacs.git / blobdiff