Warning: That file was not part of the compilation database. It may have many parsing errors.

1 | /* ========================================================================== */ |
---|---|

2 | /* === BTF package ========================================================== */ |

3 | /* ========================================================================== */ |

4 | |

5 | /* BTF_MAXTRANS: find a column permutation Q to give A*Q a zero-free diagonal |

6 | * BTF_STRONGCOMP: find a symmetric permutation P to put P*A*P' into block |

7 | * upper triangular form. |

8 | * BTF_ORDER: do both of the above (btf_maxtrans then btf_strongcomp). |

9 | * |

10 | * By Tim Davis. Copyright (c) 2004-2007, University of Florida. |

11 | * with support from Sandia National Laboratories. All Rights Reserved. |

12 | */ |

13 | |

14 | |

15 | /* ========================================================================== */ |

16 | /* === BTF_MAXTRANS ========================================================= */ |

17 | /* ========================================================================== */ |

18 | |

19 | /* BTF_MAXTRANS: finds a permutation of the columns of a matrix so that it has a |

20 | * zero-free diagonal. The input is an m-by-n sparse matrix in compressed |

21 | * column form. The array Ap of size n+1 gives the starting and ending |

22 | * positions of the columns in the array Ai. Ap[0] must be zero. The array Ai |

23 | * contains the row indices of the nonzeros of the matrix A, and is of size |

24 | * Ap[n]. The row indices of column j are located in Ai[Ap[j] ... Ap[j+1]-1]. |

25 | * Row indices must be in the range 0 to m-1. Duplicate entries may be present |

26 | * in any given column. The input matrix is not checked for validity (row |

27 | * indices out of the range 0 to m-1 will lead to an undeterminate result - |

28 | * possibly a core dump, for example). Row indices in any given column need |

29 | * not be in sorted order. However, if they are sorted and the matrix already |

30 | * has a zero-free diagonal, then the identity permutation is returned. |

31 | * |

32 | * The output of btf_maxtrans is an array Match of size n. If row i is matched |

33 | * with column j, then A(i,j) is nonzero, and then Match[i] = j. If the matrix |

34 | * is structurally nonsingular, all entries in the Match array are unique, and |

35 | * Match can be viewed as a column permutation if A is square. That is, column |

36 | * k of the original matrix becomes column Match[k] of the permuted matrix. In |

37 | * MATLAB, this can be expressed as (for non-structurally singular matrices): |

38 | * |

39 | * Match = maxtrans (A) ; |

40 | * B = A (:, Match) ; |

41 | * |

42 | * except of course here the A matrix and Match vector are all 0-based (rows |

43 | * and columns in the range 0 to n-1), not 1-based (rows/cols in range 1 to n). |

44 | * The MATLAB dmperm routine returns a row permutation. See the maxtrans |

45 | * mexFunction for more details. |

46 | * |

47 | * If row i is not matched to any column, then Match[i] is == -1. The |

48 | * btf_maxtrans routine returns the number of nonzeros on diagonal of the |

49 | * permuted matrix. |

50 | * |

51 | * In the MATLAB mexFunction interface to btf_maxtrans, 1 is added to the Match |

52 | * array to obtain a 1-based permutation. Thus, in MATLAB where A is m-by-n: |

53 | * |

54 | * q = maxtrans (A) ; % has entries in the range 0:n |

55 | * q % a column permutation (only if sprank(A)==n) |

56 | * B = A (:, q) ; % permuted matrix (only if sprank(A)==n) |

57 | * sum (q > 0) ; % same as "sprank (A)" |

58 | * |

59 | * This behaviour differs from p = dmperm (A) in MATLAB, which returns the |

60 | * matching as p(j)=i if row i and column j are matched, and p(j)=0 if column j |

61 | * is unmatched. |

62 | * |

63 | * p = dmperm (A) ; % has entries in the range 0:m |

64 | * p % a row permutation (only if sprank(A)==m) |

65 | * B = A (p, :) ; % permuted matrix (only if sprank(A)==m) |

66 | * sum (p > 0) ; % definition of sprank (A) |

67 | * |

68 | * This algorithm is based on the paper "On Algorithms for obtaining a maximum |

69 | * transversal" by Iain Duff, ACM Trans. Mathematical Software, vol 7, no. 1, |

70 | * pp. 315-330, and "Algorithm 575: Permutations for a zero-free diagonal", |

71 | * same issue, pp. 387-390. Algorithm 575 is MC21A in the Harwell Subroutine |

72 | * Library. This code is not merely a translation of the Fortran code into C. |

73 | * It is a completely new implementation of the basic underlying method (depth |

74 | * first search over a subgraph with nodes corresponding to columns matched so |

75 | * far, and cheap matching). This code was written with minimal observation of |

76 | * the MC21A/B code itself. See comments below for a comparison between the |

77 | * maxtrans and MC21A/B codes. |

78 | * |

79 | * This routine operates on a column-form matrix and produces a column |

80 | * permutation. MC21A uses a row-form matrix and produces a row permutation. |

81 | * The difference is merely one of convention in the comments and interpretation |

82 | * of the inputs and outputs. If you want a row permutation, simply pass a |

83 | * compressed-row sparse matrix to this routine and you will get a row |

84 | * permutation (just like MC21A). Similarly, you can pass a column-oriented |

85 | * matrix to MC21A and it will happily return a column permutation. |

86 | */ |

87 | |

88 | #ifndef _BTF_H |

89 | #define _BTF_H |

90 | |

91 | /* make it easy for C++ programs to include BTF */ |

92 | #ifdef __cplusplus |

93 | extern "C"{ |

94 | #endif |

95 | |

96 | #include "SuiteSparse_config.h" |

97 | |

98 | int btf_maxtrans /* returns # of columns matched */ |

99 | ( |

100 | /* --- input, not modified: --- */ |

101 | int nrow, /* A is nrow-by-ncol in compressed column form */ |

102 | int ncol, |

103 | int Ap [ ], /* size ncol+1 */ |

104 | int Ai [ ], /* size nz = Ap [ncol] */ |

105 | double maxwork, /* maximum amount of work to do is maxwork*nnz(A); no limit |

106 | * if <= 0 */ |

107 | |

108 | /* --- output, not defined on input --- */ |

109 | double *work, /* work = -1 if maxwork > 0 and the total work performed |

110 | * reached the maximum of maxwork*nnz(A). |

111 | * Otherwise, work = the total work performed. */ |

112 | |

113 | int Match [ ], /* size nrow. Match [i] = j if column j matched to row i |

114 | * (see above for the singular-matrix case) */ |

115 | |

116 | /* --- workspace, not defined on input or output --- */ |

117 | int Work [ ] /* size 5*ncol */ |

118 | ) ; |

119 | |

120 | /* long integer version (all "int" parameters become "SuiteSparse_long") */ |

121 | SuiteSparse_long btf_l_maxtrans (SuiteSparse_long, SuiteSparse_long, |

122 | SuiteSparse_long *, SuiteSparse_long *, double, double *, |

123 | SuiteSparse_long *, SuiteSparse_long *) ; |

124 | |

125 | |

126 | /* ========================================================================== */ |

127 | /* === BTF_STRONGCOMP ======================================================= */ |

128 | /* ========================================================================== */ |

129 | |

130 | /* BTF_STRONGCOMP finds the strongly connected components of a graph, returning |

131 | * a symmetric permutation. The matrix A must be square, and is provided on |

132 | * input in compressed-column form (see BTF_MAXTRANS, above). The diagonal of |

133 | * the input matrix A (or A*Q if Q is provided on input) is ignored. |

134 | * |

135 | * If Q is not NULL on input, then the strongly connected components of A*Q are |

136 | * found. Q may be flagged on input, where Q[k] < 0 denotes a flagged column k. |

137 | * The permutation is j = BTF_UNFLIP (Q [k]). On output, Q is modified (the |

138 | * flags are preserved) so that P*A*Q is in block upper triangular form. |

139 | * |

140 | * If Q is NULL, then the permutation P is returned so that P*A*P' is in upper |

141 | * block triangular form. |

142 | * |

143 | * The vector R gives the block boundaries, where block b is in rows/columns |

144 | * R[b] to R[b+1]-1 of the permuted matrix, and where b ranges from 1 to the |

145 | * number of strongly connected components found. |

146 | */ |

147 | |

148 | int btf_strongcomp /* return # of strongly connected components */ |

149 | ( |

150 | /* input, not modified: */ |

151 | int n, /* A is n-by-n in compressed column form */ |

152 | int Ap [ ], /* size n+1 */ |

153 | int Ai [ ], /* size nz = Ap [n] */ |

154 | |

155 | /* optional input, modified (if present) on output: */ |

156 | int Q [ ], /* size n, input column permutation */ |

157 | |

158 | /* output, not defined on input */ |

159 | int P [ ], /* size n. P [k] = j if row and column j are kth row/col |

160 | * in permuted matrix. */ |

161 | |

162 | int R [ ], /* size n+1. block b is in rows/cols R[b] ... R[b+1]-1 */ |

163 | |

164 | /* workspace, not defined on input or output */ |

165 | int Work [ ] /* size 4n */ |

166 | ) ; |

167 | |

168 | SuiteSparse_long btf_l_strongcomp (SuiteSparse_long, SuiteSparse_long *, |

169 | SuiteSparse_long *, SuiteSparse_long *, SuiteSparse_long *, |

170 | SuiteSparse_long *, SuiteSparse_long *) ; |

171 | |

172 | |

173 | /* ========================================================================== */ |

174 | /* === BTF_ORDER ============================================================ */ |

175 | /* ========================================================================== */ |

176 | |

177 | /* BTF_ORDER permutes a square matrix into upper block triangular form. It |

178 | * does this by first finding a maximum matching (or perhaps a limited matching |

179 | * if the work is limited), via the btf_maxtrans function. If a complete |

180 | * matching is not found, BTF_ORDER completes the permutation, but flags the |

181 | * columns of P*A*Q to denote which columns are not matched. If the matrix is |

182 | * structurally rank deficient, some of the entries on the diagonal of the |

183 | * permuted matrix will be zero. BTF_ORDER then calls btf_strongcomp to find |

184 | * the strongly-connected components. |

185 | * |

186 | * On output, P and Q are the row and column permutations, where i = P[k] if |

187 | * row i of A is the kth row of P*A*Q, and j = BTF_UNFLIP(Q[k]) if column j of |

188 | * A is the kth column of P*A*Q. If Q[k] < 0, then the (k,k)th entry in P*A*Q |

189 | * is structurally zero. |

190 | * |

191 | * The vector R gives the block boundaries, where block b is in rows/columns |

192 | * R[b] to R[b+1]-1 of the permuted matrix, and where b ranges from 1 to the |

193 | * number of strongly connected components found. |

194 | */ |

195 | |

196 | int btf_order /* returns number of blocks found */ |

197 | ( |

198 | /* --- input, not modified: --- */ |

199 | int n, /* A is n-by-n in compressed column form */ |

200 | int Ap [ ], /* size n+1 */ |

201 | int Ai [ ], /* size nz = Ap [n] */ |

202 | double maxwork, /* do at most maxwork*nnz(A) work in the maximum |

203 | * transversal; no limit if <= 0 */ |

204 | |

205 | /* --- output, not defined on input --- */ |

206 | double *work, /* return value from btf_maxtrans */ |

207 | int P [ ], /* size n, row permutation */ |

208 | int Q [ ], /* size n, column permutation */ |

209 | int R [ ], /* size n+1. block b is in rows/cols R[b] ... R[b+1]-1 */ |

210 | int *nmatch, /* # nonzeros on diagonal of P*A*Q */ |

211 | |

212 | /* --- workspace, not defined on input or output --- */ |

213 | int Work [ ] /* size 5n */ |

214 | ) ; |

215 | |

216 | SuiteSparse_long btf_l_order (SuiteSparse_long, SuiteSparse_long *, |

217 | SuiteSparse_long *, double , double *, SuiteSparse_long *, |

218 | SuiteSparse_long *, SuiteSparse_long *, SuiteSparse_long *, |

219 | SuiteSparse_long *) ; |

220 | |

221 | |

222 | /* ========================================================================== */ |

223 | /* === BTF marking of singular columns ====================================== */ |

224 | /* ========================================================================== */ |

225 | |

226 | /* BTF_FLIP is a "negation about -1", and is used to mark an integer j |

227 | * that is normally non-negative. BTF_FLIP (-1) is -1. BTF_FLIP of |

228 | * a number > -1 is negative, and BTF_FLIP of a number < -1 is positive. |

229 | * BTF_FLIP (BTF_FLIP (j)) = j for all integers j. UNFLIP (j) acts |

230 | * like an "absolute value" operation, and is always >= -1. You can test |

231 | * whether or not an integer j is "flipped" with the BTF_ISFLIPPED (j) |

232 | * macro. |

233 | */ |

234 | |

235 | #define BTF_FLIP(j) (-(j)-2) |

236 | #define BTF_ISFLIPPED(j) ((j) < -1) |

237 | #define BTF_UNFLIP(j) ((BTF_ISFLIPPED (j)) ? BTF_FLIP (j) : (j)) |

238 | |

239 | /* ========================================================================== */ |

240 | /* === BTF version ========================================================== */ |

241 | /* ========================================================================== */ |

242 | |

243 | /* All versions of BTF include these definitions. |

244 | * As an example, to test if the version you are using is 1.2 or later: |

245 | * |

246 | * if (BTF_VERSION >= BTF_VERSION_CODE (1,2)) ... |

247 | * |

248 | * This also works during compile-time: |

249 | * |

250 | * #if (BTF >= BTF_VERSION_CODE (1,2)) |

251 | * printf ("This is version 1.2 or later\n") ; |

252 | * #else |

253 | * printf ("This is an early version\n") ; |

254 | * #endif |

255 | */ |

256 | |

257 | #define BTF_DATE "May 4, 2016" |

258 | #define BTF_VERSION_CODE(main,sub) ((main) * 1000 + (sub)) |

259 | #define BTF_MAIN_VERSION 1 |

260 | #define BTF_SUB_VERSION 2 |

261 | #define BTF_SUBSUB_VERSION 6 |

262 | #define BTF_VERSION BTF_VERSION_CODE(BTF_MAIN_VERSION,BTF_SUB_VERSION) |

263 | |

264 | #ifdef __cplusplus |

265 | } |

266 | #endif |

267 | #endif |

268 |

Warning: That file was not part of the compilation database. It may have many parsing errors.