1 // ============================================================================= 2 // === SuiteSparseQR.hpp ======================================================= 3 // ============================================================================= 4 5 // User include file for C++ programs. 6 7 #ifndef SUITESPARSEQR_H 8 #define SUITESPARSEQR_H 9 10 // ----------------------------------------------------------------------------- 11 // include files 12 // ----------------------------------------------------------------------------- 13 14 #ifdef GPU_BLAS 15 #include <cublas_v2.h> 16 #endif 17 #define SUITESPARSE_GPU_EXTERN_ON 18 extern "C" 19 { 20 #include "SuiteSparseQR_definitions.h" 21 #include "cholmod.h" 22 } 23 #undef SUITESPARSE_GPU_EXTERN_ON 24 25 // ============================================================================= 26 // === spqr_gpu ================================================================ 27 // ============================================================================= 28 29 struct spqr_gpu 30 { 31 SuiteSparse_long *RimapOffsets; // Stores front offsets into Rimap 32 SuiteSparse_long RimapSize; // Allocated space for in Rimap 33 34 SuiteSparse_long *RjmapOffsets; // Stores front offsets into Rjmap 35 SuiteSparse_long RjmapSize; // Allocated space for Rjmap 36 37 SuiteSparse_long numStages; // # of Stages required to factorize 38 SuiteSparse_long *Stagingp; // Pointers into Post for boundaries 39 SuiteSparse_long *StageMap; // Mapping of front to stage # 40 size_t *FSize; // Total size of fronts in a stage 41 size_t *RSize; // Total size of R+C for a stage 42 size_t *SSize; // Total size of S for a stage 43 SuiteSparse_long *FOffsets; // F Offsets relative to a base 44 SuiteSparse_long *ROffsets; // R Offsets relative to a base 45 SuiteSparse_long *SOffsets; // S Offsets relative to a base 46 }; 47 48 // ============================================================================= 49 // === spqr_symbolic =========================================================== 50 // ============================================================================= 51 52 // The contents of this object do not change during numeric factorization. The 53 // Symbolic object depends only on the pattern of the input matrix, and not its 54 // values. These contents also do not change with column pivoting for rank 55 // detection. This makes parallelism easier to manage, since all threads can 56 // have access to this object without synchronization. 57 // 58 // The total size of the Symbolic object is (10 + 2*m + anz + 2*n + 5*nf + rnz) 59 // Long's, where the user's input A matrix is m-by-n with anz nonzeros, nf <= 60 // MIN(m,n) is the number of frontal matrices, and rnz <= nnz(R) is the number 61 // of column indices used to represent the supernodal form of R (one Long per 62 // non-pivotal column index in the leading row of each block of R). 63 64 struct spqr_symbolic 65 { 66 67 // ------------------------------------------------------------------------- 68 // row-form of the input matrix and its permutations 69 // ------------------------------------------------------------------------- 70 71 // During symbolic analysis, the nonzero pattern of S = A(P,Q) is 72 // constructed, where A is the user's input matrix. Its numerical values 73 // are also constructed, but they do not become part of the Symbolic 74 // object. The matrix S is stored in row-oriented form. The rows of S are 75 // sorted according to their leftmost column index (via PLinv). Column 76 // indices in each row of S are in strictly ascending order, even though 77 // the input matrix A need not be sorted. 78 79 SuiteSparse_long m, n, anz ; // S is m-by-n with anz entries 80 81 SuiteSparse_long *Sp ; // size m+1, row pointers of S 82 83 SuiteSparse_long *Sj ; // size anz = Sp [n], column indices of S 84 85 SuiteSparse_long *Qfill ; // size n, fill-reducing column permutation. 86 // Qfill [k] = j if column k of A is column j of S. 87 88 SuiteSparse_long *PLinv ; // size m, inverse row permutation that places 89 // S=A(P,Q) in increasing order of leftmost column 90 // index. PLinv [i] = k if row i of A is row k of S. 91 92 SuiteSparse_long *Sleft ; // size n+2. The list of rows of S whose 93 // leftmost column index is j is given by 94 // Sleft [j] ... Sleft [j+1]-1. This can be empty (that is, Sleft 95 // [j] can equal Sleft [j+1]). Sleft [n] is the number of 96 // non-empty rows of S, and Sleft [n+1] == m. That is, Sleft [n] 97 // ... Sleft [n+1]-1 gives the empty rows of S, if any. 98 99 // ------------------------------------------------------------------------- 100 // frontal matrices: pattern and tree 101 // ------------------------------------------------------------------------- 102 103 // Each frontal matrix is fm-by-fn, with fnpiv pivot columns. The fn 104 // column indices are given by a set of size fnpiv pivot columns, defined 105 // by Super, followed by the pattern Rj [ Rp[f] ... Rp[f+1]-1 ]. 106 107 // The row indices of the front are not kept. If the Householder vectors 108 // are not kept, the row indices are not needed. If the Householder 109 // vectors are kept, the row indices are computed dynamically during 110 // numerical factorization. 111 112 SuiteSparse_long nf ; // number of frontal matrices; nf <= MIN (m,n) 113 SuiteSparse_long maxfn ; // max # of columns in any front 114 115 // parent, child, and childp define the row merge tree or etree (A'A) 116 SuiteSparse_long *Parent ; // size nf+1 117 SuiteSparse_long *Child ; // size nf+1 118 SuiteSparse_long *Childp ; // size nf+2 119 120 // The parent of a front f is Parent [f], or EMPTY if f=nf. 121 // A list of children of f can be obtained in the list 122 // Child [Childp [f] ... Childp [f+1]-1]. 123 124 // Node nf in the tree is a placeholder; it does not represent a frontal 125 // matrix. All roots of the frontal "tree" (may be a forest) have the 126 // placeholder node nf as their parent. Thus, the tree of nodes 0:nf is 127 // truly a tree, with just one parent (node nf). 128 129 SuiteSparse_long *Super ; // size nf+1. Super [f] gives the first 130 // pivot column in the front F. This refers to a column of S. The 131 // number of expected pivot columns in F is thus 132 // Super [f+1] - Super [f]. 133 134 SuiteSparse_long *Rp ; // size nf+1 135 SuiteSparse_long *Rj ; // size rjsize; compressed supernodal form of R 136 137 SuiteSparse_long *Post ; // size nf+1, post ordering of frontal tree. 138 // f=Post[k] gives the kth node in the postordered tree 139 140 SuiteSparse_long rjsize ; // size of Rj 141 142 SuiteSparse_long do_rank_detection ; // TRUE: allow for tol >= 0. 143 // FALSE: ignore tol 144 145 // the rest depends on whether or not rank-detection is allowed: 146 SuiteSparse_long maxstack ; // max stack size (sequential case) 147 SuiteSparse_long hisize ; // size of Hii 148 149 SuiteSparse_long keepH ; // TRUE if H is present 150 151 SuiteSparse_long *Hip ; // size nf+1. If H is kept, the row indices 152 // of frontal matrix f are in Hii [Hip [f] ... Hip [f] + Hm [f]], 153 // where Hii and Hm are stored in the numeric object. 154 155 // There is one block row of R per frontal matrix. 156 // The fn column indices of R are given by Rj [Rp [f] ... Rp [f+1]-1], 157 // where the first fp column indices are Super [f] ... Super [f+1]-1. 158 // The remaining column indices in Rj [...] are non-pivotal, and are 159 // in the range Super [f+1] to n. The number of rows of R is at 160 // most fp, but can be less if dead columns appear in the matrix. 161 // The number of columns in the contribution block C is always 162 // cn = fn - fp, where fn = Rp [f+1] - Rp [f]. 163 164 SuiteSparse_long ntasks ; // number of tasks in task graph 165 SuiteSparse_long ns ; // number of stacks 166 167 // ------------------------------------------------------------------------- 168 // the rest of the QR symbolic object is present only if ntasks > 1 169 // ------------------------------------------------------------------------- 170 171 // Task tree (nodes 0:ntasks), including placeholder node 172 SuiteSparse_long *TaskChildp ; // size ntasks+2 173 SuiteSparse_long *TaskChild ; // size ntasks+1 174 175 SuiteSparse_long *TaskStack ; // size ntasks+1 176 177 // list of fronts for each task 178 SuiteSparse_long *TaskFront ; // size nf+1 179 SuiteSparse_long *TaskFrontp ; // size ntasks+2 180 181 SuiteSparse_long *On_stack ; // size nf+1, front f is on 182 // stack On_stack [f] 183 184 // size of each stack 185 SuiteSparse_long *Stack_maxstack ; // size ns+2 186 187 // number of rows for each front 188 SuiteSparse_long *Fm ; // size nf+1 189 190 // number of rows in the contribution block of each front 191 SuiteSparse_long *Cm ; // size nf+1 192 193 // from CHOLMOD's supernodal analysis, needed for GPU factorization 194 size_t maxcsize ; 195 size_t maxesize ; 196 SuiteSparse_long *ColCount ; 197 // SuiteSparse_long *px ; 198 199 // ------------------------------------------------------------------------- 200 // GPU structure 201 // ------------------------------------------------------------------------- 202 203 // This is NULL if the GPU is not in use. The GPU must be enabled at 204 // compile time (-DGPU_BLAS enables the GPU). If the Householder vectors 205 // are requested, if TBB is used (Common->SPQR_grain > 1), or if rank 206 // detection is requested, then the GPU is disabled. 207 208 spqr_gpu *QRgpu ; 209 210 } ; 211 212 213 // ============================================================================= 214 // === spqr_numeric ============================================================ 215 // ============================================================================= 216 217 // The Numeric object contains the numerical values of the triangular/ 218 // trapezoidal factor R, and optionally the Householder vectors H if they 219 // are kept. 220 221 template <typename Entry> struct spqr_numeric 222 { 223 224 // ------------------------------------------------------------------------- 225 // Numeric R factor 226 // ------------------------------------------------------------------------- 227 228 Entry **Rblock ; // size nf. R [f] is an (Entry *) pointer to the 229 // R block for front F. It is an upper trapezoidal 230 // of size Rm(f)-by-Rn(f), but only the upper 231 // triangular part is stored in column-packed format. 232 233 Entry **Stacks ; // size ns; an array of stacks holding the R and H 234 // factors and the current frontal matrix F at the head. 235 // This is followed by empty space, then the C blocks of 236 // prior frontal matrices at the bottom. When the 237 // factorization is complete, only the R and H part at 238 // the head of each stack is left. 239 240 SuiteSparse_long *Stack_size ; // size ns; Stack_size [s] is the size 241 // of Stacks [s] 242 243 SuiteSparse_long hisize ; // size of Hii 244 245 SuiteSparse_long n ; // A is m-by-n 246 SuiteSparse_long m ; 247 SuiteSparse_long nf ; // number of frontal matrices 248 SuiteSparse_long ntasks ; // # of tasks in task graph actually used 249 SuiteSparse_long ns ; // number of stacks actually used 250 SuiteSparse_long maxstack ; // size of sequential stack, if used 251 252 // ------------------------------------------------------------------------- 253 // for rank detection and m < n case 254 // ------------------------------------------------------------------------- 255 256 char *Rdead ; // size n, Rdead [k] = 1 if k is a dead pivot column, 257 // Rdead [k] = 0 otherwise. If no columns are dead, 258 // this is NULL. If m < n, then at least m-n columns 259 // will be dead. 260 261 SuiteSparse_long rank ; // number of live pivot columns 262 SuiteSparse_long rank1 ; // number of live pivot columns in first ntol 263 // columns of A 264 265 SuiteSparse_long maxfrank ; // max number of rows in any R block 266 267 double norm_E_fro ; // 2-norm of w, the vector of dead column 2-norms 268 269 // ------------------------------------------------------------------------- 270 // for keeping Householder vectors 271 // ------------------------------------------------------------------------- 272 273 // The factorization is R = (H_s * ... * H_2 * H_1) * P_H 274 // where P_H is the permutation HPinv, and H_1, ... H_s are the Householder 275 // vectors (s = rjsize). 276 277 SuiteSparse_long keepH ; // TRUE if H is present 278 279 SuiteSparse_long rjsize ; // size of Hstair and HTau 280 281 SuiteSparse_long *HStair ; // size rjsize. The list Hstair [Rp [f] ... 282 // Rp [f+1]-1] gives the staircase for front F 283 284 Entry *HTau ; // size rjsize. The list HTau [Rp [f] ... Rp [f+1]-1] 285 // gives the Householder coefficients for front F 286 287 SuiteSparse_long *Hii ; // size hisize, row indices of H. 288 289 SuiteSparse_long *HPinv ; // size m. HPinv [i] = k if row i of A and H 290 // is row k of R. This permutation includes 291 // QRsym->PLinv, and the permutation constructed via 292 // pivotal row ordering during factorization. 293 294 SuiteSparse_long *Hm ; // size nf, Hm [f] = # of rows in front F 295 SuiteSparse_long *Hr ; // size nf, Hr [f] = # of rows in R block of 296 // front F 297 SuiteSparse_long maxfm ; // max (Hm [0:nf-1]), computed only if H kept 298 299 } ; 300 301 302 // ============================================================================= 303 // === SuiteSparseQR_factorization ============================================= 304 // ============================================================================= 305 306 // A combined symbolic+numeric QR factorization of A or [A B], 307 // with singletons 308 309 template <typename Entry> struct SuiteSparseQR_factorization 310 { 311 312 // QR factorization of A or [A Binput] after singletons have been removed 313 double tol ; // tol used 314 spqr_symbolic *QRsym ; 315 spqr_numeric <Entry> *QRnum ; 316 317 // singletons, in compressed-row form; R is n1rows-by-n 318 SuiteSparse_long *R1p ; // size n1rows+1 319 SuiteSparse_long *R1j ; 320 Entry *R1x ; 321 SuiteSparse_long r1nz ; // nnz (R1) 322 323 // combined singleton and fill-reducing permutation 324 SuiteSparse_long *Q1fill ; 325 SuiteSparse_long *P1inv ; 326 SuiteSparse_long *HP1inv ; // NULL if n1cols == 0, in which case 327 // QRnum->HPinv serves in its place. 328 329 // Rmap and RmapInv are NULL if QR->rank == A->ncol 330 SuiteSparse_long *Rmap ; // size n. Rmap [j] = k if column j of R is 331 // the kth live column and where k < QR->rank; 332 // otherwise, if j is a dead column, then 333 // k >= QR->rank. 334 335 SuiteSparse_long *RmapInv ; 336 337 SuiteSparse_long n1rows ; // number of singleton rows of [A B] 338 SuiteSparse_long n1cols ; // number of singleton columns of [A B] 339 340 SuiteSparse_long narows ; // number of rows of A 341 SuiteSparse_long nacols ; // number of columns of A 342 SuiteSparse_long bncols ; // number of columns of B 343 SuiteSparse_long rank ; // rank estimate of A (n1rows + QRnum->rank1), 344 // ranges from 0 to min(m,n) 345 346 int allow_tol ; // if TRUE, do rank detection 347 } ; 348 349 350 // ============================================================================= 351 // === Simple user-callable SuiteSparseQR functions ============================ 352 // ============================================================================= 353 354 // SuiteSparseQR Sparse QR factorization and solve 355 // SuiteSparseQR_qmult Q'*X, Q*X, X*Q', or X*Q for X full or sparse 356 357 // returns rank(A) estimate, or EMPTY on failure 358 template <typename Entry> SuiteSparse_long SuiteSparseQR 359 ( 360 // inputs, not modified 361 int ordering, // all, except 3:given treated as 0:fixed 362 double tol, // only accept singletons above tol 363 364 SuiteSparse_long econ, // number of rows of C and R to return; a value 365 // less than the rank r of A is treated as r, and 366 // a value greater than m is treated as m. 367 368 int getCTX, // if 0: return Z = C of size econ-by-bncols 369 // if 1: return Z = C' of size bncols-by-econ 370 // if 2: return Z = X of size econ-by-bncols 371 372 cholmod_sparse *A, // m-by-n sparse matrix 373 374 // B is either sparse or dense. If Bsparse is non-NULL, B is sparse and 375 // Bdense is ignored. If Bsparse is NULL and Bdense is non-NULL, then B is 376 // dense. B is not present if both are NULL. 377 cholmod_sparse *Bsparse, 378 cholmod_dense *Bdense, 379 380 // output arrays, neither allocated nor defined on input. 381 382 // Z is the matrix C, C', or X 383 cholmod_sparse **Zsparse, 384 cholmod_dense **Zdense, 385 cholmod_sparse **R, // the R factor 386 SuiteSparse_long **E, // size n; fill-reducing ordering of A. 387 cholmod_sparse **H, // the Householder vectors (m-by-nh) 388 SuiteSparse_long **HPinv,// size m; row permutation for H 389 cholmod_dense **HTau, // size nh, Householder coefficients 390 391 // workspace and parameters 392 cholmod_common *cc 393 ) ; 394 395 // X = A\dense(B) 396 template <typename Entry> cholmod_dense *SuiteSparseQR 397 ( 398 int ordering, // all, except 3:given treated as 0:fixed 399 double tol, 400 cholmod_sparse *A, // m-by-n sparse matrix 401 cholmod_dense *B, // m-by-nrhs 402 cholmod_common *cc // workspace and parameters 403 ) ; 404 405 // X = A\dense(B) using default ordering and tolerance 406 template <typename Entry> cholmod_dense *SuiteSparseQR 407 ( 408 cholmod_sparse *A, // m-by-n sparse matrix 409 cholmod_dense *B, // m-by-nrhs 410 cholmod_common *cc // workspace and parameters 411 ) ; 412 413 // X = A\sparse(B) 414 template <typename Entry> cholmod_sparse *SuiteSparseQR 415 ( 416 int ordering, // all, except 3:given treated as 0:fixed 417 double tol, 418 cholmod_sparse *A, // m-by-n sparse matrix 419 cholmod_sparse *B, // m-by-nrhs 420 cholmod_common *cc // workspace and parameters 421 ) ; 422 423 // [Q,R,E] = qr(A), returning Q as a sparse matrix 424 template <typename Entry> SuiteSparse_long SuiteSparseQR 425 // returns rank(A) estimate 426 ( 427 int ordering, // all, except 3:given treated as 0:fixed 428 double tol, 429 SuiteSparse_long econ, 430 cholmod_sparse *A, // m-by-n sparse matrix 431 // outputs 432 cholmod_sparse **Q, // m-by-e sparse matrix where e=max(econ,rank(A)) 433 cholmod_sparse **R, // e-by-n sparse matrix 434 SuiteSparse_long **E, // permutation of 0:n-1, NULL if identity 435 cholmod_common *cc // workspace and parameters 436 ) ; 437 438 // [Q,R,E] = qr(A), discarding Q 439 template <typename Entry> SuiteSparse_long SuiteSparseQR 440 // returns rank(A) estimate 441 ( 442 int ordering, // all, except 3:given treated as 0:fixed 443 double tol, 444 SuiteSparse_long econ, 445 cholmod_sparse *A, // m-by-n sparse matrix 446 // outputs 447 cholmod_sparse **R, // e-by-n sparse matrix 448 SuiteSparse_long **E, // permutation of 0:n-1, NULL if identity 449 cholmod_common *cc // workspace and parameters 450 ) ; 451 452 // [C,R,E] = qr(A,B), where C and B are dense 453 template <typename Entry> SuiteSparse_long SuiteSparseQR 454 ( 455 // inputs, not modified 456 int ordering, // all, except 3:given treated as 0:fixed 457 double tol, // only accept singletons above tol 458 SuiteSparse_long econ, // number of rows of C and R to return 459 cholmod_sparse *A, // m-by-n sparse matrix 460 cholmod_dense *B, // m-by-nrhs dense matrix 461 // outputs 462 cholmod_dense **C, // C = Q'*B, an e-by-nrhs dense matrix 463 cholmod_sparse **R, // e-by-n sparse matrix where e=max(econ,rank(A)) 464 SuiteSparse_long **E, // permutation of 0:n-1, NULL if identity 465 cholmod_common *cc // workspace and parameters 466 ) ; 467 468 // [C,R,E] = qr(A,B), where C and B are sparse 469 template <typename Entry> SuiteSparse_long SuiteSparseQR 470 ( 471 // inputs, not modified 472 int ordering, // all, except 3:given treated as 0:fixed 473 double tol, // only accept singletons above tol 474 SuiteSparse_long econ, // number of rows of C and R to return 475 cholmod_sparse *A, // m-by-n sparse matrix 476 cholmod_sparse *B, // m-by-nrhs sparse matrix 477 // outputs 478 cholmod_sparse **C, // C = Q'*B, an e-by-nrhs sparse matrix 479 cholmod_sparse **R, // e-by-n sparse matrix where e=max(econ,rank(A)) 480 SuiteSparse_long **E, // permutation of 0:n-1, NULL if identity 481 cholmod_common *cc // workspace and parameters 482 ) ; 483 484 // [Q,R,E] = qr(A) where Q is returned in Householder form 485 template <typename Entry> SuiteSparse_long SuiteSparseQR 486 ( 487 // inputs, not modified 488 int ordering, // all, except 3:given treated as 0:fixed 489 double tol, // only accept singletons above tol 490 SuiteSparse_long econ, // number of rows of C and R to return 491 cholmod_sparse *A, // m-by-n sparse matrix 492 // outputs 493 cholmod_sparse **R, // the R factor 494 SuiteSparse_long **E, // permutation of 0:n-1, NULL if identity 495 cholmod_sparse **H, // the Householder vectors (m-by-nh) 496 SuiteSparse_long **HPinv,// size m; row permutation for H 497 cholmod_dense **HTau, // size nh, Householder coefficients 498 cholmod_common *cc // workspace and parameters 499 ) ; 500 501 // ============================================================================= 502 // === SuiteSparseQR_qmult ===================================================== 503 // ============================================================================= 504 505 // This function takes as input the matrix Q in Householder form, as returned 506 // by SuiteSparseQR (... H, HPinv, HTau, cc) above. 507 508 // returns Y of size m-by-n (NULL on failure) 509 template <typename Entry> cholmod_dense *SuiteSparseQR_qmult 510 ( 511 // inputs, no modified 512 int method, // 0,1,2,3 513 cholmod_sparse *H, // either m-by-nh or n-by-nh 514 cholmod_dense *HTau, // size 1-by-nh 515 SuiteSparse_long *HPinv,// size mh 516 cholmod_dense *Xdense, // size m-by-n 517 518 // workspace and parameters 519 cholmod_common *cc 520 ) ; 521 522 template <typename Entry> cholmod_sparse *SuiteSparseQR_qmult 523 ( 524 // inputs, no modified 525 int method, // 0,1,2,3 526 cholmod_sparse *H, // either m-by-nh or n-by-nh 527 cholmod_dense *HTau, // size 1-by-nh 528 SuiteSparse_long *HPinv,// size mh 529 cholmod_sparse *X, 530 531 // workspace and parameters 532 cholmod_common *cc 533 ) ; 534 535 // ============================================================================= 536 // === Expert user-callable SuiteSparseQR functions ============================ 537 // ============================================================================= 538 539 #ifndef NEXPERT 540 541 // These functions are "expert" routines, allowing reuse of the QR 542 // factorization for different right-hand-sides. They also allow the user to 543 // find the minimum 2-norm solution to an undertermined system of equations. 544 545 template <typename Entry> 546 SuiteSparseQR_factorization <Entry> *SuiteSparseQR_factorize 547 ( 548 // inputs, not modified: 549 int ordering, // all, except 3:given treated as 0:fixed 550 double tol, // treat columns with 2-norm <= tol as zero 551 cholmod_sparse *A, // sparse matrix to factorize 552 // workspace and parameters 553 cholmod_common *cc 554 ) ; 555 556 template <typename Entry> cholmod_dense *SuiteSparseQR_solve // returns X 557 ( 558 // inputs, not modified: 559 int system, // which system to solve 560 SuiteSparseQR_factorization <Entry> *QR, // of an m-by-n sparse matrix A 561 cholmod_dense *B, // right-hand-side, m-by-nrhs or n-by-nrhs 562 // workspace and parameters 563 cholmod_common *cc 564 ) ; 565 566 template <typename Entry> cholmod_sparse *SuiteSparseQR_solve // returns X 567 ( 568 // inputs, not modified: 569 int system, // which system to solve (0,1,2,3) 570 SuiteSparseQR_factorization <Entry> *QR, // of an m-by-n sparse matrix A 571 cholmod_sparse *Bsparse, // right-hand-side, m-by-nrhs or n-by-nrhs 572 // workspace and parameters 573 cholmod_common *cc 574 ) ; 575 576 // returns Y of size m-by-n, or NULL on failure 577 template <typename Entry> cholmod_dense *SuiteSparseQR_qmult 578 ( 579 // inputs, not modified 580 int method, // 0,1,2,3 (same as SuiteSparseQR_qmult) 581 SuiteSparseQR_factorization <Entry> *QR, // of an m-by-n sparse matrix A 582 cholmod_dense *Xdense, // size m-by-n with leading dimension ldx 583 // workspace and parameters 584 cholmod_common *cc 585 ) ; 586 587 // returns Y of size m-by-n, or NULL on failure 588 template <typename Entry> cholmod_sparse *SuiteSparseQR_qmult 589 ( 590 // inputs, not modified 591 int method, // 0,1,2,3 592 SuiteSparseQR_factorization <Entry> *QR, // of an m-by-n sparse matrix A 593 cholmod_sparse *Xsparse, // size m-by-n 594 // workspace and parameters 595 cholmod_common *cc 596 ) ; 597 598 // free the QR object 599 template <typename Entry> int SuiteSparseQR_free 600 ( 601 SuiteSparseQR_factorization <Entry> **QR, // of an m-by-n sparse matrix A 602 cholmod_common *cc 603 ) ; 604 605 // find the min 2-norm solution to a sparse linear system 606 template <typename Entry> cholmod_dense *SuiteSparseQR_min2norm 607 ( 608 int ordering, // all, except 3:given treated as 0:fixed 609 double tol, 610 cholmod_sparse *A, 611 cholmod_dense *B, 612 cholmod_common *cc 613 ) ; 614 615 template <typename Entry> cholmod_sparse *SuiteSparseQR_min2norm 616 ( 617 int ordering, // all, except 3:given treated as 0:fixed 618 double tol, 619 cholmod_sparse *A, 620 cholmod_sparse *B, 621 cholmod_common *cc 622 ) ; 623 624 // symbolic QR factorization; no singletons exploited 625 template <typename Entry> 626 SuiteSparseQR_factorization <Entry> *SuiteSparseQR_symbolic 627 ( 628 // inputs: 629 int ordering, // all, except 3:given treated as 0:fixed 630 int allow_tol, // if FALSE, tol is ignored by the numeric 631 // factorization, and no rank detection is performed 632 cholmod_sparse *A, // sparse matrix to factorize (A->x ignored) 633 cholmod_common *cc // workspace and parameters 634 ) ; 635 636 // numeric QR factorization; 637 template <typename Entry> int SuiteSparseQR_numeric 638 ( 639 // inputs: 640 double tol, // treat columns with 2-norm <= tol as zero 641 cholmod_sparse *A, // sparse matrix to factorize 642 // input/output 643 SuiteSparseQR_factorization <Entry> *QR, 644 cholmod_common *cc // workspace and parameters 645 ) ; 646 647 #endif 648 649 #endif 650