[gmx-developers] Virtual site out of the plane

[Mohammad Ghahremanpour]

> On Sep 28, 2017, at 10:51 PM, David van der Spoel <spoel at xray.bmc.uu.se> wrote:
> 
> 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?

In the manual, it has been described how a, b, and c are used in 3out, not how they are calculated. Pages 103 and 104. 
Are you referring to the same pages?

>> 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
> -- 
> Gromacs Developers mailing list
> 
> * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-developers_List before posting!
> 
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> 
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-developers or send a mail to gmx-developers-request at gromacs.org.