graph.py 8.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279
  1. from sympy.utilities.iterables import \
  2. flatten, connected_components, strongly_connected_components
  3. from .common import NonSquareMatrixError
  4. def _connected_components(M):
  5. """Returns the list of connected vertices of the graph when
  6. a square matrix is viewed as a weighted graph.
  7. Examples
  8. ========
  9. >>> from sympy import Matrix
  10. >>> A = Matrix([
  11. ... [66, 0, 0, 68, 0, 0, 0, 0, 67],
  12. ... [0, 55, 0, 0, 0, 0, 54, 53, 0],
  13. ... [0, 0, 0, 0, 1, 2, 0, 0, 0],
  14. ... [86, 0, 0, 88, 0, 0, 0, 0, 87],
  15. ... [0, 0, 10, 0, 11, 12, 0, 0, 0],
  16. ... [0, 0, 20, 0, 21, 22, 0, 0, 0],
  17. ... [0, 45, 0, 0, 0, 0, 44, 43, 0],
  18. ... [0, 35, 0, 0, 0, 0, 34, 33, 0],
  19. ... [76, 0, 0, 78, 0, 0, 0, 0, 77]])
  20. >>> A.connected_components()
  21. [[0, 3, 8], [1, 6, 7], [2, 4, 5]]
  22. Notes
  23. =====
  24. Even if any symbolic elements of the matrix can be indeterminate
  25. to be zero mathematically, this only takes the account of the
  26. structural aspect of the matrix, so they will considered to be
  27. nonzero.
  28. """
  29. if not M.is_square:
  30. raise NonSquareMatrixError
  31. V = range(M.rows)
  32. E = sorted(M.todok().keys())
  33. return connected_components((V, E))
  34. def _strongly_connected_components(M):
  35. """Returns the list of strongly connected vertices of the graph when
  36. a square matrix is viewed as a weighted graph.
  37. Examples
  38. ========
  39. >>> from sympy import Matrix
  40. >>> A = Matrix([
  41. ... [44, 0, 0, 0, 43, 0, 45, 0, 0],
  42. ... [0, 66, 62, 61, 0, 68, 0, 60, 67],
  43. ... [0, 0, 22, 21, 0, 0, 0, 20, 0],
  44. ... [0, 0, 12, 11, 0, 0, 0, 10, 0],
  45. ... [34, 0, 0, 0, 33, 0, 35, 0, 0],
  46. ... [0, 86, 82, 81, 0, 88, 0, 80, 87],
  47. ... [54, 0, 0, 0, 53, 0, 55, 0, 0],
  48. ... [0, 0, 2, 1, 0, 0, 0, 0, 0],
  49. ... [0, 76, 72, 71, 0, 78, 0, 70, 77]])
  50. >>> A.strongly_connected_components()
  51. [[0, 4, 6], [2, 3, 7], [1, 5, 8]]
  52. """
  53. if not M.is_square:
  54. raise NonSquareMatrixError
  55. # RepMatrix uses the more efficient DomainMatrix.scc() method
  56. rep = getattr(M, '_rep', None)
  57. if rep is not None:
  58. return rep.scc()
  59. V = range(M.rows)
  60. E = sorted(M.todok().keys())
  61. return strongly_connected_components((V, E))
  62. def _connected_components_decomposition(M):
  63. """Decomposes a square matrix into block diagonal form only
  64. using the permutations.
  65. Explanation
  66. ===========
  67. The decomposition is in a form of $A = P^{-1} B P$ where $P$ is a
  68. permutation matrix and $B$ is a block diagonal matrix.
  69. Returns
  70. =======
  71. P, B : PermutationMatrix, BlockDiagMatrix
  72. *P* is a permutation matrix for the similarity transform
  73. as in the explanation. And *B* is the block diagonal matrix of
  74. the result of the permutation.
  75. If you would like to get the diagonal blocks from the
  76. BlockDiagMatrix, see
  77. :meth:`~sympy.matrices.expressions.blockmatrix.BlockDiagMatrix.get_diag_blocks`.
  78. Examples
  79. ========
  80. >>> from sympy import Matrix, pprint
  81. >>> A = Matrix([
  82. ... [66, 0, 0, 68, 0, 0, 0, 0, 67],
  83. ... [0, 55, 0, 0, 0, 0, 54, 53, 0],
  84. ... [0, 0, 0, 0, 1, 2, 0, 0, 0],
  85. ... [86, 0, 0, 88, 0, 0, 0, 0, 87],
  86. ... [0, 0, 10, 0, 11, 12, 0, 0, 0],
  87. ... [0, 0, 20, 0, 21, 22, 0, 0, 0],
  88. ... [0, 45, 0, 0, 0, 0, 44, 43, 0],
  89. ... [0, 35, 0, 0, 0, 0, 34, 33, 0],
  90. ... [76, 0, 0, 78, 0, 0, 0, 0, 77]])
  91. >>> P, B = A.connected_components_decomposition()
  92. >>> pprint(P)
  93. PermutationMatrix((1 3)(2 8 5 7 4 6))
  94. >>> pprint(B)
  95. [[66 68 67] ]
  96. [[ ] ]
  97. [[86 88 87] 0 0 ]
  98. [[ ] ]
  99. [[76 78 77] ]
  100. [ ]
  101. [ [55 54 53] ]
  102. [ [ ] ]
  103. [ 0 [45 44 43] 0 ]
  104. [ [ ] ]
  105. [ [35 34 33] ]
  106. [ ]
  107. [ [0 1 2 ]]
  108. [ [ ]]
  109. [ 0 0 [10 11 12]]
  110. [ [ ]]
  111. [ [20 21 22]]
  112. >>> P = P.as_explicit()
  113. >>> B = B.as_explicit()
  114. >>> P.T*B*P == A
  115. True
  116. Notes
  117. =====
  118. This problem corresponds to the finding of the connected components
  119. of a graph, when a matrix is viewed as a weighted graph.
  120. """
  121. from sympy.combinatorics.permutations import Permutation
  122. from sympy.matrices.expressions.blockmatrix import BlockDiagMatrix
  123. from sympy.matrices.expressions.permutation import PermutationMatrix
  124. iblocks = M.connected_components()
  125. p = Permutation(flatten(iblocks))
  126. P = PermutationMatrix(p)
  127. blocks = []
  128. for b in iblocks:
  129. blocks.append(M[b, b])
  130. B = BlockDiagMatrix(*blocks)
  131. return P, B
  132. def _strongly_connected_components_decomposition(M, lower=True):
  133. """Decomposes a square matrix into block triangular form only
  134. using the permutations.
  135. Explanation
  136. ===========
  137. The decomposition is in a form of $A = P^{-1} B P$ where $P$ is a
  138. permutation matrix and $B$ is a block diagonal matrix.
  139. Parameters
  140. ==========
  141. lower : bool
  142. Makes $B$ lower block triangular when ``True``.
  143. Otherwise, makes $B$ upper block triangular.
  144. Returns
  145. =======
  146. P, B : PermutationMatrix, BlockMatrix
  147. *P* is a permutation matrix for the similarity transform
  148. as in the explanation. And *B* is the block triangular matrix of
  149. the result of the permutation.
  150. Examples
  151. ========
  152. >>> from sympy import Matrix, pprint
  153. >>> A = Matrix([
  154. ... [44, 0, 0, 0, 43, 0, 45, 0, 0],
  155. ... [0, 66, 62, 61, 0, 68, 0, 60, 67],
  156. ... [0, 0, 22, 21, 0, 0, 0, 20, 0],
  157. ... [0, 0, 12, 11, 0, 0, 0, 10, 0],
  158. ... [34, 0, 0, 0, 33, 0, 35, 0, 0],
  159. ... [0, 86, 82, 81, 0, 88, 0, 80, 87],
  160. ... [54, 0, 0, 0, 53, 0, 55, 0, 0],
  161. ... [0, 0, 2, 1, 0, 0, 0, 0, 0],
  162. ... [0, 76, 72, 71, 0, 78, 0, 70, 77]])
  163. A lower block triangular decomposition:
  164. >>> P, B = A.strongly_connected_components_decomposition()
  165. >>> pprint(P)
  166. PermutationMatrix((8)(1 4 3 2 6)(5 7))
  167. >>> pprint(B)
  168. [[44 43 45] [0 0 0] [0 0 0] ]
  169. [[ ] [ ] [ ] ]
  170. [[34 33 35] [0 0 0] [0 0 0] ]
  171. [[ ] [ ] [ ] ]
  172. [[54 53 55] [0 0 0] [0 0 0] ]
  173. [ ]
  174. [ [0 0 0] [22 21 20] [0 0 0] ]
  175. [ [ ] [ ] [ ] ]
  176. [ [0 0 0] [12 11 10] [0 0 0] ]
  177. [ [ ] [ ] [ ] ]
  178. [ [0 0 0] [2 1 0 ] [0 0 0] ]
  179. [ ]
  180. [ [0 0 0] [62 61 60] [66 68 67]]
  181. [ [ ] [ ] [ ]]
  182. [ [0 0 0] [82 81 80] [86 88 87]]
  183. [ [ ] [ ] [ ]]
  184. [ [0 0 0] [72 71 70] [76 78 77]]
  185. >>> P = P.as_explicit()
  186. >>> B = B.as_explicit()
  187. >>> P.T * B * P == A
  188. True
  189. An upper block triangular decomposition:
  190. >>> P, B = A.strongly_connected_components_decomposition(lower=False)
  191. >>> pprint(P)
  192. PermutationMatrix((0 1 5 7 4 3 2 8 6))
  193. >>> pprint(B)
  194. [[66 68 67] [62 61 60] [0 0 0] ]
  195. [[ ] [ ] [ ] ]
  196. [[86 88 87] [82 81 80] [0 0 0] ]
  197. [[ ] [ ] [ ] ]
  198. [[76 78 77] [72 71 70] [0 0 0] ]
  199. [ ]
  200. [ [0 0 0] [22 21 20] [0 0 0] ]
  201. [ [ ] [ ] [ ] ]
  202. [ [0 0 0] [12 11 10] [0 0 0] ]
  203. [ [ ] [ ] [ ] ]
  204. [ [0 0 0] [2 1 0 ] [0 0 0] ]
  205. [ ]
  206. [ [0 0 0] [0 0 0] [44 43 45]]
  207. [ [ ] [ ] [ ]]
  208. [ [0 0 0] [0 0 0] [34 33 35]]
  209. [ [ ] [ ] [ ]]
  210. [ [0 0 0] [0 0 0] [54 53 55]]
  211. >>> P = P.as_explicit()
  212. >>> B = B.as_explicit()
  213. >>> P.T * B * P == A
  214. True
  215. """
  216. from sympy.combinatorics.permutations import Permutation
  217. from sympy.matrices.expressions.blockmatrix import BlockMatrix
  218. from sympy.matrices.expressions.permutation import PermutationMatrix
  219. iblocks = M.strongly_connected_components()
  220. if not lower:
  221. iblocks = list(reversed(iblocks))
  222. p = Permutation(flatten(iblocks))
  223. P = PermutationMatrix(p)
  224. rows = []
  225. for a in iblocks:
  226. cols = []
  227. for b in iblocks:
  228. cols.append(M[a, b])
  229. rows.append(cols)
  230. B = BlockMatrix(rows)
  231. return P, B