I have a very strange situation: I calculate angle between two vectors in spherical coordinate system. First I calculate scalar product of two vectors. Calculation using numpy array gives result with relative precision 10^(-16), as it should be for 64-bit floating numbers. But then I call acos and the angle from numpy‘s scalar product gives relative precision 10^(-11), i.e. 5 orders or magnitude worse compared to direct calculation (which is the correct result). Any ideas why?

#!/usr/bin/python
import numpy as np
import math

def testPrecision(theta1, phi1, theta2, phi2):

    # We define two vectors as numpy arrays
    t1 = np.array([math.sin(theta1) * math.cos(phi1), math.sin(theta1) * math.sin(phi1), math.cos(theta1)])
    t2 = np.array([math.sin(theta2) * math.cos(phi2), math.sin(theta2) * math.sin(phi2), math.cos(theta2)])

    # scalar product of vectors using numpy
    scalarProd1 = np.einsum('i,i', t1, t2)

    # scalar product of vectors without numpy
    scalarProd2 = math.sin(theta1) * math.cos(phi1) * math.sin(theta2) * math.cos(phi2) + math.sin(theta1) * math.sin(phi1) * math.sin(theta2) * math.sin(phi2) + math.cos(theta1) * math.cos(theta2)

    print("Scalar product:nnumpy = {:.15e} vs direct = {:.15e}".format(scalarProd1, scalarProd2))

    alpha1 = math.acos(scalarProd1)
    alpha2 = math.acos(scalarProd2)

    print("Angle:nnumpy = {:.15e} vs direct = {:.15e}".format(alpha1, alpha2))

testPrecision(theta1 = 1.745329251994330e-02, theta2 = 2.154759397163566e-02, phi1 = 1.745329251994329e-01, phi2 = 2.078144368653071e-01)