Blue-moon ensemble: Difference between revisions

Revision as of 15:53, 18 April 2022

In general, constrained molecular dynamics generates biased statistical averages. It can be shown that the correct average for a quantity $a(\xi )$ can be obtained using the formula:

a(\xi )={\frac {\langle |\mathbf {Z} |^{-1/2}a(\xi ^{*})\rangle _{\xi ^{*}}}{\langle |\mathbf {Z} |^{-1/2}\rangle _{\xi ^{*}}}},

where $\langle ...\rangle _{\xi ^{*}}$ stands for the statistical average of the quantity enclosed in angular parentheses computed for a constrained ensemble and $Z$ is a mass metric tensor defined as:

Z_{\alpha ,\beta }={\sum }_{i=1}^{3N}m_{i}^{-1}\nabla _{i}\xi _{\alpha }\cdot \nabla _{i}\xi _{\beta },\,\alpha =1,...,r,\,\beta =1,...,r,

It can be shown that the free energy gradient can be computed using the equation:^[1]^[2]^[3]^[4]

{\Bigl (}{\frac {\partial A}{\partial \xi _{k}}}{\Bigr )}_{\xi ^{*}}={\frac {1}{\langle |Z|^{-1/2}\rangle _{\xi ^{*}}}}\langle |Z|^{-1/2}[\lambda _{k}+{\frac {k_{B}T}{2|Z|}}\sum _{j=1}^{r}(Z^{-1})_{kj}\sum _{i=1}^{3N}m_{i}^{-1}\nabla _{i}\xi _{j}\cdot \nabla _{i}|Z|]\rangle _{\xi ^{*}},

where $\lambda _{\xi _{k}}$ is the Lagrange multiplier associated with the parameter ${\xi _{k}}$ used in the SHAKE algorithm.^[5]

The free-energy difference between states (1) and (2) can be computed by integrating the free-energy gradients over a connecting path:

{\Delta }A_{1\rightarrow 2}=\int _{\xi (1)}^{\xi (2)}{\Bigl (}{\frac {\partial {A}}{\partial \xi }}{\Bigr )}_{\xi ^{*}}\cdot d{\xi }.

Note that as the free-energy is a state quantity, the choice of path connecting (1) with (2) is irrelevant.

How to

The information needed to determine the blue moon ensemble averages within a Constrained molecular dynamics can be obtained by setting LBLUEOUT=.TRUE. The following output is written for each MD step in the file REPORT:

>Blue_moon

       lambda        |z|^(-1/2)    GkT           |z|^(-1/2)*(lambda+GkT)
  b_m>  0.585916E+01  0.215200E+02 -0.117679E+00  0.123556E+03

with the four numerical terms indicating $\lambda _{\xi _{k}}$ , $|Z|^{-1/2}$ , $\left({\frac {k_{B}T}{2|Z|}}\sum _{j=1}^{r}(Z^{-1})_{kj}\sum _{i=1}^{3N}m_{i}^{-1}\nabla _{i}\xi _{j}\cdot \nabla _{i}|Z|\right)$ , and $\left(|Z|^{-1/2}[\lambda _{k}+{\frac {k_{B}T}{2|Z|}}\sum _{j=1}^{r}(Z^{-1})_{kj}\sum _{i=1}^{3N}m_{i}^{-1}\nabla _{i}\xi _{j}\cdot \nabla _{i}|Z|]\right)$ , respectively. With this output, free energy gradients can conveniently be determined by equation given above as ratio between averages of the last and second numerical terms, which can be performed as follows:

grep b_m REPORT |awk 'BEGIN {a=0.;b=0.} {a+=$5;b+=$3} END {print a/b}'

References

[Carter89-1] E. A. Carter, G. Ciccotti, J. T. Hynes, and R. Kapral, Chem. Phys. Lett. 156, 472 (1989).

[Otter00-2] W. K. Den Otter and W. J. Briels, Mol. Phys. 98, 773 (2000).

[Darve02-3] E. Darve, M. A. Wilson, and A. Pohorille, Mol. Simul. 28, 113 (2002).

[Fleurat05-4] P. Fleurat-Lessard and T. Ziegler, J. Chem. Phys. 123, 084101 (2005).

[Ryckaert77-5] J. P. Ryckaert, G. Ciccotti, and H. J. C. Berendsen, J. Comp. Phys. 23, 327 (1977).

[1]

[2]

[3]

[4]

[5]

@@ Line 29: / Line 29: @@
 with the four numerical terms indicating <math>\lambda_{\xi_k}</math>, <math>|Z|^{-1/2}</math>, <math>\left ( \frac{k_B T}{2 |Z|} \sum_{j=1}^{r}(Z^{-1})_{kj} \sum_{i=1}^{3N} m_i^{-1}\nabla_i \xi_j \cdot \nabla_i |Z| \right ) </math>, and <math>\left ( |Z|^{-1/2} [\lambda_k +\frac{k_B T}{2 |Z|} \sum_{j=1}^{r}(Z^{-1})_{kj} \sum_{i=1}^{3N} m_i^{-1}\nabla_i \xi_j \cdot \nabla_i |Z|] \right ) </math>, respectively. With this output, free energy gradients can conveniently be determined by equation given above as ratio between averages of the last and second numerical terms, which can be performed as follows:
+ grep b_m REPORT |awk 'BEGIN {a=0.;b=0.} {a+=$5;b+=$3} END {print a/b}'
 ==References==