_misc.py 6.1 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191
  1. import numpy as np
  2. from numpy.linalg import LinAlgError
  3. from .blas import get_blas_funcs
  4. from .lapack import get_lapack_funcs
  5. __all__ = ['LinAlgError', 'LinAlgWarning', 'norm']
  6. class LinAlgWarning(RuntimeWarning):
  7. """
  8. The warning emitted when a linear algebra related operation is close
  9. to fail conditions of the algorithm or loss of accuracy is expected.
  10. """
  11. pass
  12. def norm(a, ord=None, axis=None, keepdims=False, check_finite=True):
  13. """
  14. Matrix or vector norm.
  15. This function is able to return one of eight different matrix norms,
  16. or one of an infinite number of vector norms (described below), depending
  17. on the value of the ``ord`` parameter. For tensors with rank different from
  18. 1 or 2, only `ord=None` is supported.
  19. Parameters
  20. ----------
  21. a : array_like
  22. Input array. If `axis` is None, `a` must be 1-D or 2-D, unless `ord`
  23. is None. If both `axis` and `ord` are None, the 2-norm of
  24. ``a.ravel`` will be returned.
  25. ord : {int, inf, -inf, 'fro', 'nuc', None}, optional
  26. Order of the norm (see table under ``Notes``). inf means NumPy's
  27. `inf` object.
  28. axis : {int, 2-tuple of ints, None}, optional
  29. If `axis` is an integer, it specifies the axis of `a` along which to
  30. compute the vector norms. If `axis` is a 2-tuple, it specifies the
  31. axes that hold 2-D matrices, and the matrix norms of these matrices
  32. are computed. If `axis` is None then either a vector norm (when `a`
  33. is 1-D) or a matrix norm (when `a` is 2-D) is returned.
  34. keepdims : bool, optional
  35. If this is set to True, the axes which are normed over are left in the
  36. result as dimensions with size one. With this option the result will
  37. broadcast correctly against the original `a`.
  38. check_finite : bool, optional
  39. Whether to check that the input matrix contains only finite numbers.
  40. Disabling may give a performance gain, but may result in problems
  41. (crashes, non-termination) if the inputs do contain infinities or NaNs.
  42. Returns
  43. -------
  44. n : float or ndarray
  45. Norm of the matrix or vector(s).
  46. Notes
  47. -----
  48. For values of ``ord <= 0``, the result is, strictly speaking, not a
  49. mathematical 'norm', but it may still be useful for various numerical
  50. purposes.
  51. The following norms can be calculated:
  52. ===== ============================ ==========================
  53. ord norm for matrices norm for vectors
  54. ===== ============================ ==========================
  55. None Frobenius norm 2-norm
  56. 'fro' Frobenius norm --
  57. 'nuc' nuclear norm --
  58. inf max(sum(abs(a), axis=1)) max(abs(a))
  59. -inf min(sum(abs(a), axis=1)) min(abs(a))
  60. 0 -- sum(a != 0)
  61. 1 max(sum(abs(a), axis=0)) as below
  62. -1 min(sum(abs(a), axis=0)) as below
  63. 2 2-norm (largest sing. value) as below
  64. -2 smallest singular value as below
  65. other -- sum(abs(a)**ord)**(1./ord)
  66. ===== ============================ ==========================
  67. The Frobenius norm is given by [1]_:
  68. :math:`||A||_F = [\\sum_{i,j} abs(a_{i,j})^2]^{1/2}`
  69. The nuclear norm is the sum of the singular values.
  70. Both the Frobenius and nuclear norm orders are only defined for
  71. matrices.
  72. References
  73. ----------
  74. .. [1] G. H. Golub and C. F. Van Loan, *Matrix Computations*,
  75. Baltimore, MD, Johns Hopkins University Press, 1985, pg. 15
  76. Examples
  77. --------
  78. >>> import numpy as np
  79. >>> from scipy.linalg import norm
  80. >>> a = np.arange(9) - 4.0
  81. >>> a
  82. array([-4., -3., -2., -1., 0., 1., 2., 3., 4.])
  83. >>> b = a.reshape((3, 3))
  84. >>> b
  85. array([[-4., -3., -2.],
  86. [-1., 0., 1.],
  87. [ 2., 3., 4.]])
  88. >>> norm(a)
  89. 7.745966692414834
  90. >>> norm(b)
  91. 7.745966692414834
  92. >>> norm(b, 'fro')
  93. 7.745966692414834
  94. >>> norm(a, np.inf)
  95. 4
  96. >>> norm(b, np.inf)
  97. 9
  98. >>> norm(a, -np.inf)
  99. 0
  100. >>> norm(b, -np.inf)
  101. 2
  102. >>> norm(a, 1)
  103. 20
  104. >>> norm(b, 1)
  105. 7
  106. >>> norm(a, -1)
  107. -4.6566128774142013e-010
  108. >>> norm(b, -1)
  109. 6
  110. >>> norm(a, 2)
  111. 7.745966692414834
  112. >>> norm(b, 2)
  113. 7.3484692283495345
  114. >>> norm(a, -2)
  115. 0
  116. >>> norm(b, -2)
  117. 1.8570331885190563e-016
  118. >>> norm(a, 3)
  119. 5.8480354764257312
  120. >>> norm(a, -3)
  121. 0
  122. """
  123. # Differs from numpy only in non-finite handling and the use of blas.
  124. if check_finite:
  125. a = np.asarray_chkfinite(a)
  126. else:
  127. a = np.asarray(a)
  128. if a.size and a.dtype.char in 'fdFD' and axis is None and not keepdims:
  129. if ord in (None, 2) and (a.ndim == 1):
  130. # use blas for fast and stable euclidean norm
  131. nrm2 = get_blas_funcs('nrm2', dtype=a.dtype, ilp64='preferred')
  132. return nrm2(a)
  133. if a.ndim == 2:
  134. # Use lapack for a couple fast matrix norms.
  135. # For some reason the *lange frobenius norm is slow.
  136. lange_args = None
  137. # Make sure this works if the user uses the axis keywords
  138. # to apply the norm to the transpose.
  139. if ord == 1:
  140. if np.isfortran(a):
  141. lange_args = '1', a
  142. elif np.isfortran(a.T):
  143. lange_args = 'i', a.T
  144. elif ord == np.inf:
  145. if np.isfortran(a):
  146. lange_args = 'i', a
  147. elif np.isfortran(a.T):
  148. lange_args = '1', a.T
  149. if lange_args:
  150. lange = get_lapack_funcs('lange', dtype=a.dtype, ilp64='preferred')
  151. return lange(*lange_args)
  152. # fall back to numpy in every other case
  153. return np.linalg.norm(a, ord=ord, axis=axis, keepdims=keepdims)
  154. def _datacopied(arr, original):
  155. """
  156. Strict check for `arr` not sharing any data with `original`,
  157. under the assumption that arr = asarray(original)
  158. """
  159. if arr is original:
  160. return False
  161. if not isinstance(original, np.ndarray) and hasattr(original, '__array__'):
  162. return False
  163. return arr.base is None