Numerical error caused math domain error when using Python math.acos()

Question

I use math.acos() to calculate the angle between two vectors. The idea is to calculate the dot product of the two normalized vectors, and use arccos of the dot product to return the angle of the two vectors (in the range of 0 to pi).

Occasionally the two vectors, when normalized, have a same direction and their dot product should be 1. However, because of the numerical errors, it is actually 0.999999999998, ... or sometimes 1.0000000000002. The latter kills the math.acos() with an error ValueError: math domain error

I accidentally found this question asking about a similar problem, but was closed. Here I am re-asking, hope to get a better idea on how to avoid such errors.

In my case, I have to check if the two vectors have the same direction before doing dot production and arccos. This helps, but I am still wondering if there are better ways to do this.

Answer 1

This is an old numerical problem encountered by all who worked with vectors and floating point arithmetic. A simple fix is to wrap the code that computes the angle in a function with some clamping logic:

import numpy as np
import math

def angle_between(vhat1, vhat2, tol = 5e-6):
    cosang = np.dot(vhat1, vhat2)
    if math.fabs(cosang) > 1.0:
        if math.fabs(cosang) - 1.0 < tol:
            cosang = math.modf(cosang)[1]
        else:
            raise ValueError('Invalid arguments (vectors not normalized?)')
    return math.acos(cosang)



w = np.array([ 1., 0., 0.])
z = np.array([-1., 0., 0.])
v = np.array([ math.sqrt(2.)/2., math.sqrt(2.)/2., 0.])
z1 = np.array([-1.0000001, 0., 0.])
w1 = np.array([9., 3., -5.])

print "{0:>5} deg".format(angle_between(v, w) * 180.0 / math.pi)
print "{0:>5} deg".format(angle_between(w, z) * 180.0 / math.pi)
print "{0:>5} deg".format(angle_between(w, z1) * 180.0 / math.pi)
# this last one will raise ValueError
print "{0:>5} deg".format(angle_between(w1, z1) * 180.0 / math.pi)

Output:

 45.0 deg
180.0 deg
180.0 deg

...
ValueError: Invalid arguments (vectors not normalized?)