[gmx-developers] Virtual site out of the plane

[David van der Spoel]

On 28/09/17 14:59, Mohammad Ghahremanpour wrote:
> Hello developers,
> 
> The calc_vsite3out_param() function in vsite_parm.cpp calculates 
> the weights (called a, b, c in the code) to locate the virtual site (i) 
> out of the k-j-l plane later on in vsite.cpp.
> 
> I copied some part of the code here:
> 
> pijk = std::cos(aijk)*bij;
> pijl = std::cos(aijl)*bij;
> a    = ( pijk + (pijk*std::cos(akjl)-pijl) * std::cos(akjl) / 
> gmx::square(std::sin(akjl)) ) / bjk;
> b    = ( pijl + (pijl*std::cos(akjl)-pijk) * std::cos(akjl) / 
> gmx::square(std::sin(akjl)) ) / bjl;
> c    = -std::sqrt( gmx::square(bij) -
>                             ( gmx::square(pijk) - 
> 2*pijk*pijl*std::cos(akjl) + gmx::square(pijl) )
>                             / gmx::square(std::sin(akjl)) )
>              / ( bjk*bjl*std::sin(akjl) );
> 
> There are a lot of projections going on here which make the code 
> confusing without any documentation, unfortunately.
> Is there some reference for the way that a, b, and c are calculated here?
Hate to suggest it but have you checked the manual?
> 
> BTW, what if the number under the SQRT is negative? Then c will be NAN. 
> There is no checking to make sure that we are not calculating the SQRT 
> of a negative number.
> 
> Best,
> Mohammad
> 
> 


-- 
David van der Spoel, Ph.D., Professor of Biology
Head of Department, Cell & Molecular Biology, Uppsala University.
Box 596, SE-75124 Uppsala, Sweden. Phone: +46184714205.
http://www.icm.uu.se