From 2b2ea9aaa2b7b9e34fe55e40d05bf144cc4cd449 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Thu, 14 Nov 2024 07:59:17 -0500 Subject: [PATCH] avoid crashes with 180 degree angles --- src/MANYBODY/pair_sw_angle_table.cpp | 10 +++++++--- 1 file changed, 7 insertions(+), 3 deletions(-) diff --git a/src/MANYBODY/pair_sw_angle_table.cpp b/src/MANYBODY/pair_sw_angle_table.cpp index 6e691afeb6..a985c31d53 100644 --- a/src/MANYBODY/pair_sw_angle_table.cpp +++ b/src/MANYBODY/pair_sw_angle_table.cpp @@ -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;