by Julien Langou » Mon Dec 12, 2016 9:49 am
Hi Ilya,
Your test code looks good to me.
Please remember that DPOCON and DGECON are condition number estimators. What matters (in general) for the condition number is the order of magnitude. Is the condition number 1, 10, 100, 1000, 10000, etc. ? Whether the condition number is 14 or 16 is not that important (in general). In any case, DPOCON returns an estimate of the condition number. It does not return the condition number. I quote Higham and Tisseur [2]: "The [...] LAPACK estimators both produce estimates that in practice are almost always within a factor [...] 3 [...], of the quantities they are estimating [3], [4], [5]." So in your case an answer between 5 and 48 would be deemed to be acceptable. The original algorithm is from Hager [1].
Note 1: If you want to compute the exact condition number, you have to compute inv(A) explicitly, and then computing the 1-norm. This is (much) more expensive than POCON.
Note 2: I have no idea why using POCON on inv(A) returns a very accurate answer, while POCON on A returns a not so accurate answer.
Julien.
[1] William W. Hager, "Condition Estimates," SIAM J. Sci. Stat. Comput. 5, 1984, 311-316, 1984.
[2] Nicholas J. Higham and Françoise Tisseur, "A Block Algorithm for Matrix 1-Norm Estimation with an Application to 1-Norm Pseudospectra, "SIAM J. Matrix Anal. Appl., Vol. 21, 1185-1201, 2000.
[3] N. J. Higham, A survey of condition number estimation for triangular matrices, SIAM Rev., 29 (1987), pp. 575–596.
[4] N. J. Higham, FORTRAN codes for estimating the one-norm of a real or complex matrix, with applications to condition estimation (Algorithm 674), ACM Trans. Math. Software, 14 (1988), pp. 381–396.
[5] N. J. Higham, Experience with a matrix norm estimator, SIAM J. Sci. Stat. Comput., 11 (1990), pp. 804–809.