1 /* ========================================================================== */ 2 /* === Include/cholmod_cholesky.h =========================================== */ 3 /* ========================================================================== */ 4 5 /* ----------------------------------------------------------------------------- 6 * CHOLMOD/Include/cholmod_cholesky.h. Copyright (C) 2005-2013, Timothy A. Davis 7 * http://www.suitesparse.com 8 * -------------------------------------------------------------------------- */ 9 10 /* CHOLMOD Cholesky module. 11 * 12 * Sparse Cholesky routines: analysis, factorization, and solve. 13 * 14 * The primary routines are all that a user requires to order, analyze, and 15 * factorize a sparse symmetric positive definite matrix A (or A*A'), and 16 * to solve Ax=b (or A*A'x=b). The primary routines rely on the secondary 17 * routines, the CHOLMOD Core module, and the AMD and COLAMD packages. They 18 * make optional use of the CHOLMOD Supernodal and Partition modules, the 19 * METIS package, and the CCOLAMD package. 20 * 21 * Primary routines: 22 * ----------------- 23 * 24 * cholmod_analyze order and analyze (simplicial or supernodal) 25 * cholmod_factorize simplicial or supernodal Cholesky factorization 26 * cholmod_solve solve a linear system (simplicial or supernodal) 27 * cholmod_solve2 like cholmod_solve, but reuse workspace 28 * cholmod_spsolve solve a linear system (sparse x and b) 29 * 30 * Secondary routines: 31 * ------------------ 32 * 33 * cholmod_analyze_p analyze, with user-provided permutation or f set 34 * cholmod_factorize_p factorize, with user-provided permutation or f 35 * cholmod_analyze_ordering analyze a fill-reducing ordering 36 * cholmod_etree find the elimination tree 37 * cholmod_rowcolcounts compute the row/column counts of L 38 * cholmod_amd order using AMD 39 * cholmod_colamd order using COLAMD 40 * cholmod_rowfac incremental simplicial factorization 41 * cholmod_rowfac_mask rowfac, specific to LPDASA 42 * cholmod_rowfac_mask2 rowfac, specific to LPDASA 43 * cholmod_row_subtree find the nonzero pattern of a row of L 44 * cholmod_resymbol recompute the symbolic pattern of L 45 * cholmod_resymbol_noperm recompute the symbolic pattern of L, no L->Perm 46 * cholmod_postorder postorder a tree 47 * 48 * Requires the Core module, and two packages: AMD and COLAMD. 49 * Optionally uses the Supernodal and Partition modules. 50 * Required by the Partition module. 51 */ 52 53 #ifndef CHOLMOD_CHOLESKY_H 54 #define CHOLMOD_CHOLESKY_H 55 56 #include "cholmod_config.h" 57 #include "cholmod_core.h" 58 59 #ifndef NPARTITION 60 #include "cholmod_partition.h" 61 #endif 62 63 #ifndef NSUPERNODAL 64 #include "cholmod_supernodal.h" 65 #endif 66 67 /* -------------------------------------------------------------------------- */ 68 /* cholmod_analyze: order and analyze (simplicial or supernodal) */ 69 /* -------------------------------------------------------------------------- */ 70 71 /* Orders and analyzes A, AA', PAP', or PAA'P' and returns a symbolic factor 72 * that can later be passed to cholmod_factorize. */ 73 74 cholmod_factor *cholmod_analyze 75 ( 76 /* ---- input ---- */ 77 cholmod_sparse *A, /* matrix to order and analyze */ 78 /* --------------- */ 79 cholmod_common *Common 80 ) ; 81 82 cholmod_factor *cholmod_l_analyze (cholmod_sparse *, cholmod_common *) ; 83 84 /* -------------------------------------------------------------------------- */ 85 /* cholmod_analyze_p: analyze, with user-provided permutation or f set */ 86 /* -------------------------------------------------------------------------- */ 87 88 /* Orders and analyzes A, AA', PAP', PAA'P', FF', or PFF'P and returns a 89 * symbolic factor that can later be passed to cholmod_factorize, where 90 * F = A(:,fset) if fset is not NULL and A->stype is zero. 91 * UserPerm is tried if non-NULL. */ 92 93 cholmod_factor *cholmod_analyze_p 94 ( 95 /* ---- input ---- */ 96 cholmod_sparse *A, /* matrix to order and analyze */ 97 int *UserPerm, /* user-provided permutation, size A->nrow */ 98 int *fset, /* subset of 0:(A->ncol)-1 */ 99 size_t fsize, /* size of fset */ 100 /* --------------- */ 101 cholmod_common *Common 102 ) ; 103 104 cholmod_factor *cholmod_l_analyze_p (cholmod_sparse *, SuiteSparse_long *, 105 SuiteSparse_long *, size_t, cholmod_common *) ; 106 107 /* -------------------------------------------------------------------------- */ 108 /* cholmod_analyze_p2: analyze for sparse Cholesky or sparse QR */ 109 /* -------------------------------------------------------------------------- */ 110 111 cholmod_factor *cholmod_analyze_p2 112 ( 113 /* ---- input ---- */ 114 int for_whom, /* FOR_SPQR (0): for SPQR but not GPU-accelerated 115 FOR_CHOLESKY (1): for Cholesky (GPU or not) 116 FOR_SPQRGPU (2): for SPQR with GPU acceleration */ 117 cholmod_sparse *A, /* matrix to order and analyze */ 118 int *UserPerm, /* user-provided permutation, size A->nrow */ 119 int *fset, /* subset of 0:(A->ncol)-1 */ 120 size_t fsize, /* size of fset */ 121 /* --------------- */ 122 cholmod_common *Common 123 ) ; 124 125 cholmod_factor *cholmod_l_analyze_p2 (int, cholmod_sparse *, SuiteSparse_long *, 126 SuiteSparse_long *, size_t, cholmod_common *) ; 127 128 /* -------------------------------------------------------------------------- */ 129 /* cholmod_factorize: simplicial or supernodal Cholesky factorization */ 130 /* -------------------------------------------------------------------------- */ 131 132 /* Factorizes PAP' (or PAA'P' if A->stype is 0), using a factor obtained 133 * from cholmod_analyze. The analysis can be re-used simply by calling this 134 * routine a second time with another matrix. A must have the same nonzero 135 * pattern as that passed to cholmod_analyze. */ 136 137 int cholmod_factorize 138 ( 139 /* ---- input ---- */ 140 cholmod_sparse *A, /* matrix to factorize */ 141 /* ---- in/out --- */ 142 cholmod_factor *L, /* resulting factorization */ 143 /* --------------- */ 144 cholmod_common *Common 145 ) ; 146 147 int cholmod_l_factorize (cholmod_sparse *, cholmod_factor *, cholmod_common *) ; 148 149 /* -------------------------------------------------------------------------- */ 150 /* cholmod_factorize_p: factorize, with user-provided permutation or fset */ 151 /* -------------------------------------------------------------------------- */ 152 153 /* Same as cholmod_factorize, but with more options. */ 154 155 int cholmod_factorize_p 156 ( 157 /* ---- input ---- */ 158 cholmod_sparse *A, /* matrix to factorize */ 159 double beta [2], /* factorize beta*I+A or beta*I+A'*A */ 160 int *fset, /* subset of 0:(A->ncol)-1 */ 161 size_t fsize, /* size of fset */ 162 /* ---- in/out --- */ 163 cholmod_factor *L, /* resulting factorization */ 164 /* --------------- */ 165 cholmod_common *Common 166 ) ; 167 168 int cholmod_l_factorize_p (cholmod_sparse *, double *, SuiteSparse_long *, 169 size_t, cholmod_factor *, cholmod_common *) ; 170 171 /* -------------------------------------------------------------------------- */ 172 /* cholmod_solve: solve a linear system (simplicial or supernodal) */ 173 /* -------------------------------------------------------------------------- */ 174 175 /* Solves one of many linear systems with a dense right-hand-side, using the 176 * factorization from cholmod_factorize (or as modified by any other CHOLMOD 177 * routine). D is identity for LL' factorizations. */ 178 179 #define CHOLMOD_A 0 /* solve Ax=b */ 180 #define CHOLMOD_LDLt 1 /* solve LDL'x=b */ 181 #define CHOLMOD_LD 2 /* solve LDx=b */ 182 #define CHOLMOD_DLt 3 /* solve DL'x=b */ 183 #define CHOLMOD_L 4 /* solve Lx=b */ 184 #define CHOLMOD_Lt 5 /* solve L'x=b */ 185 #define CHOLMOD_D 6 /* solve Dx=b */ 186 #define CHOLMOD_P 7 /* permute x=Px */ 187 #define CHOLMOD_Pt 8 /* permute x=P'x */ 188 189 cholmod_dense *cholmod_solve /* returns the solution X */ 190 ( 191 /* ---- input ---- */ 192 int sys, /* system to solve */ 193 cholmod_factor *L, /* factorization to use */ 194 cholmod_dense *B, /* right-hand-side */ 195 /* --------------- */ 196 cholmod_common *Common 197 ) ; 198 199 cholmod_dense *cholmod_l_solve (int, cholmod_factor *, cholmod_dense *, 200 cholmod_common *) ; 201 202 /* -------------------------------------------------------------------------- */ 203 /* cholmod_solve2: like cholmod_solve, but with reusable workspace */ 204 /* -------------------------------------------------------------------------- */ 205 206 int cholmod_solve2 /* returns TRUE on success, FALSE on failure */ 207 ( 208 /* ---- input ---- */ 209 int sys, /* system to solve */ 210 cholmod_factor *L, /* factorization to use */ 211 cholmod_dense *B, /* right-hand-side */ 212 cholmod_sparse *Bset, 213 /* ---- output --- */ 214 cholmod_dense **X_Handle, /* solution, allocated if need be */ 215 cholmod_sparse **Xset_Handle, 216 /* ---- workspace */ 217 cholmod_dense **Y_Handle, /* workspace, or NULL */ 218 cholmod_dense **E_Handle, /* workspace, or NULL */ 219 /* --------------- */ 220 cholmod_common *Common 221 ) ; 222 223 int cholmod_l_solve2 (int, cholmod_factor *, cholmod_dense *, cholmod_sparse *, 224 cholmod_dense **, cholmod_sparse **, cholmod_dense **, cholmod_dense **, 225 cholmod_common *) ; 226 227 /* -------------------------------------------------------------------------- */ 228 /* cholmod_spsolve: solve a linear system with a sparse right-hand-side */ 229 /* -------------------------------------------------------------------------- */ 230 231 cholmod_sparse *cholmod_spsolve 232 ( 233 /* ---- input ---- */ 234 int sys, /* system to solve */ 235 cholmod_factor *L, /* factorization to use */ 236 cholmod_sparse *B, /* right-hand-side */ 237 /* --------------- */ 238 cholmod_common *Common 239 ) ; 240 241 cholmod_sparse *cholmod_l_spsolve (int, cholmod_factor *, cholmod_sparse *, 242 cholmod_common *) ; 243 244 /* -------------------------------------------------------------------------- */ 245 /* cholmod_etree: find the elimination tree of A or A'*A */ 246 /* -------------------------------------------------------------------------- */ 247 248 int cholmod_etree 249 ( 250 /* ---- input ---- */ 251 cholmod_sparse *A, 252 /* ---- output --- */ 253 int *Parent, /* size ncol. Parent [j] = p if p is the parent of j */ 254 /* --------------- */ 255 cholmod_common *Common 256 ) ; 257 258 int cholmod_l_etree (cholmod_sparse *, SuiteSparse_long *, cholmod_common *) ; 259 260 /* -------------------------------------------------------------------------- */ 261 /* cholmod_rowcolcounts: compute the row/column counts of L */ 262 /* -------------------------------------------------------------------------- */ 263 264 int cholmod_rowcolcounts 265 ( 266 /* ---- input ---- */ 267 cholmod_sparse *A, /* matrix to analyze */ 268 int *fset, /* subset of 0:(A->ncol)-1 */ 269 size_t fsize, /* size of fset */ 270 int *Parent, /* size nrow. Parent [i] = p if p is the parent of i */ 271 int *Post, /* size nrow. Post [k] = i if i is the kth node in 272 * the postordered etree. */ 273 /* ---- output --- */ 274 int *RowCount, /* size nrow. RowCount [i] = # entries in the ith row of 275 * L, including the diagonal. */ 276 int *ColCount, /* size nrow. ColCount [i] = # entries in the ith 277 * column of L, including the diagonal. */ 278 int *First, /* size nrow. First [i] = k is the least postordering 279 * of any descendant of i. */ 280 int *Level, /* size nrow. Level [i] is the length of the path from 281 * i to the root, with Level [root] = 0. */ 282 /* --------------- */ 283 cholmod_common *Common 284 ) ; 285 286 int cholmod_l_rowcolcounts (cholmod_sparse *, SuiteSparse_long *, size_t, 287 SuiteSparse_long *, SuiteSparse_long *, SuiteSparse_long *, 288 SuiteSparse_long *, SuiteSparse_long *, SuiteSparse_long *, 289 cholmod_common *) ; 290 291 /* -------------------------------------------------------------------------- */ 292 /* cholmod_analyze_ordering: analyze a fill-reducing ordering */ 293 /* -------------------------------------------------------------------------- */ 294 295 int cholmod_analyze_ordering 296 ( 297 /* ---- input ---- */ 298 cholmod_sparse *A, /* matrix to analyze */ 299 int ordering, /* ordering method used */ 300 int *Perm, /* size n, fill-reducing permutation to analyze */ 301 int *fset, /* subset of 0:(A->ncol)-1 */ 302 size_t fsize, /* size of fset */ 303 /* ---- output --- */ 304 int *Parent, /* size n, elimination tree */ 305 int *Post, /* size n, postordering of elimination tree */ 306 int *ColCount, /* size n, nnz in each column of L */ 307 /* ---- workspace */ 308 int *First, /* size nworkspace for cholmod_postorder */ 309 int *Level, /* size n workspace for cholmod_postorder */ 310 /* --------------- */ 311 cholmod_common *Common 312 ) ; 313 314 int cholmod_l_analyze_ordering (cholmod_sparse *, int, SuiteSparse_long *, 315 SuiteSparse_long *, size_t, SuiteSparse_long *, SuiteSparse_long *, 316 SuiteSparse_long *, SuiteSparse_long *, SuiteSparse_long *, 317 cholmod_common *) ; 318 319 /* -------------------------------------------------------------------------- */ 320 /* cholmod_amd: order using AMD */ 321 /* -------------------------------------------------------------------------- */ 322 323 /* Finds a permutation P to reduce fill-in in the factorization of P*A*P' 324 * or P*A*A'P' */ 325 326 int cholmod_amd 327 ( 328 /* ---- input ---- */ 329 cholmod_sparse *A, /* matrix to order */ 330 int *fset, /* subset of 0:(A->ncol)-1 */ 331 size_t fsize, /* size of fset */ 332 /* ---- output --- */ 333 int *Perm, /* size A->nrow, output permutation */ 334 /* --------------- */ 335 cholmod_common *Common 336 ) ; 337 338 int cholmod_l_amd (cholmod_sparse *, SuiteSparse_long *, size_t, 339 SuiteSparse_long *, cholmod_common *) ; 340 341 /* -------------------------------------------------------------------------- */ 342 /* cholmod_colamd: order using COLAMD */ 343 /* -------------------------------------------------------------------------- */ 344 345 /* Finds a permutation P to reduce fill-in in the factorization of P*A*A'*P'. 346 * Orders F*F' where F = A (:,fset) if fset is not NULL */ 347 348 int cholmod_colamd 349 ( 350 /* ---- input ---- */ 351 cholmod_sparse *A, /* matrix to order */ 352 int *fset, /* subset of 0:(A->ncol)-1 */ 353 size_t fsize, /* size of fset */ 354 int postorder, /* if TRUE, follow with a coletree postorder */ 355 /* ---- output --- */ 356 int *Perm, /* size A->nrow, output permutation */ 357 /* --------------- */ 358 cholmod_common *Common 359 ) ; 360 361 int cholmod_l_colamd (cholmod_sparse *, SuiteSparse_long *, size_t, int, 362 SuiteSparse_long *, cholmod_common *) ; 363 364 /* -------------------------------------------------------------------------- */ 365 /* cholmod_rowfac: incremental simplicial factorization */ 366 /* -------------------------------------------------------------------------- */ 367 368 /* Partial or complete simplicial factorization. Rows and columns kstart:kend-1 369 * of L and D must be initially equal to rows/columns kstart:kend-1 of the 370 * identity matrix. Row k can only be factorized if all descendants of node 371 * k in the elimination tree have been factorized. */ 372 373 int cholmod_rowfac 374 ( 375 /* ---- input ---- */ 376 cholmod_sparse *A, /* matrix to factorize */ 377 cholmod_sparse *F, /* used for A*A' case only. F=A' or A(:,fset)' */ 378 double beta [2], /* factorize beta*I+A or beta*I+A'*A */ 379 size_t kstart, /* first row to factorize */ 380 size_t kend, /* last row to factorize is kend-1 */ 381 /* ---- in/out --- */ 382 cholmod_factor *L, 383 /* --------------- */ 384 cholmod_common *Common 385 ) ; 386 387 int cholmod_l_rowfac (cholmod_sparse *, cholmod_sparse *, double *, size_t, 388 size_t, cholmod_factor *, cholmod_common *) ; 389 390 /* -------------------------------------------------------------------------- */ 391 /* cholmod_rowfac_mask: incremental simplicial factorization */ 392 /* -------------------------------------------------------------------------- */ 393 394 /* cholmod_rowfac_mask is a version of cholmod_rowfac that is specific to 395 * LPDASA. It is unlikely to be needed by any other application. */ 396 397 int cholmod_rowfac_mask 398 ( 399 /* ---- input ---- */ 400 cholmod_sparse *A, /* matrix to factorize */ 401 cholmod_sparse *F, /* used for A*A' case only. F=A' or A(:,fset)' */ 402 double beta [2], /* factorize beta*I+A or beta*I+A'*A */ 403 size_t kstart, /* first row to factorize */ 404 size_t kend, /* last row to factorize is kend-1 */ 405 int *mask, /* if mask[i] >= 0, then set row i to zero */ 406 int *RLinkUp, /* link list of rows to compute */ 407 /* ---- in/out --- */ 408 cholmod_factor *L, 409 /* --------------- */ 410 cholmod_common *Common 411 ) ; 412 413 int cholmod_l_rowfac_mask (cholmod_sparse *, cholmod_sparse *, double *, size_t, 414 size_t, SuiteSparse_long *, SuiteSparse_long *, cholmod_factor *, 415 cholmod_common *) ; 416 417 int cholmod_rowfac_mask2 418 ( 419 /* ---- input ---- */ 420 cholmod_sparse *A, /* matrix to factorize */ 421 cholmod_sparse *F, /* used for A*A' case only. F=A' or A(:,fset)' */ 422 double beta [2], /* factorize beta*I+A or beta*I+A'*A */ 423 size_t kstart, /* first row to factorize */ 424 size_t kend, /* last row to factorize is kend-1 */ 425 int *mask, /* if mask[i] >= maskmark, then set row i to zero */ 426 int maskmark, 427 int *RLinkUp, /* link list of rows to compute */ 428 /* ---- in/out --- */ 429 cholmod_factor *L, 430 /* --------------- */ 431 cholmod_common *Common 432 ) ; 433 434 int cholmod_l_rowfac_mask2 (cholmod_sparse *, cholmod_sparse *, double *, 435 size_t, size_t, SuiteSparse_long *, SuiteSparse_long, SuiteSparse_long *, 436 cholmod_factor *, cholmod_common *) ; 437 438 /* -------------------------------------------------------------------------- */ 439 /* cholmod_row_subtree: find the nonzero pattern of a row of L */ 440 /* -------------------------------------------------------------------------- */ 441 442 /* Find the nonzero pattern of x for the system Lx=b where L = (0:k-1,0:k-1) 443 * and b = kth column of A or A*A' (rows 0 to k-1 only) */ 444 445 int cholmod_row_subtree 446 ( 447 /* ---- input ---- */ 448 cholmod_sparse *A, /* matrix to analyze */ 449 cholmod_sparse *F, /* used for A*A' case only. F=A' or A(:,fset)' */ 450 size_t k, /* row k of L */ 451 int *Parent, /* elimination tree */ 452 /* ---- output --- */ 453 cholmod_sparse *R, /* pattern of L(k,:), n-by-1 with R->nzmax >= n */ 454 /* --------------- */ 455 cholmod_common *Common 456 ) ; 457 458 int cholmod_l_row_subtree (cholmod_sparse *, cholmod_sparse *, size_t, 459 SuiteSparse_long *, cholmod_sparse *, cholmod_common *) ; 460 461 /* -------------------------------------------------------------------------- */ 462 /* cholmod_lsolve_pattern: find the nonzero pattern of x=L\b */ 463 /* -------------------------------------------------------------------------- */ 464 465 int cholmod_lsolve_pattern 466 ( 467 /* ---- input ---- */ 468 cholmod_sparse *B, /* sparse right-hand-side (a single sparse column) */ 469 cholmod_factor *L, /* the factor L from which parent(i) is derived */ 470 /* ---- output --- */ 471 cholmod_sparse *X, /* pattern of X=L\B, n-by-1 with X->nzmax >= n */ 472 /* --------------- */ 473 cholmod_common *Common 474 ) ; 475 476 int cholmod_l_lsolve_pattern (cholmod_sparse *, cholmod_factor *, 477 cholmod_sparse *, cholmod_common *) ; 478 479 /* -------------------------------------------------------------------------- */ 480 /* cholmod_row_lsubtree: find the nonzero pattern of a row of L */ 481 /* -------------------------------------------------------------------------- */ 482 483 /* Identical to cholmod_row_subtree, except that it finds the elimination tree 484 * from L itself. */ 485 486 int cholmod_row_lsubtree 487 ( 488 /* ---- input ---- */ 489 cholmod_sparse *A, /* matrix to analyze */ 490 int *Fi, size_t fnz, /* nonzero pattern of kth row of A', not required 491 * for the symmetric case. Need not be sorted. */ 492 size_t k, /* row k of L */ 493 cholmod_factor *L, /* the factor L from which parent(i) is derived */ 494 /* ---- output --- */ 495 cholmod_sparse *R, /* pattern of L(k,:), n-by-1 with R->nzmax >= n */ 496 /* --------------- */ 497 cholmod_common *Common 498 ) ; 499 500 int cholmod_l_row_lsubtree (cholmod_sparse *, SuiteSparse_long *, size_t, 501 size_t, cholmod_factor *, cholmod_sparse *, cholmod_common *) ; 502 503 /* -------------------------------------------------------------------------- */ 504 /* cholmod_resymbol: recompute the symbolic pattern of L */ 505 /* -------------------------------------------------------------------------- */ 506 507 /* Remove entries from L that are not in the factorization of P*A*P', P*A*A'*P', 508 * or P*F*F'*P' (depending on A->stype and whether fset is NULL or not). 509 * 510 * cholmod_resymbol is the same as cholmod_resymbol_noperm, except that it 511 * first permutes A according to L->Perm. A can be upper/lower/unsymmetric, 512 * in contrast to cholmod_resymbol_noperm (which can be lower or unsym). */ 513 514 int cholmod_resymbol 515 ( 516 /* ---- input ---- */ 517 cholmod_sparse *A, /* matrix to analyze */ 518 int *fset, /* subset of 0:(A->ncol)-1 */ 519 size_t fsize, /* size of fset */ 520 int pack, /* if TRUE, pack the columns of L */ 521 /* ---- in/out --- */ 522 cholmod_factor *L, /* factorization, entries pruned on output */ 523 /* --------------- */ 524 cholmod_common *Common 525 ) ; 526 527 int cholmod_l_resymbol (cholmod_sparse *, SuiteSparse_long *, size_t, int, 528 cholmod_factor *, cholmod_common *) ; 529 530 /* -------------------------------------------------------------------------- */ 531 /* cholmod_resymbol_noperm: recompute the symbolic pattern of L, no L->Perm */ 532 /* -------------------------------------------------------------------------- */ 533 534 /* Remove entries from L that are not in the factorization of A, A*A', 535 * or F*F' (depending on A->stype and whether fset is NULL or not). */ 536 537 int cholmod_resymbol_noperm 538 ( 539 /* ---- input ---- */ 540 cholmod_sparse *A, /* matrix to analyze */ 541 int *fset, /* subset of 0:(A->ncol)-1 */ 542 size_t fsize, /* size of fset */ 543 int pack, /* if TRUE, pack the columns of L */ 544 /* ---- in/out --- */ 545 cholmod_factor *L, /* factorization, entries pruned on output */ 546 /* --------------- */ 547 cholmod_common *Common 548 ) ; 549 550 int cholmod_l_resymbol_noperm (cholmod_sparse *, SuiteSparse_long *, size_t, int, 551 cholmod_factor *, cholmod_common *) ; 552 553 /* -------------------------------------------------------------------------- */ 554 /* cholmod_rcond: compute rough estimate of reciprocal of condition number */ 555 /* -------------------------------------------------------------------------- */ 556 557 double cholmod_rcond /* return min(diag(L)) / max(diag(L)) */ 558 ( 559 /* ---- input ---- */ 560 cholmod_factor *L, 561 /* --------------- */ 562 cholmod_common *Common 563 ) ; 564 565 double cholmod_l_rcond (cholmod_factor *, cholmod_common *) ; 566 567 /* -------------------------------------------------------------------------- */ 568 /* cholmod_postorder: Compute the postorder of a tree */ 569 /* -------------------------------------------------------------------------- */ 570 571 SuiteSparse_long cholmod_postorder /* return # of nodes postordered */ 572 ( 573 /* ---- input ---- */ 574 int *Parent, /* size n. Parent [j] = p if p is the parent of j */ 575 size_t n, 576 int *Weight_p, /* size n, optional. Weight [j] is weight of node j */ 577 /* ---- output --- */ 578 int *Post, /* size n. Post [k] = j is kth in postordered tree */ 579 /* --------------- */ 580 cholmod_common *Common 581 ) ; 582 583 SuiteSparse_long cholmod_l_postorder (SuiteSparse_long *, size_t, 584 SuiteSparse_long *, SuiteSparse_long *, cholmod_common *) ; 585 586 #endif 587