Volatile results for numpy / scipy eigendecompositions

I find that scipy.linalg.eig sometimes gives inconsistent results. But not every time.

>>> import numpy as np
>>> import scipy.linalg as lin
>>> modmat=np.random.random((150,150))
>>> modmat=modmat+modmat.T  # the data i am interested in is described by real symmetric matrices
>>> d,v=lin.eig(modmat)
>>> dx=d.copy()
>>> vx=v.copy()
>>> d,v=lin.eig(modmat)
>>> np.all(d==dx)
False
>>> np.all(v==vx)
False
>>> e,w=lin.eigh(modmat)
>>> ex=e.copy()
>>> wx=w.copy()
>>> e,w=lin.eigh(modmat)
>>> np.all(e==ex)
True
>>> e,w=lin.eigh(modmat)
>>> np.all(e==ex)
False

Although I am not the greatest master of linear algebra, I understand that eigendecomposition is inherently prone to strange rounding errors, but I don’t understand why repeating the calculation will result in a different value. But my results and reproducibility are changing.

What exactly is the problem's problem - well, sometimes the results are acceptably different, and sometimes not. Here are some examples:

>>> d[1]
(9.8986888573772465+0j)
>>> dx[1]
(9.8986888573772092+0j)

The difference above ~ 3e-13 does not seem extremely complicated. Instead, the real problem (at least for my current project) is that some of the eigenvalues ​​cannot match their own sign.

>>> np.all(np.sign(d)==np.sign(dx))
False
>>> np.nonzero(np.sign(d)!=np.sign(dx))
(array([ 38,  39,  40,  41,  42,  45,  46,  47,  79,  80,  81,  82,  83,
    84, 109, 112]),)
>>> d[38]
(-6.4011617320002525+0j)
>>> dx[38]
(6.1888785138080209+0j)

MATLAB, , .