From 300e381633d0add1072ebf2fe0e7684489565a3f Mon Sep 17 00:00:00 2001 From: Noel Naughton Date: Thu, 4 Jan 2024 22:46:58 -0500 Subject: [PATCH] Update _rotations.py Make sure that arccos is never given a value outside of its allowable range. The previous version would occasionally yield a nan for stiff rods. --- elastica/_rotations.py | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/elastica/_rotations.py b/elastica/_rotations.py index ab8882848..81ee2d82e 100644 --- a/elastica/_rotations.py +++ b/elastica/_rotations.py @@ -146,9 +146,12 @@ def _inv_rotate(director_collection): ) ) - # TODO HARDCODED bugfix has to be changed. Remove 1e-14 tolerance - theta = arccos(0.5 * trace - 0.5 - 1e-10) + # Clip the trace to between -3 and 3. It any deviation beyond this is numerical error + trace = min(trace,3.0) + trace = max(trace,-3.0) + theta = arccos(0.5 * trace - 0.5 ) + # TODO HARDCODED bugfix has to be changed. Remove 1e-14 tolerance vector_collection[0, k] *= -0.5 * theta / sin(theta + 1e-14) vector_collection[1, k] *= -0.5 * theta / sin(theta + 1e-14) vector_collection[2, k] *= -0.5 * theta / sin(theta + 1e-14)