[gmx-developers] Virtual site out of the plane

[Mohammad Ghahremanpour]

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? 

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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-developers/attachments/20170928/4d11b6df/attachment.html>