I've been trying to use the output of numpy.linalg.eigh as a
rotation matrix, but the matrix formed from the eigenvectors it
returns is improper; that is, it's not a pure rotation, but a
combination of a rotation and a reflection. You can see that as
follows.
GG = np.array([[2.0, 0.146, 0.0064],
[0.146, -1.0, 0.0003],
[0.0064, 0.0003, -1.0]])
vals,vecs = np.linalg.eigh(GG)
np.linalg.det(vecs)
This returns -0.99999999999. A pure rotation would return 1.0.
Am I missing something about the way this solver works? How can I
get eigh or some other function to return a proper rotation
matrix? Perhaps the fact that this matrix is traceless causes
problems?
Consider a square matrix A, and let us assume v is an eigenvector with e the corresponding eigenvalue, that is
A*v = e*v
Now consider what happens if we scale v by some scalar s, lets set w = v*s:
A*w = A*(v*s) = (A*v)*s = (e*v) = e*(v*s) = e*w
It turns out w is also an eigenvector for eigenvalue e! Actually eigenvectors are only defined up to scalar multiples, so the convention is that we normalize them to a length 1. So if you scale one eigenvector by -1 you still have a valid eigensystem, this time with determinant +0.999...9.
Finally, why is it 0.999...9 instead of 1? We are using floating point arithmetic, which is not exact. If you are trying to solve numerical problems I'd recommend familiarizing yourself with floating point numbers, their problems and limitations and the best practices.
Evidently