numarray.linear_algebra.eigenvectors bug ?

[Eugene Druker]

Hi all,
I need eigen values and vectors for symmetric matrices (like VCV). 
Solving with numarray and testing the results, I've got 
strange results - input and output matrices of some sizes are 
very different. That's an example of differencies:

Size   MaxDif(E.Val)    MaxDif(E.Vect)
   2   0                1.110223e-016
   3   1.776357e-015    2.053913e-015
  10   1.221245e-015    7.255307e-014
  11   3.552714e-015    0.8455322
  12   1.526557e-016    1.345035e-013
  13   8.881784e-016    0.7796273
  14   2.775558e-016    0.7504375
  15   2.331468e-015    2.045031e-013

This is for numarray 0.7 on W2K.
I hoped to have all the differencies below 1e-15 or so.
What's wrong ? Program text is below.

Thanks,
Eugene Druker


import numarray, numarray.linear_algebra

EigenVectors = numarray.linear_algebra.eigenvectors
Multiply = numarray.matrixmultiply
Zeros = numarray.zeros
Transpose = numarray.transpose

def makeCovMat(size):
    covmat = Zeros((size,size),'f8')
    d,t = 0.9,0.8
    for i in xrange(size):
        v = covmat[i,i] = d
        for j in xrange(i+1,size):
            v *= t
            covmat[i,j] = covmat[j,i] = v
    return covmat

print 'Size   MaxDif(E.Val)    MaxDif(E.Vect)'
for size in [2,3]+range(10,24):
    covmat = makeCovMat(size)
    eval,evec = EigenVectors(covmat)
    # Test S = V L V'
    mdiag = Zeros((size,size),'f8') # L
    for j in xrange(size):
        mdiag[j,j] = eval[j]
    testmat = Multiply(Transpose(evec),Multiply(mdiag,evec))
    ival,ivec = EigenVectors(testmat)
    print ' %3d   %-14.7g   %-14.7g' % \
          (size, (eval-ival).max(), (evec-ivec).max())