[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

bug in EIGENQL ?



Fellow IDLers:

I encountered the following strange behavior when using the IDL (v 5.2
under HP-UX 10.20) routine EIGENQL, which is supposed to compute the
eigenvectors and eigenvalues of a real, symmetric matrix.  The error
is reproducible on my system and seems to occur the first time EIGENQL
is called during an IDL session but not thereafter.


-----------------------------------------------------------------------------------
First I read in my data file (contents appears at bottom of this article)...

IDL> openr, lunr, "t", /get_lun
IDL> means = fltarr(7)
IDL> covar = fltarr(7,7)
IDL> readf, lunr, means
IDL> readf, lunr, covar
IDL> free_lun, lunr

-----------------------------------------------------------------------------------
I verify that I read what I wanted to... 

IDL> print, means
      196.000      162.000      194.300      191.700      165.400      193.200      174.900

IDL> print, covar
      277.700      294.900      284.200      288.300      296.000      242.200      236.800
      294.900      373.000      304.300      314.000      359.400      267.700      281.600
      284.200      304.300      293.500      304.200      312.400      259.300      253.300
      288.300      314.000      304.200      340.600      348.000      307.700      297.700
      296.000      359.400      312.400      348.000      385.500      314.300      320.400
      242.200      267.700      259.300      307.700      314.300      321.900      306.100
      236.800      281.600      253.300      297.700      320.400      306.100      305.300

-----------------------------------------------------------------------------------
Now I attempt to compute the eigenvectors of the matrix.... 

IDL> evals = eigenql(covar,eigenvectors=evecs)

-----------------------------------------------------------------------------------
The following diagnostic message occurs only the first time EIGENQL is
called.  I don't know what it means.

% Restored file: EIGENQL.

-----------------------------------------------------------------------------------
Now I print the results.  The eigenvalues are correct, but notice the spurious
zeros in the first row and last column of the eigenvector matrix.  In
fact, the whole matrix is erroneous.

IDL> print, evals
      2106.65      118.913      51.9952      14.3094      3.72573      1.56117     0.345064

IDL> print, evecs
      0.00000      0.00000      0.00000      0.00000      0.00000      0.00000      1.00000
   -0.0203527    -0.142684    -0.475245 -4.46735e-05     0.748321    -0.439755      0.00000
    0.0638729     0.732782     0.405909     0.256386     0.128762    -0.460299      0.00000
     0.133166    -0.413454     0.573429    -0.546051    0.0375968    -0.427687      0.00000
    -0.761143     0.134580    -0.242504    -0.282152    -0.362916    -0.363901      0.00000
   -0.0365828    -0.482871    0.0594889     0.715859    -0.292982    -0.404557      0.00000
     0.630155     0.143010    -0.467131    -0.209839    -0.452190    -0.340196      0.00000

-----------------------------------------------------------------------------------
Now I run EIGENQL a second time with identical parameters.  This time,
the results for both evals and evecs are correct.....


IDL> evals = eigenql(covar,eigenvectors=evecs)

IDL> print, evals                             
      2106.65      118.913      51.9952      14.3094      3.72573      1.56117     0.345064

IDL> print, evecs                             
     0.344896     0.395337     0.361339     0.395519     0.420594     0.362555     0.359734
     0.366553     0.450309     0.285925   -0.0748501    0.0173607    -0.588532    -0.478365
    -0.426247     0.580345    -0.400150    -0.314331     0.376371    -0.157386     0.236990
    -0.293417    -0.320451   -0.0798885     0.555279     0.577378    -0.307122    -0.262319
     0.256999    -0.395998     0.163642    -0.310046     0.235583    -0.498243     0.592017
     0.227300    -0.211230  -0.00957581    -0.541312     0.540564     0.390869    -0.406961
     0.600993   0.00285527    -0.770914     0.204148  -0.00832092   -0.0454444    0.0260860


Why the erroneous results on the first call?

		- Grant

==================== contents of input file "t" ==========================


 196.0     162.0     194.3     191.7     165.4     193.2     174.9    

 277.7     294.9     284.2     288.3     296.0     242.2     236.8    
 294.9     373.0     304.3     314.0     359.4     267.7     281.6    
 284.2     304.3     293.5     304.2     312.4     259.3     253.3    
 288.3     314.0     304.2     340.6     348.0     307.7     297.7    
 296.0     359.4     312.4     348.0     385.5     314.3     320.4    
 242.2     267.7     259.3     307.7     314.3     321.9     306.1    
 236.8     281.6     253.3     297.7     320.4     306.1     305.3    
-- 
**** Reply to: gpetty@purdue.edu, NOT address in header! *****
Grant W. Petty                        |Assoc. Prof.,  Atmospheric Science
Dept. of Earth & Atmospheric Sciences |Voice:  (765)-494-2544
Purdue University, West Lafayette IN  |Fax:    (765)-496-1210