Coverage for src/cvx/linalg/core/condition.py: 100%

1"""Condition-number helpers: the NaN-aware ``cond`` and the ill-conditioning warning checks.""" 

2 

3from __future__ import annotations 

4 

5import warnings 

6from typing import Literal 

7 

8import numpy as np 

9 

10from .exceptions import IllConditionedMatrixWarning 

11from .types import Matrix 

12 

13try: # SciPy lets check_and_warn_condition skip the SVD for well-conditioned SPD input. 

14 from scipy.linalg import cho_factor as _cho_factor # type: ignore[import-untyped] 

15 from scipy.linalg.lapack import get_lapack_funcs as _get_lapack_funcs # type: ignore[import-untyped] 

16 

17 _HAVE_SCIPY = True 

18except ImportError: # pragma: no cover - depends on the environment; the fallback is tested by patching the flag 

19 _HAVE_SCIPY = False 

20 

21_SCREEN_MARGIN: float = 10.0 

22"""Safety factor on the LAPACK condition estimate before it may stand in for the exact 2-norm value.""" 

23 

24DEFAULT_COND_THRESHOLD: float = 1e12 

25"""Default condition-number threshold above which an IllConditionedMatrixWarning is emitted.""" 

26 

27 

28def cond(matrix: Matrix, p: int | float | Literal["fro", "nuc"] | None = None) -> float: 

29 """Return the condition number of a matrix. 

30 

31 Returns ``nan`` if the matrix contains any non-finite (NaN or inf) entries. 

32 Otherwise delegates to :func:`numpy.linalg.cond`. 

33 

34 Args: 

35 matrix: Input matrix. 

36 p: Order of the norm used to compute the condition number. 

37 Accepts the same values as :func:`numpy.linalg.cond` 

38 (``None``, ``1``, ``-1``, ``2``, ``-2``, ``numpy.inf``, 

39 ``-numpy.inf``, ``'fro'``). Defaults to ``None`` which 

40 corresponds to the 2-norm (largest singular value divided by 

41 the smallest). 

42 

43 Returns: 

44 The condition number as a ``float``, or ``nan`` when the matrix 

45 contains non-finite entries. 

46 

47 Examples: 

48 >>> import numpy as np 

49 >>> cond(np.eye(3)) 

50 1.0 

51 >>> import math 

52 >>> math.isnan(cond(np.array([[float('nan'), 1.0], [1.0, 2.0]]))) 

53 True 

54 >>> cond(np.diag([1.0, 1e10]), p=1) 

55 10000000000.0 

56 """ 

57 if not np.all(np.isfinite(matrix)): 

58 return float("nan") 

59 return float(np.linalg.cond(matrix, p=p)) 

60 

61 

62def warn_ill_conditioned(cond_value: float, threshold: float | None, stacklevel: int = 3) -> None: 

63 """Emit IllConditionedMatrixWarning when *cond_value* exceeds *threshold*. 

64 

65 Args: 

66 cond_value: Condition number to compare against the threshold. 

67 threshold: Upper bound before a warning is issued, or ``None`` to never warn. 

68 stacklevel: Stack level passed to :func:`warnings.warn` so the warning 

69 points at the caller of the public API. Defaults to ``3``. 

70 

71 Example: 

72 >>> import warnings 

73 >>> with warnings.catch_warnings(record=True) as w: 

74 ... warnings.simplefilter("always") 

75 ... warn_ill_conditioned(2.0, 0.5) 

76 ... len(w) 

77 1 

78 """ 

79 if threshold is not None and cond_value > threshold: 

80 warnings.warn( 

81 f"Matrix condition number {cond_value:.3e} exceeds threshold {threshold:.3e}; " 

82 "results may be numerically unreliable.", 

83 IllConditionedMatrixWarning, 

84 stacklevel=stacklevel, 

85 ) 

86 

87 

88def check_and_warn_condition(matrix: Matrix, threshold: float | None) -> None: 

89 """Emit IllConditionedMatrixWarning when the condition number exceeds threshold. 

90 

91 Args: 

92 matrix: Square matrix whose condition number is checked. 

93 threshold: Upper bound before a warning is issued, or ``None`` to skip 

94 the check entirely. 

95 

96 The condition number is the 2-norm one, as from :func:`cond`. With SciPy 

97 installed, a symmetric positive-definite matrix is first screened by a 

98 LAPACK estimate from its Cholesky factor; when that proves the matrix well 

99 within *threshold*, the full SVD is skipped. 

100 

101 Example: 

102 >>> import numpy as np 

103 >>> import warnings 

104 >>> with warnings.catch_warnings(record=True) as w: 

105 ... warnings.simplefilter("always") 

106 ... check_and_warn_condition(np.eye(2), 0.5) 

107 ... len(w) 

108 1 

109 """ 

110 if threshold is None or _certainly_within(matrix, threshold): 

111 return 

112 warn_ill_conditioned(cond(matrix), threshold, stacklevel=4) 

113 

114 

115def _certainly_within(matrix: Matrix, threshold: float) -> bool: 

116 """Return True when a cheap estimate proves the 2-norm condition number is at most *threshold*. 

117 

118 For a symmetric matrix ``||A||_2 <= ||A||_1``, so the 1-norm condition 

119 number bounds the 2-norm one from above. LAPACK ``pocon`` estimates the 

120 former from a Cholesky factor for about the cost of one more solve, whereas 

121 the exact 2-norm value needs a full SVD. The estimate can fall short of the 

122 true 1-norm value, hence the :data:`_SCREEN_MARGIN` safety factor. 

123 

124 ``False`` means only "not proven": SciPy is missing, or the matrix is not 

125 symmetric positive-definite, or the estimate is too close to the threshold. 

126 The caller then computes the exact condition number. 

127 

128 Args: 

129 matrix: Square matrix whose condition number is screened. 

130 threshold: Upper bound the condition number is tested against. 

131 

132 Returns: 

133 ``True`` if the 2-norm condition number is certainly at most *threshold*. 

134 """ 

135 if not _HAVE_SCIPY or not np.array_equal(matrix, matrix.T): 

136 return False 

137 try: 

138 factor, lower = _cho_factor(matrix) 

139 except (np.linalg.LinAlgError, ValueError): # not positive-definite, or not finite 

140 return False 

141 (pocon,) = _get_lapack_funcs(("pocon",), (factor,)) 

142 rcond, _info = pocon(factor, np.linalg.norm(matrix, 1), uplo="L" if lower else "U") 

143 return bool(rcond * threshold >= _SCREEN_MARGIN)