from math import sqrt
from cvxopt import matrix
-from cvxopt.lapack import gees, syevr
+from cvxopt.lapack import gees, gesdd, syevr
from . import options
"""
return matrix(mat, (len(mat), 1))
+
+
+def condition_number(mat):
+ """
+ Return the condition number of the given matrix.
+
+ The condition number of a matrix quantifies how hard it is to do
+ numerical computation with that matrix. It is usually defined as
+ the ratio of the norm of the matrix to the norm of its inverse, and
+ therefore depends on the norm used. One way to compute the condition
+ number with respect to the 2-norm is as the ratio of the matrix's
+ largest and smallest singular values. Since we have easy access to
+ those singular values, that is the algorithm we use.
+
+ The larger the condition number is, the worse the matrix is.
+
+ Parameters
+ ----------
+ mat : matrix
+ The matrix whose condition number you want.
+
+ Returns
+ -------
+
+ float
+ The nonnegative condition number of ``mat``.
+
+ Examples
+ --------
+
+ >>> condition_number(identity(1, typecode='d'))
+ 1.0
+ >>> condition_number(identity(2, typecode='d'))
+ 1.0
+ >>> condition_number(identity(3, typecode='d'))
+ 1.0
+
+ >>> A = matrix([[2,1],[1,2]], tc='d')
+ >>> abs(condition_number(A) - 3.0) < options.ABS_TOL
+ True
+
+ >>> A = matrix([[2,1j],[-1j,2]], tc='z')
+ >>> abs(condition_number(A) - 3.0) < options.ABS_TOL
+ True
+
+ """
+ num_eigs = min(mat.size)
+ eigs = matrix(0, (num_eigs,1), tc='d')
+ gesdd(mat, eigs)
+
+ if len(eigs) > 0:
+ return eigs[0]/eigs[-1]
+ else:
+ return 0