commit 79a554ffcf33156b81f3098e4e67821da867e2d6
author Sameer Agarwal <sameeragarwal@google.com> Fri Jan 13 11:37:27 2023 -0800
committer Sameer Agarwal <sameeragarwal@google.com> Fri Jan 13 11:37:27 2023 -0800
tree 8f026bfcb7487d4e13d8b09ca69b576dcedae707
parent f1113c08abec23462182a5a8235f8ab47277ec4d

Fix a bug in QuaternionRotatePoint.

In https://ceres-solver-review.git.corp.google.com/c/ceres-solver/+/23802

the computation of the norm of a quaternion

scale = 1/sqrt(q[0] * q[0] + q[1] * q[1] + q[2] * q[2] + q[3] * q[3]);

was replaced by

scale = 1/hypot(q[0], q[1], hypot(q[2], q[3]));

while this appear to be a more accurate computation because of the
use of hypot which can handle over and underflow it introduces a
bug for the case where q[2] = q[3] = 0.

While the hypot(q[2], q[3]) == 0 as scalars, if q[2] and q[3] are
jets, then the derivative will be NaN. Which means that even though
q[0] or q[1] is non-zero and the norm of the quaternion is non-zero,
and the resulting derivative is finite, this way of computing the
scale will produce nans in the derivative of scale.

The following quaternion will replicate the problem described above.

using Jet = ceres::Jet<double, 4>;
std::array<Jet, 4> quaternion = {Jet(1.0, 0), Jet(0.0, 1), Jet(0.0, 2), Jet(0.0, 3)};

This CL reverts the change to QuaternionRotatePoint and
adds a test for it.

Thanks to Jonathan Taylor for reproducing this bug.

Change-Id: I0fbbcc77d6945a38563d82efba4429f4b5278cd5