avoid crashes with 180 degree angles

This commit is contained in:
Axel Kohlmeyer
2024-11-14 07:59:17 -05:00
parent a3cd14bb82
commit 2b2ea9aaa2

View File

@ -402,14 +402,18 @@ void PairSWAngleTable::threebody_table(Param *paramij, Param *paramik, ParamTabl
rinv12 = 1.0/(r1*r2);
cs = (delr1[0]*delr2[0] + delr1[1]*delr2[1] + delr1[2]*delr2[2]) * rinv12;
var = acos(cs);
cs = MIN(-1.0,MAX(cs,1.0));
// look up energy (f(theta), ftheta) and force (df(theta)/dtheta, fprimetheta) at
// angle theta (var) in angle table belonging to parameter set paramijk
var = acos(cs);
uf_lookup(table_paramijk, var, ftheta, fprimetheta);
acosprime = 1.0 / (sqrt(1 - cs*cs ) );
if ((cs*cs - 1.0) > TINY)
acosprime = 1.0 / (sqrt(1 - cs*cs ) );
else
acosprime = 0.0;
facradtable = facexp*ftheta;
frad1table = facradtable*gsrainvsq1;