[gmx-developers] Question on SHAKE equations

[David van der Spoel]

On 2013-06-07 19:26, "Pablo García Risueño" wrote:
> Dear Gromacs developers
>
> I am studying the Gromacs implementation of the SHAKE algorithm. As far as
> I know, it corresponds to the cshake procedure, which lies in the shakef.c
> file, and the equation (5.6) of the SHAKE paper
> (see http://www.sciencedirect.com/science/article/pii/0021999177900985) is
> essentially in
>
> acor = omega*diff*m2[ll]/rrpr
> lagr[ll] += acor;
>
> where the relationship between variables is:
>
> SHAKE var               Gromacs var
>                            omega=1
> d^2 - (r')^2               diff
> 1/(2(1/m_i+1/mj))          m2[ll]
> \vec{r}_{ij}\cdot\vec{r}'_{ij}  rrpr
> g                         lagr[ll]
>
> being ll the constraint index. However, I think that the term proportional
> to g^2 in the equation (5.6) of SHAKE paper is not in the Gromacs code.
> This term is relatively small, but I was assuming that it is not
> negligible. Could some developer explain me whether this term is
> neglected, or it is included in the calculations elsewhere?
>

I don't have the complete overview, so I will assume your observation is 
correct.
Since the algorithm is iterative the second term may not be needed, 
although it might be converge in fewer iterations with it.

A practical issue is that not many people use shake, since LINCS works 
in parallel too.

> Thank you very much. Best regards
>
>
>
> --
>
> Dr. Pablo García Risueño
>
> Institut für Physik und IRIS Adlershof, Humboldt Universität zu Berlin,
> Zum Grossen Windkanal 6, 12489 Berlin, Germany
>
> Tel. +49 030 209366369
>


-- 
David van der Spoel, Ph.D., Professor of Biology
Dept. of Cell & Molec. Biol., Uppsala University.
Box 596, 75124 Uppsala, Sweden. Phone:	+46184714205.
spoel at xray.bmc.uu.se    http://folding.bmc.uu.se