Eigenvector ambiguity - How to enforce a specific sign convention

Question

I am writing a program in Python that uses numpy.linalg.eigh to diagonalize a Hermitian matrix (a Hamiltonian). I diagonalize many such matrices and use the resultant eigenvector matrices for multiple unitary transformations of some other matrix. By "eigenvector matrix", I mean a matrix whose columns are the eigenvectors of the original matrix.

Unfortunately, I am hitting a potential problem because of the eigenvector sign ambiguity (i.e., eigenvectors are only defined up to a constant and normalization still does not fix the sign of an eigenvector). Specifically, the result I am calculating depends on the interference patterns produced by the successive unitary transformations. Thus, I anticipate that the sign ambiguity will become a problem.

My question:

What is the best way (or the industry standard) to enforce a particular sign convention for the eigenvectors?

I have thought of/come across the following:

1. Ensure the first coefficient of each eigenvector is positive. Problem: some of these coefficients are zero or within numerical error of zero.
2. Ensure the first coefficient of largest magnitude is positive. Problem: some of the eigenvectors have multiple coefficients with the same magnitude within numerical error. Numerical error then "randomly" determines which coefficient is "bigger."
3. Ensure the sum of the coefficients is positive. Problem: some coefficients are equal in magnitude but opposite in sign, leaving the sign still ambiguous/determined by numerical error. (I also see other problems with this approach).
4. Add a small number (such as 1E-16) to the eigenvector, ensure that the first coefficient is positive, then subtract the number. Problem: Maybe none important for me, but this makes me uneasy as I am not sure what problems it may cause.
5. (Inspired by Eigenshuffle and Sign correction in SVD and PCA) Pick a reference vector and ensure that the dot product of every eigenvector with this vector is positive. Problem: How to pick the vector? A random vector increases the likelihood that no eigenvectors are orthogonal to it (within numerical error), but there is no guarantee. Alternatively, one could choose a set of random vectors (all with positive coefficients) to increase the likelihood that the vector space is "spanned" well-enough.

I have tried to find what is the "standard" convention but I have a hard time finding anything particularly useful, particularly in Python. There is a solution for SVD (Sign correction in SVD and PCA), but I don't have any data vectors to compare to. There is Eigenshuffle (which is for Matlab and I am using Python), but my matrices are not usually successive small modifications of each other (though some are).

I am leaning toward solution 5 at it seems pretty intuitive; we are simply ensuring that all eigenvectors are in the same high-dimensional "quadrant". Also, having two or three random reference vectors with positive coefficients should cover almost all eigenvectors with very high probability, assuming the dimensionality of the system is not too big (my system has a dimensionality of 9).

Acea

If 1. and 4. Truely solves your problem I don't see why adding 1E-16 or treating 0 as true is a problem. By default numpy uses float64 and any changes to that will be insignificant if your algorithm is numerically stable, on the flip side, if your algorithm is not numerically stable, you'll run into problems no matter what you do. – Conclusive

What do you mean by "treating 0 as true"? Also, solution 1 does not work as many of my eigenvectors have initial coefficients of zero and I can't check in Python (ASFAIK) whether it is positive or negative zero. I presume my code is numerically stable, though I don't have a good intuition for when it may become unstable. It is basically a matrix multiplication of 20 to 30 unitary matrices. True, mapping any coefficient within 1E-16 of zero to +1E-16 should have a negligible effect. – Acea

what I mean is to Python, positive and negative zero doesn't matter at all. And if your code is numerically not stable, it won't matter what you do it will always blow up due to how float works. – Conclusive

Are you confident that all the eigenvalues of your Hamiltonian are different? If they are not you have even worse non-uniquess of the matrix of eigenvectors. – Barilla

@Barilla That is a good point. I am confident that some of the Hamiltonians have degenerate eigenvalues. (If some don't, then I have an even bigger problem than the non-uniqueness of the eigenvectors as the physics demands degeneracies in some cases.) Hmmm... I will need to test carefully whether my final result changes if I modify the eigenvector signs or take different linear combinations of degenerate eigenvectors. I have hope (but not confidence) that it won't matter. – Acea

@RockyLi Ok. I guess I can always just add a small number and see how much the result changes. – Acea

Also, even if I test the final code and find that the non-uniqueness of the matrix of eigenvectors doesn't matter, I think it would be useful to know about any standard or useful approaches to the eigenvector sign ambiguity problem. I had a hard time finding solid resources and I think collecting some knowledge together on one page would be useful to the community. – Acea