Hi Everyone,
This is a strange problem which may not be about Lapack, although it's showing up when using Lapack to invert a double precision real matrix. I am using:
Fedora Linux 2.6.30.9-96 on an I64 box
Lapack compiled with g77 3.4.6 20060404
Lapack routines DGETRF followed by DGETRI
Here are the calls to the Lapack routines:
call dgetrf(irank,irank,sscpD,irank,iScratch,ier)
call dgetri(irank,sscpD,irank,iScratch,Scrtch2,irank*16,ier)
The matrix has rank 640. It's a sum-of-squares and cross-products matrix for a statistics problem I know these routines are not as efficient for handling a symmetric positive definite matrix but that's not the problem. Here's the upper left corner of the product of the matrix with its inverse when I simply pass the matrix directly to the Lapack routines:
130157. 235105. 13758.4 -92677.5 35805.9 -13174.7
363250. -163466. 501326. 117770. 557132. -44305.9
-67196.9 -14119.3 -116038. 32273.7 -112083. 20780.7
98882.4 136193. 2502.16 39098.6 21843.1 -6672.41
54807.9 -91053.6 50633.2 32585.4 72267.4 13974.0
-214997. -43332.4 -198698. -130761. -252292. 5121.31
But if, instead, I cast the elements to single precision and then back to double precision, I get this very different and reasonable result, although it's not great. If I generate the matrix with a single precision routine, the inversion is even better behaved and with no funny business. But I really need to use the double precision version.
1.00009 0.714234E-09 -0.487233E-04 0.140948E-08 -0.159189E-03 0.289091E-08
0.424778E-06 1.00018 -0.536447E-07 0.195398E-03 0.217082E-07 0.541233E-03
0.531085E-03 0.567480E-08 0.999560 0.985579E-08 -0.915230E-03 0.231129E-07
-0.206659E-07 0.244289E-04 0.113778E-08 0.999989 0.512776E-07 -0.189541E-05
-0.215070E-03 -0.208492E-08 0.172455E-03 -0.373849E-08 1.00036 -0.870138E-08
0.182251E-07 -0.304040E-04 0.108860E-07 0.470498E-04 0.246131E-07 1.00009
Maybe there's a NaN in the double precision array or something, I don't know, wrong with some of the numbers. Any thoughts would be most welcome. Thanks.
Don

