Fix precision issues and a bug in vec_roll_to_mat3_normalized.

When the input vector gets close to -Y, y and theta becomes totally
unreliable. It is thus necessary to compute the result in a different
way based on x and z. The code already had a special case, but:

- The threshold for using the special case was way too low.
- The special case was not precise enough to extend the threshold.
- The special case math had a sign error, resulting in a jump.

This adds tests for the computation precision and fixes the issues
by adjusting the threshold, and replacing the special case with one
based on a quadratic Taylor expansion of sqrt instead of linear.

Replacing the special case fixes the bug and results in a compatibility
break, requiring versioning for the roll of affected bones.

Differential Revision: https://developer.blender.org/D9551

Split and extend unit tests for vec_roll_to_mat3_normalized.

Separate the huge test into huge logical parts and add more cases
to check. Also add a utility to check that the matrix is orthogonal,
with arbitrary epsilon values and calculations in double.

A couple of tests deliberately fail, to be fixed in following commits.

Ref D9551

Fix T82455: vec_roll_to_mat3_normalized returns NaN when nor close to -Y.

In this case theta is completely unsafe to use, so a different
threshold based on x and z has to be used to avoid division by zero.

Ref D9551