1 /* ========================================================================== */ 2 /* === Include/cholmod_core.h =============================================== */ 3 /* ========================================================================== */ 4 5 /* ----------------------------------------------------------------------------- 6 * CHOLMOD/Include/cholmod_core.h. 7 * Copyright (C) 2005-2019, Univ. of Florida. Author: Timothy A. Davis 8 * -------------------------------------------------------------------------- */ 9 10 /* CHOLMOD Core module: basic CHOLMOD objects and routines. 11 * Required by all CHOLMOD modules. Requires no other module or package. 12 * 13 * The CHOLMOD modules are: 14 * 15 * Core basic data structures and definitions 16 * Check check/print the 5 CHOLMOD objects, & 3 types of integer vectors 17 * Cholesky sparse Cholesky factorization 18 * Modify sparse Cholesky update/downdate/row-add/row-delete 19 * MatrixOps sparse matrix functions (add, multiply, norm, ...) 20 * Supernodal supernodal sparse Cholesky factorization 21 * Partition graph-partitioning based orderings 22 * 23 * The CHOLMOD objects: 24 * -------------------- 25 * 26 * cholmod_common parameters, statistics, and workspace 27 * cholmod_sparse a sparse matrix in compressed column form 28 * cholmod_factor an LL' or LDL' factorization 29 * cholmod_dense a dense matrix 30 * cholmod_triplet a sparse matrix in "triplet" form 31 * 32 * The Core module described here defines the CHOLMOD data structures, and 33 * basic operations on them. To create and solve a sparse linear system Ax=b, 34 * the user must create A and b, populate them with values, and then pass them 35 * to the routines in the CHOLMOD Cholesky module. There are two primary 36 * methods for creating A: (1) allocate space for a column-oriented sparse 37 * matrix and fill it with pattern and values, or (2) create a triplet form 38 * matrix and convert it to a sparse matrix. The latter option is simpler. 39 * 40 * The matrices b and x are typically dense matrices, but can also be sparse. 41 * You can allocate and free them as dense matrices with the 42 * cholmod_allocate_dense and cholmod_free_dense routines. 43 * 44 * The cholmod_factor object contains the symbolic and numeric LL' or LDL' 45 * factorization of sparse symmetric matrix. The matrix must be positive 46 * definite for an LL' factorization. It need only be symmetric and have well- 47 * conditioned leading submatrices for it to have an LDL' factorization 48 * (CHOLMOD does not pivot for numerical stability). It is typically created 49 * with the cholmod_factorize routine in the Cholesky module, but can also 50 * be initialized to L=D=I in the Core module and then modified by the Modify 51 * module. It must be freed with cholmod_free_factor, defined below. 52 * 53 * The Core routines for each object are described below. Each list is split 54 * into two parts: the primary routines and secondary routines. 55 * 56 * ============================================================================ 57 * === cholmod_common ========================================================= 58 * ============================================================================ 59 * 60 * The Common object contains control parameters, statistics, and 61 * You must call cholmod_start before calling any other CHOLMOD routine, and 62 * must call cholmod_finish as your last call to CHOLMOD, with two exceptions: 63 * you may call cholmod_print_common and cholmod_check_common in the Check 64 * module after calling cholmod_finish. 65 * 66 * cholmod_start first call to CHOLMOD 67 * cholmod_finish last call to CHOLMOD 68 * ----------------------------- 69 * cholmod_defaults restore default parameters 70 * cholmod_maxrank maximum rank for update/downdate 71 * cholmod_allocate_work allocate workspace in Common 72 * cholmod_free_work free workspace in Common 73 * cholmod_clear_flag clear Flag workspace in Common 74 * cholmod_error called when CHOLMOD encounters an error 75 * cholmod_dbound for internal use in CHOLMOD only 76 * cholmod_hypot compute sqrt (x*x + y*y) accurately 77 * cholmod_divcomplex complex division, c = a/b 78 * 79 * ============================================================================ 80 * === cholmod_sparse ========================================================= 81 * ============================================================================ 82 * 83 * A sparse matrix is held in compressed column form. In the basic type 84 * ("packed", which corresponds to a MATLAB sparse matrix), an n-by-n matrix 85 * with nz entries is held in three arrays: p of size n+1, i of size nz, and x 86 * of size nz. Row indices of column j are held in i [p [j] ... p [j+1]-1] and 87 * in the same locations in x. There may be no duplicate entries in a column. 88 * Row indices in each column may be sorted or unsorted (CHOLMOD keeps track). 89 * A->stype determines the storage mode: 0 if both upper/lower parts are stored, 90 * -1 if A is symmetric and just tril(A) is stored, +1 if symmetric and triu(A) 91 * is stored. 92 * 93 * cholmod_allocate_sparse allocate a sparse matrix 94 * cholmod_free_sparse free a sparse matrix 95 * ----------------------------- 96 * cholmod_reallocate_sparse change the size (# entries) of sparse matrix 97 * cholmod_nnz number of nonzeros in a sparse matrix 98 * cholmod_speye sparse identity matrix 99 * cholmod_spzeros sparse zero matrix 100 * cholmod_transpose transpose a sparse matrix 101 * cholmod_ptranspose transpose/permute a sparse matrix 102 * cholmod_transpose_unsym transpose/permute an unsymmetric sparse matrix 103 * cholmod_transpose_sym transpose/permute a symmetric sparse matrix 104 * cholmod_sort sort row indices in each column of sparse matrix 105 * cholmod_band C = tril (triu (A,k1), k2) 106 * cholmod_band_inplace A = tril (triu (A,k1), k2) 107 * cholmod_aat C = A*A' 108 * cholmod_copy_sparse C = A, create an exact copy of a sparse matrix 109 * cholmod_copy C = A, with possible change of stype 110 * cholmod_add C = alpha*A + beta*B 111 * cholmod_sparse_xtype change the xtype of a sparse matrix 112 * 113 * ============================================================================ 114 * === cholmod_factor ========================================================= 115 * ============================================================================ 116 * 117 * The data structure for an LL' or LDL' factorization is too complex to 118 * describe in one sentence. This object can hold the symbolic analysis alone, 119 * or in combination with a "simplicial" (similar to a sparse matrix) or 120 * "supernodal" form of the numerical factorization. Only the routine to free 121 * a factor is primary, since a factor object is created by the factorization 122 * routine (cholmod_factorize). It must be freed with cholmod_free_factor. 123 * 124 * cholmod_free_factor free a factor 125 * ----------------------------- 126 * cholmod_allocate_factor allocate a factor (LL' or LDL') 127 * cholmod_reallocate_factor change the # entries in a factor 128 * cholmod_change_factor change the type of factor (e.g., LDL' to LL') 129 * cholmod_pack_factor pack the columns of a factor 130 * cholmod_reallocate_column resize a single column of a factor 131 * cholmod_factor_to_sparse create a sparse matrix copy of a factor 132 * cholmod_copy_factor create a copy of a factor 133 * cholmod_factor_xtype change the xtype of a factor 134 * 135 * Note that there is no cholmod_sparse_to_factor routine to create a factor 136 * as a copy of a sparse matrix. It could be done, after a fashion, but a 137 * lower triangular sparse matrix would not necessarily have a chordal graph, 138 * which would break the many CHOLMOD routines that rely on this property. 139 * 140 * ============================================================================ 141 * === cholmod_dense ========================================================== 142 * ============================================================================ 143 * 144 * The solve routines and some of the MatrixOps and Modify routines use dense 145 * matrices as inputs. These are held in column-major order. With a leading 146 * dimension of d, the entry in row i and column j is held in x [i+j*d]. 147 * 148 * cholmod_allocate_dense allocate a dense matrix 149 * cholmod_free_dense free a dense matrix 150 * ----------------------------- 151 * cholmod_zeros allocate a dense matrix of all zeros 152 * cholmod_ones allocate a dense matrix of all ones 153 * cholmod_eye allocate a dense identity matrix 154 * cholmod_sparse_to_dense create a dense matrix copy of a sparse matrix 155 * cholmod_dense_to_sparse create a sparse matrix copy of a dense matrix 156 * cholmod_copy_dense create a copy of a dense matrix 157 * cholmod_copy_dense2 copy a dense matrix (pre-allocated) 158 * cholmod_dense_xtype change the xtype of a dense matrix 159 * cholmod_ensure_dense ensure a dense matrix has a given size and type 160 * 161 * ============================================================================ 162 * === cholmod_triplet ======================================================== 163 * ============================================================================ 164 * 165 * A sparse matrix held in triplet form is the simplest one for a user to 166 * create. It consists of a list of nz entries in arbitrary order, held in 167 * three arrays: i, j, and x, each of length nk. The kth entry is in row i[k], 168 * column j[k], with value x[k]. There may be duplicate values; if A(i,j) 169 * appears more than once, its value is the sum of the entries with those row 170 * and column indices. 171 * 172 * cholmod_allocate_triplet allocate a triplet matrix 173 * cholmod_triplet_to_sparse create a sparse matrix copy of a triplet matrix 174 * cholmod_free_triplet free a triplet matrix 175 * ----------------------------- 176 * cholmod_reallocate_triplet change the # of entries in a triplet matrix 177 * cholmod_sparse_to_triplet create a triplet matrix copy of a sparse matrix 178 * cholmod_copy_triplet create a copy of a triplet matrix 179 * cholmod_triplet_xtype change the xtype of a triplet matrix 180 * 181 * ============================================================================ 182 * === memory management ====================================================== 183 * ============================================================================ 184 * 185 * cholmod_malloc malloc wrapper 186 * cholmod_calloc calloc wrapper 187 * cholmod_free free wrapper 188 * cholmod_realloc realloc wrapper 189 * cholmod_realloc_multiple realloc wrapper for multiple objects 190 * 191 * ============================================================================ 192 * === Core CHOLMOD prototypes ================================================ 193 * ============================================================================ 194 * 195 * All CHOLMOD routines (in all modules) use the following protocol for return 196 * values, with one exception: 197 * 198 * int TRUE (1) if successful, or FALSE (0) otherwise. 199 * (exception: cholmod_divcomplex) 200 * SuiteSparse_long a value >= 0 if successful, or -1 otherwise. 201 * double a value >= 0 if successful, or -1 otherwise. 202 * size_t a value > 0 if successful, or 0 otherwise. 203 * void * a non-NULL pointer to newly allocated memory if 204 * successful, or NULL otherwise. 205 * cholmod_sparse * a non-NULL pointer to a newly allocated matrix 206 * if successful, or NULL otherwise. 207 * cholmod_factor * a non-NULL pointer to a newly allocated factor 208 * if successful, or NULL otherwise. 209 * cholmod_triplet * a non-NULL pointer to a newly allocated triplet 210 * matrix if successful, or NULL otherwise. 211 * cholmod_dense * a non-NULL pointer to a newly allocated triplet 212 * matrix if successful, or NULL otherwise. 213 * 214 * The last parameter to all routines is always a pointer to the CHOLMOD 215 * Common object. 216 * 217 * TRUE and FALSE are not defined here, since they may conflict with the user 218 * program. A routine that described here returning TRUE or FALSE returns 1 219 * or 0, respectively. Any TRUE/FALSE parameter is true if nonzero, false if 220 * zero. 221 */ 222 223 #ifndef CHOLMOD_CORE_H 224 #define CHOLMOD_CORE_H 225 226 /* ========================================================================== */ 227 /* === CHOLMOD version ====================================================== */ 228 /* ========================================================================== */ 229 230 /* All versions of CHOLMOD will include the following definitions. 231 * As an example, to test if the version you are using is 1.3 or later: 232 * 233 * if (CHOLMOD_VERSION >= CHOLMOD_VER_CODE (1,3)) ... 234 * 235 * This also works during compile-time: 236 * 237 * #if CHOLMOD_VERSION >= CHOLMOD_VER_CODE (1,3) 238 * printf ("This is version 1.3 or later\n") ; 239 * #else 240 * printf ("This is version is earlier than 1.3\n") ; 241 * #endif 242 */ 243 244 #define CHOLMOD_HAS_VERSION_FUNCTION 245 246 #define CHOLMOD_DATE "Oct 22, 2019" 247 #define CHOLMOD_VER_CODE(main,sub) ((main) * 1000 + (sub)) 248 #define CHOLMOD_MAIN_VERSION 3 249 #define CHOLMOD_SUB_VERSION 0 250 #define CHOLMOD_SUBSUB_VERSION 14 251 #define CHOLMOD_VERSION \ 252 CHOLMOD_VER_CODE(CHOLMOD_MAIN_VERSION,CHOLMOD_SUB_VERSION) 253 254 255 /* ========================================================================== */ 256 /* === non-CHOLMOD include files ============================================ */ 257 /* ========================================================================== */ 258 259 /* This is the only non-CHOLMOD include file imposed on the user program. 260 * It required for size_t definition used here. CHOLMOD itself includes other 261 * ANSI C89 standard #include files, but does not expose them to the user. 262 * 263 * CHOLMOD assumes that your C compiler is ANSI C89 compliant. It does not make 264 * use of ANSI C99 features. 265 */ 266 267 #include <stddef.h> 268 #include <stdlib.h> 269 270 /* ========================================================================== */ 271 /* === CUDA BLAS for the GPU ================================================ */ 272 /* ========================================================================== */ 273 274 /* The number of OMP threads should typically be set to the number of cores */ 275 /* per socket inthe machine being used. This maximizes memory performance. */ 276 #ifndef CHOLMOD_OMP_NUM_THREADS 277 #define CHOLMOD_OMP_NUM_THREADS 4 278 #endif 279 280 /* Define buffering parameters for GPU processing */ 281 #ifndef SUITESPARSE_GPU_EXTERN_ON 282 #ifdef GPU_BLAS 283 #include <cublas_v2.h> 284 #endif 285 #endif 286 287 #define CHOLMOD_DEVICE_SUPERNODE_BUFFERS 6 288 #define CHOLMOD_HOST_SUPERNODE_BUFFERS 8 289 #define CHOLMOD_DEVICE_STREAMS 2 290 291 /* ========================================================================== */ 292 /* === CHOLMOD objects ====================================================== */ 293 /* ========================================================================== */ 294 295 /* Each CHOLMOD object has its own type code. */ 296 297 #define CHOLMOD_COMMON 0 298 #define CHOLMOD_SPARSE 1 299 #define CHOLMOD_FACTOR 2 300 #define CHOLMOD_DENSE 3 301 #define CHOLMOD_TRIPLET 4 302 303 /* ========================================================================== */ 304 /* === CHOLMOD Common ======================================================= */ 305 /* ========================================================================== */ 306 307 /* itype defines the types of integer used: */ 308 #define CHOLMOD_INT 0 /* all integer arrays are int */ 309 #define CHOLMOD_INTLONG 1 /* most are int, some are SuiteSparse_long */ 310 #define CHOLMOD_LONG 2 /* all integer arrays are SuiteSparse_long */ 311 312 /* The itype of all parameters for all CHOLMOD routines must match. 313 * FUTURE WORK: CHOLMOD_INTLONG is not yet supported. 314 */ 315 316 /* dtype defines what the numerical type is (double or float): */ 317 #define CHOLMOD_DOUBLE 0 /* all numerical values are double */ 318 #define CHOLMOD_SINGLE 1 /* all numerical values are float */ 319 320 /* The dtype of all parameters for all CHOLMOD routines must match. 321 * 322 * Scalar floating-point values are always passed as double arrays of size 2 323 * (for the real and imaginary parts). They are typecast to float as needed. 324 * FUTURE WORK: the float case is not supported yet. 325 */ 326 327 /* xtype defines the kind of numerical values used: */ 328 #define CHOLMOD_PATTERN 0 /* pattern only, no numerical values */ 329 #define CHOLMOD_REAL 1 /* a real matrix */ 330 #define CHOLMOD_COMPLEX 2 /* a complex matrix (ANSI C99 compatible) */ 331 #define CHOLMOD_ZOMPLEX 3 /* a complex matrix (MATLAB compatible) */ 332 333 /* The xtype of all parameters for all CHOLMOD routines must match. 334 * 335 * CHOLMOD_PATTERN: x and z are ignored. 336 * CHOLMOD_DOUBLE: x is non-null of size nzmax, z is ignored. 337 * CHOLMOD_COMPLEX: x is non-null of size 2*nzmax doubles, z is ignored. 338 * CHOLMOD_ZOMPLEX: x and z are non-null of size nzmax 339 * 340 * In the real case, z is ignored. The kth entry in the matrix is x [k]. 341 * There are two methods for the complex case. In the ANSI C99-compatible 342 * CHOLMOD_COMPLEX case, the real and imaginary parts of the kth entry 343 * are in x [2*k] and x [2*k+1], respectively. z is ignored. In the 344 * MATLAB-compatible CHOLMOD_ZOMPLEX case, the real and imaginary 345 * parts of the kth entry are in x [k] and z [k]. 346 * 347 * Scalar floating-point values are always passed as double arrays of size 2 348 * (real and imaginary parts). The imaginary part of a scalar is ignored if 349 * the routine operates on a real matrix. 350 * 351 * These Modules support complex and zomplex matrices, with a few exceptions: 352 * 353 * Check all routines 354 * Cholesky all routines 355 * Core all except cholmod_aat, add, band, copy 356 * Demo all routines 357 * Partition all routines 358 * Supernodal all routines support any real, complex, or zomplex input. 359 * There will never be a supernodal zomplex L; a complex 360 * supernodal L is created if A is zomplex. 361 * Tcov all routines 362 * Valgrind all routines 363 * 364 * These Modules provide partial support for complex and zomplex matrices: 365 * 366 * MATLAB all routines support real and zomplex only, not complex, 367 * with the exception of ldlupdate, which supports 368 * real matrices only. This is a minor constraint since 369 * MATLAB's matrices are all real or zomplex. 370 * MatrixOps only norm_dense, norm_sparse, and sdmult support complex 371 * and zomplex 372 * 373 * These Modules do not support complex and zomplex matrices at all: 374 * 375 * Modify all routines support real matrices only 376 */ 377 378 /* Definitions for cholmod_common: */ 379 #define CHOLMOD_MAXMETHODS 9 /* maximum number of different methods that */ 380 /* cholmod_analyze can try. Must be >= 9. */ 381 382 /* Common->status values. zero means success, negative means a fatal error, 383 * positive is a warning. */ 384 #define CHOLMOD_OK 0 /* success */ 385 #define CHOLMOD_NOT_INSTALLED (-1) /* failure: method not installed */ 386 #define CHOLMOD_OUT_OF_MEMORY (-2) /* failure: out of memory */ 387 #define CHOLMOD_TOO_LARGE (-3) /* failure: integer overflow occured */ 388 #define CHOLMOD_INVALID (-4) /* failure: invalid input */ 389 #define CHOLMOD_GPU_PROBLEM (-5) /* failure: GPU fatal error */ 390 #define CHOLMOD_NOT_POSDEF (1) /* warning: matrix not pos. def. */ 391 #define CHOLMOD_DSMALL (2) /* warning: D for LDL' or diag(L) or */ 392 /* LL' has tiny absolute value */ 393 394 /* ordering method (also used for L->ordering) */ 395 #define CHOLMOD_NATURAL 0 /* use natural ordering */ 396 #define CHOLMOD_GIVEN 1 /* use given permutation */ 397 #define CHOLMOD_AMD 2 /* use minimum degree (AMD) */ 398 #define CHOLMOD_METIS 3 /* use METIS' nested dissection */ 399 #define CHOLMOD_NESDIS 4 /* use CHOLMOD's version of nested dissection:*/ 400 /* node bisector applied recursively, followed 401 * by constrained minimum degree (CSYMAMD or 402 * CCOLAMD) */ 403 #define CHOLMOD_COLAMD 5 /* use AMD for A, COLAMD for A*A' */ 404 405 /* POSTORDERED is not a method, but a result of natural ordering followed by a 406 * weighted postorder. It is used for L->ordering, not method [ ].ordering. */ 407 #define CHOLMOD_POSTORDERED 6 /* natural ordering, postordered. */ 408 409 /* supernodal strategy (for Common->supernodal) */ 410 #define CHOLMOD_SIMPLICIAL 0 /* always do simplicial */ 411 #define CHOLMOD_AUTO 1 /* select simpl/super depending on matrix */ 412 #define CHOLMOD_SUPERNODAL 2 /* always do supernodal */ 413 414 typedef struct cholmod_common_struct 415 { 416 /* ---------------------------------------------------------------------- */ 417 /* parameters for symbolic/numeric factorization and update/downdate */ 418 /* ---------------------------------------------------------------------- */ 419 420 double dbound ; /* Smallest absolute value of diagonal entries of D 421 * for LDL' factorization and update/downdate/rowadd/ 422 * rowdel, or the diagonal of L for an LL' factorization. 423 * Entries in the range 0 to dbound are replaced with dbound. 424 * Entries in the range -dbound to 0 are replaced with -dbound. No 425 * changes are made to the diagonal if dbound <= 0. Default: zero */ 426 427 double grow0 ; /* For a simplicial factorization, L->i and L->x can 428 * grow if necessary. grow0 is the factor by which 429 * it grows. For the initial space, L is of size MAX (1,grow0) times 430 * the required space. If L runs out of space, the new size of L is 431 * MAX(1.2,grow0) times the new required space. If you do not plan on 432 * modifying the LDL' factorization in the Modify module, set grow0 to 433 * zero (or set grow2 to 0, see below). Default: 1.2 */ 434 435 double grow1 ; 436 437 size_t grow2 ; /* For a simplicial factorization, each column j of L 438 * is initialized with space equal to 439 * grow1*L->ColCount[j] + grow2. If grow0 < 1, grow1 < 1, or grow2 == 0, 440 * then the space allocated is exactly equal to L->ColCount[j]. If the 441 * column j runs out of space, it increases to grow1*need + grow2 in 442 * size, where need is the total # of nonzeros in that column. If you do 443 * not plan on modifying the factorization in the Modify module, set 444 * grow2 to zero. Default: grow1 = 1.2, grow2 = 5. */ 445 446 size_t maxrank ; /* rank of maximum update/downdate. Valid values: 447 * 2, 4, or 8. A value < 2 is set to 2, and a 448 * value > 8 is set to 8. It is then rounded up to the next highest 449 * power of 2, if not already a power of 2. Workspace (Xwork, below) of 450 * size nrow-by-maxrank double's is allocated for the update/downdate. 451 * If an update/downdate of rank-k is requested, with k > maxrank, 452 * it is done in steps of maxrank. Default: 8, which is fastest. 453 * Memory usage can be reduced by setting maxrank to 2 or 4. 454 */ 455 456 double supernodal_switch ; /* supernodal vs simplicial factorization */ 457 int supernodal ; /* If Common->supernodal <= CHOLMOD_SIMPLICIAL 458 * (0) then cholmod_analyze performs a 459 * simplicial analysis. If >= CHOLMOD_SUPERNODAL (2), then a supernodal 460 * analysis is performed. If == CHOLMOD_AUTO (1) and 461 * flop/nnz(L) < Common->supernodal_switch, then a simplicial analysis 462 * is done. A supernodal analysis done otherwise. 463 * Default: CHOLMOD_AUTO. Default supernodal_switch = 40 */ 464 465 int final_asis ; /* If TRUE, then ignore the other final_* parameters 466 * (except for final_pack). 467 * The factor is left as-is when done. Default: TRUE.*/ 468 469 int final_super ; /* If TRUE, leave a factor in supernodal form when 470 * supernodal factorization is finished. If FALSE, 471 * then convert to a simplicial factor when done. 472 * Default: TRUE */ 473 474 int final_ll ; /* If TRUE, leave factor in LL' form when done. 475 * Otherwise, leave in LDL' form. Default: FALSE */ 476 477 int final_pack ; /* If TRUE, pack the columns when done. If TRUE, and 478 * cholmod_factorize is called with a symbolic L, L is 479 * allocated with exactly the space required, using L->ColCount. If you 480 * plan on modifying the factorization, set Common->final_pack to FALSE, 481 * and each column will be given a little extra slack space for future 482 * growth in fill-in due to updates. Default: TRUE */ 483 484 int final_monotonic ; /* If TRUE, ensure columns are monotonic when done. 485 * Default: TRUE */ 486 487 int final_resymbol ;/* if cholmod_factorize performed a supernodal 488 * factorization, final_resymbol is true, and 489 * final_super is FALSE (convert a simplicial numeric factorization), 490 * then numerically zero entries that resulted from relaxed supernodal 491 * amalgamation are removed. This does not remove entries that are zero 492 * due to exact numeric cancellation, since doing so would break the 493 * update/downdate rowadd/rowdel routines. Default: FALSE. */ 494 495 /* supernodal relaxed amalgamation parameters: */ 496 double zrelax [3] ; 497 size_t nrelax [3] ; 498 499 /* Let ns be the total number of columns in two adjacent supernodes. 500 * Let z be the fraction of zero entries in the two supernodes if they 501 * are merged (z includes zero entries from prior amalgamations). The 502 * two supernodes are merged if: 503 * (ns <= nrelax [0]) || (no new zero entries added) || 504 * (ns <= nrelax [1] && z < zrelax [0]) || 505 * (ns <= nrelax [2] && z < zrelax [1]) || (z < zrelax [2]) 506 * 507 * Default parameters result in the following rule: 508 * (ns <= 4) || (no new zero entries added) || 509 * (ns <= 16 && z < 0.8) || (ns <= 48 && z < 0.1) || (z < 0.05) 510 */ 511 512 int prefer_zomplex ; /* X = cholmod_solve (sys, L, B, Common) computes 513 * x=A\b or solves a related system. If L and B are 514 * both real, then X is real. Otherwise, X is returned as 515 * CHOLMOD_COMPLEX if Common->prefer_zomplex is FALSE, or 516 * CHOLMOD_ZOMPLEX if Common->prefer_zomplex is TRUE. This parameter 517 * is needed because there is no supernodal zomplex L. Suppose the 518 * caller wants all complex matrices to be stored in zomplex form 519 * (MATLAB, for example). A supernodal L is returned in complex form 520 * if A is zomplex. B can be real, and thus X = cholmod_solve (L,B) 521 * should return X as zomplex. This cannot be inferred from the input 522 * arguments L and B. Default: FALSE, since all data types are 523 * supported in CHOLMOD_COMPLEX form and since this is the native type 524 * of LAPACK and the BLAS. Note that the MATLAB/cholmod.c mexFunction 525 * sets this parameter to TRUE, since MATLAB matrices are in 526 * CHOLMOD_ZOMPLEX form. 527 */ 528 529 int prefer_upper ; /* cholmod_analyze and cholmod_factorize work 530 * fastest when a symmetric matrix is stored in 531 * upper triangular form when a fill-reducing ordering is used. In 532 * MATLAB, this corresponds to how x=A\b works. When the matrix is 533 * ordered as-is, they work fastest when a symmetric matrix is in lower 534 * triangular form. In MATLAB, R=chol(A) does the opposite. This 535 * parameter affects only how cholmod_read returns a symmetric matrix. 536 * If TRUE (the default case), a symmetric matrix is always returned in 537 * upper-triangular form (A->stype = 1). */ 538 539 int quick_return_if_not_posdef ; /* if TRUE, the supernodal numeric 540 * factorization will return quickly if 541 * the matrix is not positive definite. Default: FALSE. */ 542 543 int prefer_binary ; /* cholmod_read_triplet converts a symmetric 544 * pattern-only matrix into a real matrix. If 545 * prefer_binary is FALSE, the diagonal entries are set to 1 + the degree 546 * of the row/column, and off-diagonal entries are set to -1 (resulting 547 * in a positive definite matrix if the diagonal is zero-free). Most 548 * symmetric patterns are the pattern a positive definite matrix. If 549 * this parameter is TRUE, then the matrix is returned with a 1 in each 550 * entry, instead. Default: FALSE. Added in v1.3. */ 551 552 /* ---------------------------------------------------------------------- */ 553 /* printing and error handling options */ 554 /* ---------------------------------------------------------------------- */ 555 556 int print ; /* print level. Default: 3 */ 557 int precise ; /* if TRUE, print 16 digits. Otherwise print 5 */ 558 559 /* CHOLMOD print_function replaced with SuiteSparse_config.print_func */ 560 561 int try_catch ; /* if TRUE, then ignore errors; CHOLMOD is in the middle 562 * of a try/catch block. No error message is printed 563 * and the Common->error_handler function is not called. */ 564 565 void (*error_handler) (int status, const char *file, 566 int line, const char *message) ; 567 568 /* Common->error_handler is the user's error handling routine. If not 569 * NULL, this routine is called if an error occurs in CHOLMOD. status 570 * can be CHOLMOD_OK (0), negative for a fatal error, and positive for 571 * a warning. file is a string containing the name of the source code 572 * file where the error occured, and line is the line number in that 573 * file. message is a string describing the error in more detail. */ 574 575 /* ---------------------------------------------------------------------- */ 576 /* ordering options */ 577 /* ---------------------------------------------------------------------- */ 578 579 /* The cholmod_analyze routine can try many different orderings and select 580 * the best one. It can also try one ordering method multiple times, with 581 * different parameter settings. The default is to use three orderings, 582 * the user's permutation (if provided), AMD which is the fastest ordering 583 * and generally gives good fill-in, and METIS. CHOLMOD's nested dissection 584 * (METIS with a constrained AMD) usually gives a better ordering than METIS 585 * alone (by about 5% to 10%) but it takes more time. 586 * 587 * If you know the method that is best for your matrix, set Common->nmethods 588 * to 1 and set Common->method [0] to the set of parameters for that method. 589 * If you set it to 1 and do not provide a permutation, then only AMD will 590 * be called. 591 * 592 * If METIS is not available, the default # of methods tried is 2 (the user 593 * permutation, if any, and AMD). 594 * 595 * To try other methods, set Common->nmethods to the number of methods you 596 * want to try. The suite of default methods and their parameters is 597 * described in the cholmod_defaults routine, and summarized here: 598 * 599 * Common->method [i]: 600 * i = 0: user-provided ordering (cholmod_analyze_p only) 601 * i = 1: AMD (for both A and A*A') 602 * i = 2: METIS 603 * i = 3: CHOLMOD's nested dissection (NESDIS), default parameters 604 * i = 4: natural 605 * i = 5: NESDIS with nd_small = 20000 606 * i = 6: NESDIS with nd_small = 4, no constrained minimum degree 607 * i = 7: NESDIS with no dense node removal 608 * i = 8: AMD for A, COLAMD for A*A' 609 * 610 * You can modify the suite of methods you wish to try by modifying 611 * Common.method [...] after calling cholmod_start or cholmod_defaults. 612 * 613 * For example, to use AMD, followed by a weighted postordering: 614 * 615 * Common->nmethods = 1 ; 616 * Common->method [0].ordering = CHOLMOD_AMD ; 617 * Common->postorder = TRUE ; 618 * 619 * To use the natural ordering (with no postordering): 620 * 621 * Common->nmethods = 1 ; 622 * Common->method [0].ordering = CHOLMOD_NATURAL ; 623 * Common->postorder = FALSE ; 624 * 625 * If you are going to factorize hundreds or more matrices with the same 626 * nonzero pattern, you may wish to spend a great deal of time finding a 627 * good permutation. In this case, try setting Common->nmethods to 9. 628 * The time spent in cholmod_analysis will be very high, but you need to 629 * call it only once. 630 * 631 * cholmod_analyze sets Common->current to a value between 0 and nmethods-1. 632 * Each ordering method uses the set of options defined by this parameter. 633 */ 634 635 int nmethods ; /* The number of ordering methods to try. Default: 0. 636 * nmethods = 0 is a special case. cholmod_analyze 637 * will try the user-provided ordering (if given) and AMD. Let fl and 638 * lnz be the flop count and nonzeros in L from AMD's ordering. Let 639 * anz be the number of nonzeros in the upper or lower triangular part 640 * of the symmetric matrix A. If fl/lnz < 500 or lnz/anz < 5, then this 641 * is a good ordering, and METIS is not attempted. Otherwise, METIS is 642 * tried. The best ordering found is used. If nmethods > 0, the 643 * methods used are given in the method[ ] array, below. The first 644 * three methods in the default suite of orderings is (1) use the given 645 * permutation (if provided), (2) use AMD, and (3) use METIS. Maximum 646 * allowed value is CHOLMOD_MAXMETHODS. */ 647 648 int current ; /* The current method being tried. Default: 0. Valid 649 * range is 0 to nmethods-1. */ 650 651 int selected ; /* The best method found. */ 652 653 /* The suite of ordering methods and parameters: */ 654 655 struct cholmod_method_struct 656 { 657 /* statistics for this method */ 658 double lnz ; /* nnz(L) excl. zeros from supernodal amalgamation, 659 * for a "pure" L */ 660 661 double fl ; /* flop count for a "pure", real simplicial LL' 662 * factorization, with no extra work due to 663 * amalgamation. Subtract n to get the LDL' flop count. Multiply 664 * by about 4 if the matrix is complex or zomplex. */ 665 666 /* ordering method parameters */ 667 double prune_dense ;/* dense row/col control for AMD, SYMAMD, CSYMAMD, 668 * and NESDIS (cholmod_nested_dissection). For a 669 * symmetric n-by-n matrix, rows/columns with more than 670 * MAX (16, prune_dense * sqrt (n)) entries are removed prior to 671 * ordering. They appear at the end of the re-ordered matrix. 672 * 673 * If prune_dense < 0, only completely dense rows/cols are removed. 674 * 675 * This paramater is also the dense column control for COLAMD and 676 * CCOLAMD. For an m-by-n matrix, columns with more than 677 * MAX (16, prune_dense * sqrt (MIN (m,n))) entries are removed prior 678 * to ordering. They appear at the end of the re-ordered matrix. 679 * CHOLMOD factorizes A*A', so it calls COLAMD and CCOLAMD with A', 680 * not A. Thus, this parameter affects the dense *row* control for 681 * CHOLMOD's matrix, and the dense *column* control for COLAMD and 682 * CCOLAMD. 683 * 684 * Removing dense rows and columns improves the run-time of the 685 * ordering methods. It has some impact on ordering quality 686 * (usually minimal, sometimes good, sometimes bad). 687 * 688 * Default: 10. */ 689 690 double prune_dense2 ;/* dense row control for COLAMD and CCOLAMD. 691 * Rows with more than MAX (16, dense2 * sqrt (n)) 692 * for an m-by-n matrix are removed prior to ordering. CHOLMOD's 693 * matrix is transposed before ordering it with COLAMD or CCOLAMD, 694 * so this controls the dense *columns* of CHOLMOD's matrix, and 695 * the dense *rows* of COLAMD's or CCOLAMD's matrix. 696 * 697 * If prune_dense2 < 0, only completely dense rows/cols are removed. 698 * 699 * Default: -1. Note that this is not the default for COLAMD and 700 * CCOLAMD. -1 is best for Cholesky. 10 is best for LU. */ 701 702 double nd_oksep ; /* in NESDIS, when a node separator is computed, it 703 * discarded if nsep >= nd_oksep*n, where nsep is 704 * the number of nodes in the separator, and n is the size of the 705 * graph being cut. Valid range is 0 to 1. If 1 or greater, the 706 * separator is discarded if it consists of the entire graph. 707 * Default: 1 */ 708 709 double other_1 [4] ; /* future expansion */ 710 711 size_t nd_small ; /* do not partition graphs with fewer nodes than 712 * nd_small, in NESDIS. Default: 200 (same as 713 * METIS) */ 714 715 size_t other_2 [4] ; /* future expansion */ 716 717 int aggressive ; /* Aggresive absorption in AMD, COLAMD, SYMAMD, 718 * CCOLAMD, and CSYMAMD. Default: TRUE */ 719 720 int order_for_lu ; /* CCOLAMD can be optimized to produce an ordering 721 * for LU or Cholesky factorization. CHOLMOD only 722 * performs a Cholesky factorization. However, you may wish to use 723 * CHOLMOD as an interface for CCOLAMD but use it for your own LU 724 * factorization. In this case, order_for_lu should be set to FALSE. 725 * When factorizing in CHOLMOD itself, you should *** NEVER *** set 726 * this parameter FALSE. Default: TRUE. */ 727 728 int nd_compress ; /* If TRUE, compress the graph and subgraphs before 729 * partitioning them in NESDIS. Default: TRUE */ 730 731 int nd_camd ; /* If 1, follow the nested dissection ordering 732 * with a constrained minimum degree ordering that 733 * respects the partitioning just found (using CAMD). If 2, use 734 * CSYMAMD instead. If you set nd_small very small, you may not need 735 * this ordering, and can save time by setting it to zero (no 736 * constrained minimum degree ordering). Default: 1. */ 737 738 int nd_components ; /* The nested dissection ordering finds a node 739 * separator that splits the graph into two parts, 740 * which may be unconnected. If nd_components is TRUE, each of 741 * these connected components is split independently. If FALSE, 742 * each part is split as a whole, even if it consists of more than 743 * one connected component. Default: FALSE */ 744 745 /* fill-reducing ordering to use */ 746 int ordering ; 747 748 size_t other_3 [4] ; /* future expansion */ 749 750 } method [CHOLMOD_MAXMETHODS + 1] ; 751 752 int postorder ; /* If TRUE, cholmod_analyze follows the ordering with a 753 * weighted postorder of the elimination tree. Improves 754 * supernode amalgamation. Does not affect fundamental nnz(L) and 755 * flop count. Default: TRUE. */ 756 757 int default_nesdis ; /* Default: FALSE. If FALSE, then the default 758 * ordering strategy (when Common->nmethods == 0) 759 * is to try the given ordering (if present), AMD, and then METIS if AMD 760 * reports high fill-in. If Common->default_nesdis is TRUE then NESDIS 761 * is used instead in the default strategy. */ 762 763 /* ---------------------------------------------------------------------- */ 764 /* memory management, complex divide, and hypot function pointers moved */ 765 /* ---------------------------------------------------------------------- */ 766 767 /* Function pointers moved from here (in CHOLMOD 2.2.0) to 768 SuiteSparse_config.[ch]. See CHOLMOD/Include/cholmod_back.h 769 for a set of macros that can be #include'd or copied into your 770 application to define these function pointers on any version of CHOLMOD. 771 */ 772 773 /* ---------------------------------------------------------------------- */ 774 /* METIS workarounds */ 775 /* ---------------------------------------------------------------------- */ 776 777 /* These workarounds were put into place for METIS 4.0.1. They are safe 778 to use with METIS 5.1.0, but they might not longer be necessary. */ 779 780 double metis_memory ; /* This is a parameter for CHOLMOD's interface to 781 * METIS, not a parameter to METIS itself. METIS 782 * uses an amount of memory that is difficult to estimate precisely 783 * beforehand. If it runs out of memory, it terminates your program. 784 * All routines in CHOLMOD except for CHOLMOD's interface to METIS 785 * return an error status and safely return to your program if they run 786 * out of memory. To mitigate this problem, the CHOLMOD interface 787 * can allocate a single block of memory equal in size to an empirical 788 * upper bound of METIS's memory usage times the Common->metis_memory 789 * parameter, and then immediately free it. It then calls METIS. If 790 * this pre-allocation fails, it is possible that METIS will fail as 791 * well, and so CHOLMOD returns with an out-of-memory condition without 792 * calling METIS. 793 * 794 * METIS_NodeND (used in the CHOLMOD_METIS ordering option) with its 795 * default parameter settings typically uses about (4*nz+40n+4096) 796 * times sizeof(int) memory, where nz is equal to the number of entries 797 * in A for the symmetric case or AA' if an unsymmetric matrix is 798 * being ordered (where nz includes both the upper and lower parts 799 * of A or AA'). The observed "upper bound" (with 2 exceptions), 800 * measured in an instrumented copy of METIS 4.0.1 on thousands of 801 * matrices, is (10*nz+50*n+4096) * sizeof(int). Two large matrices 802 * exceeded this bound, one by almost a factor of 2 (Gupta/gupta2). 803 * 804 * If your program is terminated by METIS, try setting metis_memory to 805 * 2.0, or even higher if needed. By default, CHOLMOD assumes that METIS 806 * does not have this problem (so that CHOLMOD will work correctly when 807 * this issue is fixed in METIS). Thus, the default value is zero. 808 * This work-around is not guaranteed anyway. 809 * 810 * If a matrix exceeds this predicted memory usage, AMD is attempted 811 * instead. It, too, may run out of memory, but if it does so it will 812 * not terminate your program. 813 */ 814 815 double metis_dswitch ; /* METIS_NodeND in METIS 4.0.1 gives a seg */ 816 size_t metis_nswitch ; /* fault with one matrix of order n = 3005 and 817 * nz = 6,036,025. This is a very dense graph. 818 * The workaround is to use AMD instead of METIS for matrices of dimension 819 * greater than Common->metis_nswitch (default 3000) or more and with 820 * density of Common->metis_dswitch (default 0.66) or more. 821 * cholmod_nested_dissection has no problems with the same matrix, even 822 * though it uses METIS_ComputeVertexSeparator on this matrix. If this 823 * seg fault does not affect you, set metis_nswitch to zero or less, 824 * and CHOLMOD will not switch to AMD based just on the density of the 825 * matrix (it will still switch to AMD if the metis_memory parameter 826 * causes the switch). 827 */ 828 829 /* ---------------------------------------------------------------------- */ 830 /* workspace */ 831 /* ---------------------------------------------------------------------- */ 832 833 /* CHOLMOD has several routines that take less time than the size of 834 * workspace they require. Allocating and initializing the workspace would 835 * dominate the run time, unless workspace is allocated and initialized 836 * just once. CHOLMOD allocates this space when needed, and holds it here 837 * between calls to CHOLMOD. cholmod_start sets these pointers to NULL 838 * (which is why it must be the first routine called in CHOLMOD). 839 * cholmod_finish frees the workspace (which is why it must be the last 840 * call to CHOLMOD). 841 */ 842 843 size_t nrow ; /* size of Flag and Head */ 844 SuiteSparse_long mark ; /* mark value for Flag array */ 845 size_t iworksize ; /* size of Iwork. Upper bound: 6*nrow+ncol */ 846 size_t xworksize ; /* size of Xwork, in bytes. 847 * maxrank*nrow*sizeof(double) for update/downdate. 848 * 2*nrow*sizeof(double) otherwise */ 849 850 /* initialized workspace: contents needed between calls to CHOLMOD */ 851 void *Flag ; /* size nrow, an integer array. Kept cleared between 852 * calls to cholmod rouines (Flag [i] < mark) */ 853 854 void *Head ; /* size nrow+1, an integer array. Kept cleared between 855 * calls to cholmod routines (Head [i] = EMPTY) */ 856 857 void *Xwork ; /* a double array. Its size varies. It is nrow for 858 * most routines (cholmod_rowfac, cholmod_add, 859 * cholmod_aat, cholmod_norm, cholmod_ssmult) for the real case, twice 860 * that when the input matrices are complex or zomplex. It is of size 861 * 2*nrow for cholmod_rowadd and cholmod_rowdel. For cholmod_updown, 862 * its size is maxrank*nrow where maxrank is 2, 4, or 8. Kept cleared 863 * between calls to cholmod (set to zero). */ 864 865 /* uninitialized workspace, contents not needed between calls to CHOLMOD */ 866 void *Iwork ; /* size iworksize, 2*nrow+ncol for most routines, 867 * up to 6*nrow+ncol for cholmod_analyze. */ 868 869 int itype ; /* If CHOLMOD_LONG, Flag, Head, and Iwork are 870 * SuiteSparse_long. Otherwise all three are int. */ 871 872 int dtype ; /* double or float */ 873 874 /* Common->itype and Common->dtype are used to define the types of all 875 * sparse matrices, triplet matrices, dense matrices, and factors 876 * created using this Common struct. The itypes and dtypes of all 877 * parameters to all CHOLMOD routines must match. */ 878 879 int no_workspace_reallocate ; /* this is an internal flag, used as a 880 * precaution by cholmod_analyze. It is normally false. If true, 881 * cholmod_allocate_work is not allowed to reallocate any workspace; 882 * they must use the existing workspace in Common (Iwork, Flag, Head, 883 * and Xwork). Added for CHOLMOD v1.1 */ 884 885 /* ---------------------------------------------------------------------- */ 886 /* statistics */ 887 /* ---------------------------------------------------------------------- */ 888 889 /* fl and lnz are set only in cholmod_analyze and cholmod_rowcolcounts, 890 * in the Cholesky modudle. modfl is set only in the Modify module. */ 891 892 int status ; /* error code */ 893 double fl ; /* LL' flop count from most recent analysis */ 894 double lnz ; /* fundamental nz in L */ 895 double anz ; /* nonzeros in tril(A) if A is symmetric/lower, 896 * triu(A) if symmetric/upper, or tril(A*A') if 897 * unsymmetric, in last call to cholmod_analyze. */ 898 double modfl ; /* flop count from most recent update/downdate/ 899 * rowadd/rowdel (excluding flops to modify the 900 * solution to Lx=b, if computed) */ 901 size_t malloc_count ; /* # of objects malloc'ed minus the # free'd*/ 902 size_t memory_usage ; /* peak memory usage in bytes */ 903 size_t memory_inuse ; /* current memory usage in bytes */ 904 905 double nrealloc_col ; /* # of column reallocations */ 906 double nrealloc_factor ;/* # of factor reallocations due to col. reallocs */ 907 double ndbounds_hit ; /* # of times diagonal modified by dbound */ 908 909 double rowfacfl ; /* # of flops in last call to cholmod_rowfac */ 910 double aatfl ; /* # of flops to compute A(:,f)*A(:,f)' */ 911 912 int called_nd ; /* TRUE if the last call to 913 * cholmod_analyze called NESDIS or METIS. */ 914 int blas_ok ; /* FALSE if BLAS int overflow; TRUE otherwise */ 915 916 /* ---------------------------------------------------------------------- */ 917 /* SuiteSparseQR control parameters: */ 918 /* ---------------------------------------------------------------------- */ 919 920 double SPQR_grain ; /* task size is >= max (total flops / grain) */ 921 double SPQR_small ; /* task size is >= small */ 922 int SPQR_shrink ; /* controls stack realloc method */ 923 int SPQR_nthreads ; /* number of TBB threads, 0 = auto */ 924 925 /* ---------------------------------------------------------------------- */ 926 /* SuiteSparseQR statistics */ 927 /* ---------------------------------------------------------------------- */ 928 929 /* was other1 [0:3] */ 930 double SPQR_flopcount ; /* flop count for SPQR */ 931 double SPQR_analyze_time ; /* analysis time in seconds for SPQR */ 932 double SPQR_factorize_time ; /* factorize time in seconds for SPQR */ 933 double SPQR_solve_time ; /* backsolve time in seconds */ 934 935 /* was SPQR_xstat [0:3] */ 936 double SPQR_flopcount_bound ; /* upper bound on flop count */ 937 double SPQR_tol_used ; /* tolerance used */ 938 double SPQR_norm_E_fro ; /* Frobenius norm of dropped entries */ 939 940 /* was SPQR_istat [0:9] */ 941 SuiteSparse_long SPQR_istat [10] ; 942 943 /* ---------------------------------------------------------------------- */ 944 /* GPU configuration and statistics */ 945 /* ---------------------------------------------------------------------- */ 946 947 /* useGPU: 1 if gpu-acceleration is requested */ 948 /* 0 if gpu-acceleration is prohibited */ 949 /* -1 if gpu-acceleration is undefined in which case the */ 950 /* environment CHOLMOD_USE_GPU will be queried and used. */ 951 /* useGPU=-1 is only used by CHOLMOD and treated as 0 by SPQR */ 952 int useGPU; 953 954 /* for CHOLMOD: */ 955 size_t maxGpuMemBytes; 956 double maxGpuMemFraction; 957 958 /* for SPQR: */ 959 size_t gpuMemorySize; /* Amount of memory in bytes on the GPU */ 960 double gpuKernelTime; /* Time taken by GPU kernels */ 961 SuiteSparse_long gpuFlops; /* Number of flops performed by the GPU */ 962 int gpuNumKernelLaunches; /* Number of GPU kernel launches */ 963 964 /* If not using the GPU, these items are not used, but they should be 965 present so that the CHOLMOD Common has the same size whether the GPU 966 is used or not. This way, all packages will agree on the size of 967 the CHOLMOD Common, regardless of whether or not they are compiled 968 with the GPU libraries or not */ 969 970 #ifdef GPU_BLAS 971 /* in CUDA, these three types are pointers */ 972 #define CHOLMOD_CUBLAS_HANDLE cublasHandle_t 973 #define CHOLMOD_CUDASTREAM cudaStream_t 974 #define CHOLMOD_CUDAEVENT cudaEvent_t 975 #else 976 /* ... so make them void * pointers if the GPU is not being used */ 977 #define CHOLMOD_CUBLAS_HANDLE void * 978 #define CHOLMOD_CUDASTREAM void * 979 #define CHOLMOD_CUDAEVENT void * 980 #endif 981 982 CHOLMOD_CUBLAS_HANDLE cublasHandle ; 983 984 /* a set of streams for general use */ 985 CHOLMOD_CUDASTREAM gpuStream[CHOLMOD_HOST_SUPERNODE_BUFFERS]; 986 987 CHOLMOD_CUDAEVENT cublasEventPotrf [3] ; 988 CHOLMOD_CUDAEVENT updateCKernelsComplete; 989 CHOLMOD_CUDAEVENT updateCBuffersFree[CHOLMOD_HOST_SUPERNODE_BUFFERS]; 990 991 void *dev_mempool; /* pointer to single allocation of device memory */ 992 size_t dev_mempool_size; 993 994 void *host_pinned_mempool; /* pointer to single allocation of pinned mem */ 995 size_t host_pinned_mempool_size; 996 997 size_t devBuffSize; 998 int ibuffer; 999 1000 double syrkStart ; /* time syrk started */ 1001 1002 /* run times of the different parts of CHOLMOD (GPU and CPU) */ 1003 double cholmod_cpu_gemm_time ; 1004 double cholmod_cpu_syrk_time ; 1005 double cholmod_cpu_trsm_time ; 1006 double cholmod_cpu_potrf_time ; 1007 double cholmod_gpu_gemm_time ; 1008 double cholmod_gpu_syrk_time ; 1009 double cholmod_gpu_trsm_time ; 1010 double cholmod_gpu_potrf_time ; 1011 double cholmod_assemble_time ; 1012 double cholmod_assemble_time2 ; 1013 1014 /* number of times the BLAS are called on the CPU and the GPU */ 1015 size_t cholmod_cpu_gemm_calls ; 1016 size_t cholmod_cpu_syrk_calls ; 1017 size_t cholmod_cpu_trsm_calls ; 1018 size_t cholmod_cpu_potrf_calls ; 1019 size_t cholmod_gpu_gemm_calls ; 1020 size_t cholmod_gpu_syrk_calls ; 1021 size_t cholmod_gpu_trsm_calls ; 1022 size_t cholmod_gpu_potrf_calls ; 1023 1024 } cholmod_common ; 1025 1026 /* size_t BLAS statistcs in Common: */ 1027 #define CHOLMOD_CPU_GEMM_CALLS cholmod_cpu_gemm_calls 1028 #define CHOLMOD_CPU_SYRK_CALLS cholmod_cpu_syrk_calls 1029 #define CHOLMOD_CPU_TRSM_CALLS cholmod_cpu_trsm_calls 1030 #define CHOLMOD_CPU_POTRF_CALLS cholmod_cpu_potrf_calls 1031 #define CHOLMOD_GPU_GEMM_CALLS cholmod_gpu_gemm_calls 1032 #define CHOLMOD_GPU_SYRK_CALLS cholmod_gpu_syrk_calls 1033 #define CHOLMOD_GPU_TRSM_CALLS cholmod_gpu_trsm_calls 1034 #define CHOLMOD_GPU_POTRF_CALLS cholmod_gpu_potrf_calls 1035 1036 /* double BLAS statistics in Common: */ 1037 #define CHOLMOD_CPU_GEMM_TIME cholmod_cpu_gemm_time 1038 #define CHOLMOD_CPU_SYRK_TIME cholmod_cpu_syrk_time 1039 #define CHOLMOD_CPU_TRSM_TIME cholmod_cpu_trsm_time 1040 #define CHOLMOD_CPU_POTRF_TIME cholmod_cpu_potrf_time 1041 #define CHOLMOD_GPU_GEMM_TIME cholmod_gpu_gemm_time 1042 #define CHOLMOD_GPU_SYRK_TIME cholmod_gpu_syrk_time 1043 #define CHOLMOD_GPU_TRSM_TIME cholmod_gpu_trsm_time 1044 #define CHOLMOD_GPU_POTRF_TIME cholmod_gpu_potrf_time 1045 #define CHOLMOD_ASSEMBLE_TIME cholmod_assemble_time 1046 #define CHOLMOD_ASSEMBLE_TIME2 cholmod_assemble_time2 1047 1048 /* for supernodal analysis */ 1049 #define CHOLMOD_ANALYZE_FOR_SPQR 0 1050 #define CHOLMOD_ANALYZE_FOR_CHOLESKY 1 1051 #define CHOLMOD_ANALYZE_FOR_SPQRGPU 2 1052 1053 /* -------------------------------------------------------------------------- */ 1054 /* cholmod_start: first call to CHOLMOD */ 1055 /* -------------------------------------------------------------------------- */ 1056 1057 int cholmod_start 1058 ( 1059 cholmod_common *Common 1060 ) ; 1061 1062 int cholmod_l_start (cholmod_common *) ; 1063 1064 /* -------------------------------------------------------------------------- */ 1065 /* cholmod_finish: last call to CHOLMOD */ 1066 /* -------------------------------------------------------------------------- */ 1067 1068 int cholmod_finish 1069 ( 1070 cholmod_common *Common 1071 ) ; 1072 1073 int cholmod_l_finish (cholmod_common *) ; 1074 1075 /* -------------------------------------------------------------------------- */ 1076 /* cholmod_defaults: restore default parameters */ 1077 /* -------------------------------------------------------------------------- */ 1078 1079 int cholmod_defaults 1080 ( 1081 cholmod_common *Common 1082 ) ; 1083 1084 int cholmod_l_defaults (cholmod_common *) ; 1085 1086 /* -------------------------------------------------------------------------- */ 1087 /* cholmod_maxrank: return valid maximum rank for update/downdate */ 1088 /* -------------------------------------------------------------------------- */ 1089 1090 size_t cholmod_maxrank /* returns validated value of Common->maxrank */ 1091 ( 1092 /* ---- input ---- */ 1093 size_t n, /* A and L will have n rows */ 1094 /* --------------- */ 1095 cholmod_common *Common 1096 ) ; 1097 1098 size_t cholmod_l_maxrank (size_t, cholmod_common *) ; 1099 1100 /* -------------------------------------------------------------------------- */ 1101 /* cholmod_allocate_work: allocate workspace in Common */ 1102 /* -------------------------------------------------------------------------- */ 1103 1104 int cholmod_allocate_work 1105 ( 1106 /* ---- input ---- */ 1107 size_t nrow, /* size: Common->Flag (nrow), Common->Head (nrow+1) */ 1108 size_t iworksize, /* size of Common->Iwork */ 1109 size_t xworksize, /* size of Common->Xwork */ 1110 /* --------------- */ 1111 cholmod_common *Common 1112 ) ; 1113 1114 int cholmod_l_allocate_work (size_t, size_t, size_t, cholmod_common *) ; 1115 1116 /* -------------------------------------------------------------------------- */ 1117 /* cholmod_free_work: free workspace in Common */ 1118 /* -------------------------------------------------------------------------- */ 1119 1120 int cholmod_free_work 1121 ( 1122 cholmod_common *Common 1123 ) ; 1124 1125 int cholmod_l_free_work (cholmod_common *) ; 1126 1127 /* -------------------------------------------------------------------------- */ 1128 /* cholmod_clear_flag: clear Flag workspace in Common */ 1129 /* -------------------------------------------------------------------------- */ 1130 1131 /* use a macro for speed */ 1132 #define CHOLMOD_CLEAR_FLAG(Common) \ 1133 { \ 1134 Common->mark++ ; \ 1135 if (Common->mark <= 0) \ 1136 { \ 1137 Common->mark = EMPTY ; \ 1138 CHOLMOD (clear_flag) (Common) ; \ 1139 } \ 1140 } 1141 1142 SuiteSparse_long cholmod_clear_flag 1143 ( 1144 cholmod_common *Common 1145 ) ; 1146 1147 SuiteSparse_long cholmod_l_clear_flag (cholmod_common *) ; 1148 1149 /* -------------------------------------------------------------------------- */ 1150 /* cholmod_error: called when CHOLMOD encounters an error */ 1151 /* -------------------------------------------------------------------------- */ 1152 1153 int cholmod_error 1154 ( 1155 /* ---- input ---- */ 1156 int status, /* error status */ 1157 const char *file, /* name of source code file where error occured */ 1158 int line, /* line number in source code file where error occured*/ 1159 const char *message,/* error message */ 1160 /* --------------- */ 1161 cholmod_common *Common 1162 ) ; 1163 1164 int cholmod_l_error (int, const char *, int, const char *, cholmod_common *) ; 1165 1166 /* -------------------------------------------------------------------------- */ 1167 /* cholmod_dbound: for internal use in CHOLMOD only */ 1168 /* -------------------------------------------------------------------------- */ 1169 1170 double cholmod_dbound /* returns modified diagonal entry of D or L */ 1171 ( 1172 /* ---- input ---- */ 1173 double dj, /* diagonal entry of D for LDL' or L for LL' */ 1174 /* --------------- */ 1175 cholmod_common *Common 1176 ) ; 1177 1178 double cholmod_l_dbound (double, cholmod_common *) ; 1179 1180 /* -------------------------------------------------------------------------- */ 1181 /* cholmod_hypot: compute sqrt (x*x + y*y) accurately */ 1182 /* -------------------------------------------------------------------------- */ 1183 1184 double cholmod_hypot 1185 ( 1186 /* ---- input ---- */ 1187 double x, double y 1188 ) ; 1189 1190 double cholmod_l_hypot (double, double) ; 1191 1192 /* -------------------------------------------------------------------------- */ 1193 /* cholmod_divcomplex: complex division, c = a/b */ 1194 /* -------------------------------------------------------------------------- */ 1195 1196 int cholmod_divcomplex /* return 1 if divide-by-zero, 0 otherise */ 1197 ( 1198 /* ---- input ---- */ 1199 double ar, double ai, /* real and imaginary parts of a */ 1200 double br, double bi, /* real and imaginary parts of b */ 1201 /* ---- output --- */ 1202 double *cr, double *ci /* real and imaginary parts of c */ 1203 ) ; 1204 1205 int cholmod_l_divcomplex (double, double, double, double, double *, double *) ; 1206 1207 1208 /* ========================================================================== */ 1209 /* === Core/cholmod_sparse ================================================== */ 1210 /* ========================================================================== */ 1211 1212 /* A sparse matrix stored in compressed-column form. */ 1213 1214 typedef struct cholmod_sparse_struct 1215 { 1216 size_t nrow ; /* the matrix is nrow-by-ncol */ 1217 size_t ncol ; 1218 size_t nzmax ; /* maximum number of entries in the matrix */ 1219 1220 /* pointers to int or SuiteSparse_long: */ 1221 void *p ; /* p [0..ncol], the column pointers */ 1222 void *i ; /* i [0..nzmax-1], the row indices */ 1223 1224 /* for unpacked matrices only: */ 1225 void *nz ; /* nz [0..ncol-1], the # of nonzeros in each col. In 1226 * packed form, the nonzero pattern of column j is in 1227 * A->i [A->p [j] ... A->p [j+1]-1]. In unpacked form, column j is in 1228 * A->i [A->p [j] ... A->p [j]+A->nz[j]-1] instead. In both cases, the 1229 * numerical values (if present) are in the corresponding locations in 1230 * the array x (or z if A->xtype is CHOLMOD_ZOMPLEX). */ 1231 1232 /* pointers to double or float: */ 1233 void *x ; /* size nzmax or 2*nzmax, if present */ 1234 void *z ; /* size nzmax, if present */ 1235 1236 int stype ; /* Describes what parts of the matrix are considered: 1237 * 1238 * 0: matrix is "unsymmetric": use both upper and lower triangular parts 1239 * (the matrix may actually be symmetric in pattern and value, but 1240 * both parts are explicitly stored and used). May be square or 1241 * rectangular. 1242 * >0: matrix is square and symmetric, use upper triangular part. 1243 * Entries in the lower triangular part are ignored. 1244 * <0: matrix is square and symmetric, use lower triangular part. 1245 * Entries in the upper triangular part are ignored. 1246 * 1247 * Note that stype>0 and stype<0 are different for cholmod_sparse and 1248 * cholmod_triplet. See the cholmod_triplet data structure for more 1249 * details. 1250 */ 1251 1252 int itype ; /* CHOLMOD_INT: p, i, and nz are int. 1253 * CHOLMOD_INTLONG: p is SuiteSparse_long, 1254 * i and nz are int. 1255 * CHOLMOD_LONG: p, i, and nz are SuiteSparse_long */ 1256 1257 int xtype ; /* pattern, real, complex, or zomplex */ 1258 int dtype ; /* x and z are double or float */ 1259 int sorted ; /* TRUE if columns are sorted, FALSE otherwise */ 1260 int packed ; /* TRUE if packed (nz ignored), FALSE if unpacked 1261 * (nz is required) */ 1262 1263 } cholmod_sparse ; 1264 1265 typedef struct cholmod_descendant_score_t { 1266 double score; 1267 SuiteSparse_long d; 1268 } descendantScore; 1269 1270 /* For sorting descendant supernodes with qsort */ 1271 int cholmod_score_comp (struct cholmod_descendant_score_t *i, 1272 struct cholmod_descendant_score_t *j); 1273 1274 int cholmod_l_score_comp (struct cholmod_descendant_score_t *i, 1275 struct cholmod_descendant_score_t *j); 1276 1277 /* -------------------------------------------------------------------------- */ 1278 /* cholmod_allocate_sparse: allocate a sparse matrix */ 1279 /* -------------------------------------------------------------------------- */ 1280 1281 cholmod_sparse *cholmod_allocate_sparse 1282 ( 1283 /* ---- input ---- */ 1284 size_t nrow, /* # of rows of A */ 1285 size_t ncol, /* # of columns of A */ 1286 size_t nzmax, /* max # of nonzeros of A */ 1287 int sorted, /* TRUE if columns of A sorted, FALSE otherwise */ 1288 int packed, /* TRUE if A will be packed, FALSE otherwise */ 1289 int stype, /* stype of A */ 1290 int xtype, /* CHOLMOD_PATTERN, _REAL, _COMPLEX, or _ZOMPLEX */ 1291 /* --------------- */ 1292 cholmod_common *Common 1293 ) ; 1294 1295 cholmod_sparse *cholmod_l_allocate_sparse (size_t, size_t, size_t, int, int, 1296 int, int, cholmod_common *) ; 1297 1298 /* -------------------------------------------------------------------------- */ 1299 /* cholmod_free_sparse: free a sparse matrix */ 1300 /* -------------------------------------------------------------------------- */ 1301 1302 int cholmod_free_sparse 1303 ( 1304 /* ---- in/out --- */ 1305 cholmod_sparse **A, /* matrix to deallocate, NULL on output */ 1306 /* --------------- */ 1307 cholmod_common *Common 1308 ) ; 1309 1310 int cholmod_l_free_sparse (cholmod_sparse **, cholmod_common *) ; 1311 1312 /* -------------------------------------------------------------------------- */ 1313 /* cholmod_reallocate_sparse: change the size (# entries) of sparse matrix */ 1314 /* -------------------------------------------------------------------------- */ 1315 1316 int cholmod_reallocate_sparse 1317 ( 1318 /* ---- input ---- */ 1319 size_t nznew, /* new # of entries in A */ 1320 /* ---- in/out --- */ 1321 cholmod_sparse *A, /* matrix to reallocate */ 1322 /* --------------- */ 1323 cholmod_common *Common 1324 ) ; 1325 1326 int cholmod_l_reallocate_sparse ( size_t, cholmod_sparse *, cholmod_common *) ; 1327 1328 /* -------------------------------------------------------------------------- */ 1329 /* cholmod_nnz: return number of nonzeros in a sparse matrix */ 1330 /* -------------------------------------------------------------------------- */ 1331 1332 SuiteSparse_long cholmod_nnz 1333 ( 1334 /* ---- input ---- */ 1335 cholmod_sparse *A, 1336 /* --------------- */ 1337 cholmod_common *Common 1338 ) ; 1339 1340 SuiteSparse_long cholmod_l_nnz (cholmod_sparse *, cholmod_common *) ; 1341 1342 /* -------------------------------------------------------------------------- */ 1343 /* cholmod_speye: sparse identity matrix */ 1344 /* -------------------------------------------------------------------------- */ 1345 1346 cholmod_sparse *cholmod_speye 1347 ( 1348 /* ---- input ---- */ 1349 size_t nrow, /* # of rows of A */ 1350 size_t ncol, /* # of columns of A */ 1351 int xtype, /* CHOLMOD_PATTERN, _REAL, _COMPLEX, or _ZOMPLEX */ 1352 /* --------------- */ 1353 cholmod_common *Common 1354 ) ; 1355 1356 cholmod_sparse *cholmod_l_speye (size_t, size_t, int, cholmod_common *) ; 1357 1358 /* -------------------------------------------------------------------------- */ 1359 /* cholmod_spzeros: sparse zero matrix */ 1360 /* -------------------------------------------------------------------------- */ 1361 1362 cholmod_sparse *cholmod_spzeros 1363 ( 1364 /* ---- input ---- */ 1365 size_t nrow, /* # of rows of A */ 1366 size_t ncol, /* # of columns of A */ 1367 size_t nzmax, /* max # of nonzeros of A */ 1368 int xtype, /* CHOLMOD_PATTERN, _REAL, _COMPLEX, or _ZOMPLEX */ 1369 /* --------------- */ 1370 cholmod_common *Common 1371 ) ; 1372 1373 cholmod_sparse *cholmod_l_spzeros (size_t, size_t, size_t, int, 1374 cholmod_common *) ; 1375 1376 /* -------------------------------------------------------------------------- */ 1377 /* cholmod_transpose: transpose a sparse matrix */ 1378 /* -------------------------------------------------------------------------- */ 1379 1380 /* Return A' or A.' The "values" parameter is 0, 1, or 2 to denote the pattern 1381 * transpose, the array transpose (A.'), and the complex conjugate transpose 1382 * (A'). 1383 */ 1384 1385 cholmod_sparse *cholmod_transpose 1386 ( 1387 /* ---- input ---- */ 1388 cholmod_sparse *A, /* matrix to transpose */ 1389 int values, /* 0: pattern, 1: array transpose, 2: conj. transpose */ 1390 /* --------------- */ 1391 cholmod_common *Common 1392 ) ; 1393 1394 cholmod_sparse *cholmod_l_transpose (cholmod_sparse *, int, cholmod_common *) ; 1395 1396 /* -------------------------------------------------------------------------- */ 1397 /* cholmod_transpose_unsym: transpose an unsymmetric sparse matrix */ 1398 /* -------------------------------------------------------------------------- */ 1399 1400 /* Compute F = A', A (:,f)', or A (p,f)', where A is unsymmetric and F is 1401 * already allocated. See cholmod_transpose for a simpler routine. */ 1402 1403 int cholmod_transpose_unsym 1404 ( 1405 /* ---- input ---- */ 1406 cholmod_sparse *A, /* matrix to transpose */ 1407 int values, /* 0: pattern, 1: array transpose, 2: conj. transpose */ 1408 int *Perm, /* size nrow, if present (can be NULL) */ 1409 int *fset, /* subset of 0:(A->ncol)-1 */ 1410 size_t fsize, /* size of fset */ 1411 /* ---- output --- */ 1412 cholmod_sparse *F, /* F = A', A(:,f)', or A(p,f)' */ 1413 /* --------------- */ 1414 cholmod_common *Common 1415 ) ; 1416 1417 int cholmod_l_transpose_unsym (cholmod_sparse *, int, SuiteSparse_long *, 1418 SuiteSparse_long *, size_t, cholmod_sparse *, cholmod_common *) ; 1419 1420 /* -------------------------------------------------------------------------- */ 1421 /* cholmod_transpose_sym: transpose a symmetric sparse matrix */ 1422 /* -------------------------------------------------------------------------- */ 1423 1424 /* Compute F = A' or A (p,p)', where A is symmetric and F is already allocated. 1425 * See cholmod_transpose for a simpler routine. */ 1426 1427 int cholmod_transpose_sym 1428 ( 1429 /* ---- input ---- */ 1430 cholmod_sparse *A, /* matrix to transpose */ 1431 int values, /* 0: pattern, 1: array transpose, 2: conj. transpose */ 1432 int *Perm, /* size nrow, if present (can be NULL) */ 1433 /* ---- output --- */ 1434 cholmod_sparse *F, /* F = A' or A(p,p)' */ 1435 /* --------------- */ 1436 cholmod_common *Common 1437 ) ; 1438 1439 int cholmod_l_transpose_sym (cholmod_sparse *, int, SuiteSparse_long *, 1440 cholmod_sparse *, cholmod_common *) ; 1441 1442 /* -------------------------------------------------------------------------- */ 1443 /* cholmod_ptranspose: transpose a sparse matrix */ 1444 /* -------------------------------------------------------------------------- */ 1445 1446 /* Return A' or A(p,p)' if A is symmetric. Return A', A(:,f)', or A(p,f)' if 1447 * A is unsymmetric. */ 1448 1449 cholmod_sparse *cholmod_ptranspose 1450 ( 1451 /* ---- input ---- */ 1452 cholmod_sparse *A, /* matrix to transpose */ 1453 int values, /* 0: pattern, 1: array transpose, 2: conj. transpose */ 1454 int *Perm, /* if non-NULL, F = A(p,f) or A(p,p) */ 1455 int *fset, /* subset of 0:(A->ncol)-1 */ 1456 size_t fsize, /* size of fset */ 1457 /* --------------- */ 1458 cholmod_common *Common 1459 ) ; 1460 1461 cholmod_sparse *cholmod_l_ptranspose (cholmod_sparse *, int, SuiteSparse_long *, 1462 SuiteSparse_long *, size_t, cholmod_common *) ; 1463 1464 /* -------------------------------------------------------------------------- */ 1465 /* cholmod_sort: sort row indices in each column of sparse matrix */ 1466 /* -------------------------------------------------------------------------- */ 1467 1468 int cholmod_sort 1469 ( 1470 /* ---- in/out --- */ 1471 cholmod_sparse *A, /* matrix to sort */ 1472 /* --------------- */ 1473 cholmod_common *Common 1474 ) ; 1475 1476 int cholmod_l_sort (cholmod_sparse *, cholmod_common *) ; 1477 1478 /* -------------------------------------------------------------------------- */ 1479 /* cholmod_band: C = tril (triu (A,k1), k2) */ 1480 /* -------------------------------------------------------------------------- */ 1481 1482 cholmod_sparse *cholmod_band 1483 ( 1484 /* ---- input ---- */ 1485 cholmod_sparse *A, /* matrix to extract band matrix from */ 1486 SuiteSparse_long k1, /* ignore entries below the k1-st diagonal */ 1487 SuiteSparse_long k2, /* ignore entries above the k2-nd diagonal */ 1488 int mode, /* >0: numerical, 0: pattern, <0: pattern (no diag) */ 1489 /* --------------- */ 1490 cholmod_common *Common 1491 ) ; 1492 1493 cholmod_sparse *cholmod_l_band (cholmod_sparse *, SuiteSparse_long, 1494 SuiteSparse_long, int, cholmod_common *) ; 1495 1496 /* -------------------------------------------------------------------------- */ 1497 /* cholmod_band_inplace: A = tril (triu (A,k1), k2) */ 1498 /* -------------------------------------------------------------------------- */ 1499 1500 int cholmod_band_inplace 1501 ( 1502 /* ---- input ---- */ 1503 SuiteSparse_long k1, /* ignore entries below the k1-st diagonal */ 1504 SuiteSparse_long k2, /* ignore entries above the k2-nd diagonal */ 1505 int mode, /* >0: numerical, 0: pattern, <0: pattern (no diag) */ 1506 /* ---- in/out --- */ 1507 cholmod_sparse *A, /* matrix from which entries not in band are removed */ 1508 /* --------------- */ 1509 cholmod_common *Common 1510 ) ; 1511 1512 int cholmod_l_band_inplace (SuiteSparse_long, SuiteSparse_long, int, 1513 cholmod_sparse *, cholmod_common *) ; 1514 1515 /* -------------------------------------------------------------------------- */ 1516 /* cholmod_aat: C = A*A' or A(:,f)*A(:,f)' */ 1517 /* -------------------------------------------------------------------------- */ 1518 1519 cholmod_sparse *cholmod_aat 1520 ( 1521 /* ---- input ---- */ 1522 cholmod_sparse *A, /* input matrix; C=A*A' is constructed */ 1523 int *fset, /* subset of 0:(A->ncol)-1 */ 1524 size_t fsize, /* size of fset */ 1525 int mode, /* >0: numerical, 0: pattern, <0: pattern (no diag), 1526 * -2: pattern only, no diagonal, add 50%+n extra 1527 * space to C */ 1528 /* --------------- */ 1529 cholmod_common *Common 1530 ) ; 1531 1532 cholmod_sparse *cholmod_l_aat (cholmod_sparse *, SuiteSparse_long *, size_t, 1533 int, cholmod_common *) ; 1534 1535 /* -------------------------------------------------------------------------- */ 1536 /* cholmod_copy_sparse: C = A, create an exact copy of a sparse matrix */ 1537 /* -------------------------------------------------------------------------- */ 1538 1539 cholmod_sparse *cholmod_copy_sparse 1540 ( 1541 /* ---- input ---- */ 1542 cholmod_sparse *A, /* matrix to copy */ 1543 /* --------------- */ 1544 cholmod_common *Common 1545 ) ; 1546 1547 cholmod_sparse *cholmod_l_copy_sparse (cholmod_sparse *, cholmod_common *) ; 1548 1549 /* -------------------------------------------------------------------------- */ 1550 /* cholmod_copy: C = A, with possible change of stype */ 1551 /* -------------------------------------------------------------------------- */ 1552 1553 cholmod_sparse *cholmod_copy 1554 ( 1555 /* ---- input ---- */ 1556 cholmod_sparse *A, /* matrix to copy */ 1557 int stype, /* requested stype of C */ 1558 int mode, /* >0: numerical, 0: pattern, <0: pattern (no diag) */ 1559 /* --------------- */ 1560 cholmod_common *Common 1561 ) ; 1562 1563 cholmod_sparse *cholmod_l_copy (cholmod_sparse *, int, int, cholmod_common *) ; 1564 1565 /* -------------------------------------------------------------------------- */ 1566 /* cholmod_add: C = alpha*A + beta*B */ 1567 /* -------------------------------------------------------------------------- */ 1568 1569 cholmod_sparse *cholmod_add 1570 ( 1571 /* ---- input ---- */ 1572 cholmod_sparse *A, /* matrix to add */ 1573 cholmod_sparse *B, /* matrix to add */ 1574 double alpha [2], /* scale factor for A */ 1575 double beta [2], /* scale factor for B */ 1576 int values, /* if TRUE compute the numerical values of C */ 1577 int sorted, /* if TRUE, sort columns of C */ 1578 /* --------------- */ 1579 cholmod_common *Common 1580 ) ; 1581 1582 cholmod_sparse *cholmod_l_add (cholmod_sparse *, cholmod_sparse *, double *, 1583 double *, int, int, cholmod_common *) ; 1584 1585 /* -------------------------------------------------------------------------- */ 1586 /* cholmod_sparse_xtype: change the xtype of a sparse matrix */ 1587 /* -------------------------------------------------------------------------- */ 1588 1589 int cholmod_sparse_xtype 1590 ( 1591 /* ---- input ---- */ 1592 int to_xtype, /* requested xtype (pattern, real, complex, zomplex) */ 1593 /* ---- in/out --- */ 1594 cholmod_sparse *A, /* sparse matrix to change */ 1595 /* --------------- */ 1596 cholmod_common *Common 1597 ) ; 1598 1599 int cholmod_l_sparse_xtype (int, cholmod_sparse *, cholmod_common *) ; 1600 1601 1602 /* ========================================================================== */ 1603 /* === Core/cholmod_factor ================================================== */ 1604 /* ========================================================================== */ 1605 1606 /* A symbolic and numeric factorization, either simplicial or supernodal. 1607 * In all cases, the row indices in the columns of L are kept sorted. */ 1608 1609 typedef struct cholmod_factor_struct 1610 { 1611 /* ---------------------------------------------------------------------- */ 1612 /* for both simplicial and supernodal factorizations */ 1613 /* ---------------------------------------------------------------------- */ 1614 1615 size_t n ; /* L is n-by-n */ 1616 1617 size_t minor ; /* If the factorization failed, L->minor is the column 1618 * at which it failed (in the range 0 to n-1). A value 1619 * of n means the factorization was successful or 1620 * the matrix has not yet been factorized. */ 1621 1622 /* ---------------------------------------------------------------------- */ 1623 /* symbolic ordering and analysis */ 1624 /* ---------------------------------------------------------------------- */ 1625 1626 void *Perm ; /* size n, permutation used */ 1627 void *ColCount ; /* size n, column counts for simplicial L */ 1628 1629 void *IPerm ; /* size n, inverse permutation. Only created by 1630 * cholmod_solve2 if Bset is used. */ 1631 1632 /* ---------------------------------------------------------------------- */ 1633 /* simplicial factorization */ 1634 /* ---------------------------------------------------------------------- */ 1635 1636 size_t nzmax ; /* size of i and x */ 1637 1638 void *p ; /* p [0..ncol], the column pointers */ 1639 void *i ; /* i [0..nzmax-1], the row indices */ 1640 void *x ; /* x [0..nzmax-1], the numerical values */ 1641 void *z ; 1642 void *nz ; /* nz [0..ncol-1], the # of nonzeros in each column. 1643 * i [p [j] ... p [j]+nz[j]-1] contains the row indices, 1644 * and the numerical values are in the same locatins 1645 * in x. The value of i [p [k]] is always k. */ 1646 1647 void *next ; /* size ncol+2. next [j] is the next column in i/x */ 1648 void *prev ; /* size ncol+2. prev [j] is the prior column in i/x. 1649 * head of the list is ncol+1, and the tail is ncol. */ 1650 1651 /* ---------------------------------------------------------------------- */ 1652 /* supernodal factorization */ 1653 /* ---------------------------------------------------------------------- */ 1654 1655 /* Note that L->x is shared with the simplicial data structure. L->x has 1656 * size L->nzmax for a simplicial factor, and size L->xsize for a supernodal 1657 * factor. */ 1658 1659 size_t nsuper ; /* number of supernodes */ 1660 size_t ssize ; /* size of s, integer part of supernodes */ 1661 size_t xsize ; /* size of x, real part of supernodes */ 1662 size_t maxcsize ; /* size of largest update matrix */ 1663 size_t maxesize ; /* max # of rows in supernodes, excl. triangular part */ 1664 1665 void *super ; /* size nsuper+1, first col in each supernode */ 1666 void *pi ; /* size nsuper+1, pointers to integer patterns */ 1667 void *px ; /* size nsuper+1, pointers to real parts */ 1668 void *s ; /* size ssize, integer part of supernodes */ 1669 1670 /* ---------------------------------------------------------------------- */ 1671 /* factorization type */ 1672 /* ---------------------------------------------------------------------- */ 1673 1674 int ordering ; /* ordering method used */ 1675 1676 int is_ll ; /* TRUE if LL', FALSE if LDL' */ 1677 int is_super ; /* TRUE if supernodal, FALSE if simplicial */ 1678 int is_monotonic ; /* TRUE if columns of L appear in order 0..n-1. 1679 * Only applicable to simplicial numeric types. */ 1680 1681 /* There are 8 types of factor objects that cholmod_factor can represent 1682 * (only 6 are used): 1683 * 1684 * Numeric types (xtype is not CHOLMOD_PATTERN) 1685 * -------------------------------------------- 1686 * 1687 * simplicial LDL': (is_ll FALSE, is_super FALSE). Stored in compressed 1688 * column form, using the simplicial components above (nzmax, p, i, 1689 * x, z, nz, next, and prev). The unit diagonal of L is not stored, 1690 * and D is stored in its place. There are no supernodes. 1691 * 1692 * simplicial LL': (is_ll TRUE, is_super FALSE). Uses the same storage 1693 * scheme as the simplicial LDL', except that D does not appear. 1694 * The first entry of each column of L is the diagonal entry of 1695 * that column of L. 1696 * 1697 * supernodal LDL': (is_ll FALSE, is_super TRUE). Not used. 1698 * FUTURE WORK: add support for supernodal LDL' 1699 * 1700 * supernodal LL': (is_ll TRUE, is_super TRUE). A supernodal factor, 1701 * using the supernodal components described above (nsuper, ssize, 1702 * xsize, maxcsize, maxesize, super, pi, px, s, x, and z). 1703 * 1704 * 1705 * Symbolic types (xtype is CHOLMOD_PATTERN) 1706 * ----------------------------------------- 1707 * 1708 * simplicial LDL': (is_ll FALSE, is_super FALSE). Nothing is present 1709 * except Perm and ColCount. 1710 * 1711 * simplicial LL': (is_ll TRUE, is_super FALSE). Identical to the 1712 * simplicial LDL', except for the is_ll flag. 1713 * 1714 * supernodal LDL': (is_ll FALSE, is_super TRUE). Not used. 1715 * FUTURE WORK: add support for supernodal LDL' 1716 * 1717 * supernodal LL': (is_ll TRUE, is_super TRUE). A supernodal symbolic 1718 * factorization. The simplicial symbolic information is present 1719 * (Perm and ColCount), as is all of the supernodal factorization 1720 * except for the numerical values (x and z). 1721 */ 1722 1723 int itype ; /* The integer arrays are Perm, ColCount, p, i, nz, 1724 * next, prev, super, pi, px, and s. If itype is 1725 * CHOLMOD_INT, all of these are int arrays. 1726 * CHOLMOD_INTLONG: p, pi, px are SuiteSparse_long, others int. 1727 * CHOLMOD_LONG: all integer arrays are SuiteSparse_long. */ 1728 int xtype ; /* pattern, real, complex, or zomplex */ 1729 int dtype ; /* x and z double or float */ 1730 1731 int useGPU; /* Indicates the symbolic factorization supports 1732 * GPU acceleration */ 1733 1734 } cholmod_factor ; 1735 1736 1737 /* -------------------------------------------------------------------------- */ 1738 /* cholmod_allocate_factor: allocate a factor (symbolic LL' or LDL') */ 1739 /* -------------------------------------------------------------------------- */ 1740 1741 cholmod_factor *cholmod_allocate_factor 1742 ( 1743 /* ---- input ---- */ 1744 size_t n, /* L is n-by-n */ 1745 /* --------------- */ 1746 cholmod_common *Common 1747 ) ; 1748 1749 cholmod_factor *cholmod_l_allocate_factor (size_t, cholmod_common *) ; 1750 1751 /* -------------------------------------------------------------------------- */ 1752 /* cholmod_free_factor: free a factor */ 1753 /* -------------------------------------------------------------------------- */ 1754 1755 int cholmod_free_factor 1756 ( 1757 /* ---- in/out --- */ 1758 cholmod_factor **L, /* factor to free, NULL on output */ 1759 /* --------------- */ 1760 cholmod_common *Common 1761 ) ; 1762 1763 int cholmod_l_free_factor (cholmod_factor **, cholmod_common *) ; 1764 1765 /* -------------------------------------------------------------------------- */ 1766 /* cholmod_reallocate_factor: change the # entries in a factor */ 1767 /* -------------------------------------------------------------------------- */ 1768 1769 int cholmod_reallocate_factor 1770 ( 1771 /* ---- input ---- */ 1772 size_t nznew, /* new # of entries in L */ 1773 /* ---- in/out --- */ 1774 cholmod_factor *L, /* factor to modify */ 1775 /* --------------- */ 1776 cholmod_common *Common 1777 ) ; 1778 1779 int cholmod_l_reallocate_factor (size_t, cholmod_factor *, cholmod_common *) ; 1780 1781 /* -------------------------------------------------------------------------- */ 1782 /* cholmod_change_factor: change the type of factor (e.g., LDL' to LL') */ 1783 /* -------------------------------------------------------------------------- */ 1784 1785 int cholmod_change_factor 1786 ( 1787 /* ---- input ---- */ 1788 int to_xtype, /* to CHOLMOD_PATTERN, _REAL, _COMPLEX, _ZOMPLEX */ 1789 int to_ll, /* TRUE: convert to LL', FALSE: LDL' */ 1790 int to_super, /* TRUE: convert to supernodal, FALSE: simplicial */ 1791 int to_packed, /* TRUE: pack simplicial columns, FALSE: do not pack */ 1792 int to_monotonic, /* TRUE: put simplicial columns in order, FALSE: not */ 1793 /* ---- in/out --- */ 1794 cholmod_factor *L, /* factor to modify */ 1795 /* --------------- */ 1796 cholmod_common *Common 1797 ) ; 1798 1799 int cholmod_l_change_factor ( int, int, int, int, int, cholmod_factor *, 1800 cholmod_common *) ; 1801 1802 /* -------------------------------------------------------------------------- */ 1803 /* cholmod_pack_factor: pack the columns of a factor */ 1804 /* -------------------------------------------------------------------------- */ 1805 1806 /* Pack the columns of a simplicial factor. Unlike cholmod_change_factor, 1807 * it can pack the columns of a factor even if they are not stored in their 1808 * natural order (non-monotonic). */ 1809 1810 int cholmod_pack_factor 1811 ( 1812 /* ---- in/out --- */ 1813 cholmod_factor *L, /* factor to modify */ 1814 /* --------------- */ 1815 cholmod_common *Common 1816 ) ; 1817 1818 int cholmod_l_pack_factor (cholmod_factor *, cholmod_common *) ; 1819 1820 /* -------------------------------------------------------------------------- */ 1821 /* cholmod_reallocate_column: resize a single column of a factor */ 1822 /* -------------------------------------------------------------------------- */ 1823 1824 int cholmod_reallocate_column 1825 ( 1826 /* ---- input ---- */ 1827 size_t j, /* the column to reallocate */ 1828 size_t need, /* required size of column j */ 1829 /* ---- in/out --- */ 1830 cholmod_factor *L, /* factor to modify */ 1831 /* --------------- */ 1832 cholmod_common *Common 1833 ) ; 1834 1835 int cholmod_l_reallocate_column (size_t, size_t, cholmod_factor *, 1836 cholmod_common *) ; 1837 1838 /* -------------------------------------------------------------------------- */ 1839 /* cholmod_factor_to_sparse: create a sparse matrix copy of a factor */ 1840 /* -------------------------------------------------------------------------- */ 1841 1842 /* Only operates on numeric factors, not symbolic ones */ 1843 1844 cholmod_sparse *cholmod_factor_to_sparse 1845 ( 1846 /* ---- in/out --- */ 1847 cholmod_factor *L, /* factor to copy, converted to symbolic on output */ 1848 /* --------------- */ 1849 cholmod_common *Common 1850 ) ; 1851 1852 cholmod_sparse *cholmod_l_factor_to_sparse (cholmod_factor *, 1853 cholmod_common *) ; 1854 1855 /* -------------------------------------------------------------------------- */ 1856 /* cholmod_copy_factor: create a copy of a factor */ 1857 /* -------------------------------------------------------------------------- */ 1858 1859 cholmod_factor *cholmod_copy_factor 1860 ( 1861 /* ---- input ---- */ 1862 cholmod_factor *L, /* factor to copy */ 1863 /* --------------- */ 1864 cholmod_common *Common 1865 ) ; 1866 1867 cholmod_factor *cholmod_l_copy_factor (cholmod_factor *, cholmod_common *) ; 1868 1869 /* -------------------------------------------------------------------------- */ 1870 /* cholmod_factor_xtype: change the xtype of a factor */ 1871 /* -------------------------------------------------------------------------- */ 1872 1873 int cholmod_factor_xtype 1874 ( 1875 /* ---- input ---- */ 1876 int to_xtype, /* requested xtype (real, complex, or zomplex) */ 1877 /* ---- in/out --- */ 1878 cholmod_factor *L, /* factor to change */ 1879 /* --------------- */ 1880 cholmod_common *Common 1881 ) ; 1882 1883 int cholmod_l_factor_xtype (int, cholmod_factor *, cholmod_common *) ; 1884 1885 1886 /* ========================================================================== */ 1887 /* === Core/cholmod_dense =================================================== */ 1888 /* ========================================================================== */ 1889 1890 /* A dense matrix in column-oriented form. It has no itype since it contains 1891 * no integers. Entry in row i and column j is located in x [i+j*d]. 1892 */ 1893 1894 typedef struct cholmod_dense_struct 1895 { 1896 size_t nrow ; /* the matrix is nrow-by-ncol */ 1897 size_t ncol ; 1898 size_t nzmax ; /* maximum number of entries in the matrix */ 1899 size_t d ; /* leading dimension (d >= nrow must hold) */ 1900 void *x ; /* size nzmax or 2*nzmax, if present */ 1901 void *z ; /* size nzmax, if present */ 1902 int xtype ; /* pattern, real, complex, or zomplex */ 1903 int dtype ; /* x and z double or float */ 1904 1905 } cholmod_dense ; 1906 1907 /* -------------------------------------------------------------------------- */ 1908 /* cholmod_allocate_dense: allocate a dense matrix (contents uninitialized) */ 1909 /* -------------------------------------------------------------------------- */ 1910 1911 cholmod_dense *cholmod_allocate_dense 1912 ( 1913 /* ---- input ---- */ 1914 size_t nrow, /* # of rows of matrix */ 1915 size_t ncol, /* # of columns of matrix */ 1916 size_t d, /* leading dimension */ 1917 int xtype, /* CHOLMOD_REAL, _COMPLEX, or _ZOMPLEX */ 1918 /* --------------- */ 1919 cholmod_common *Common 1920 ) ; 1921 1922 cholmod_dense *cholmod_l_allocate_dense (size_t, size_t, size_t, int, 1923 cholmod_common *) ; 1924 1925 /* -------------------------------------------------------------------------- */ 1926 /* cholmod_zeros: allocate a dense matrix and set it to zero */ 1927 /* -------------------------------------------------------------------------- */ 1928 1929 cholmod_dense *cholmod_zeros 1930 ( 1931 /* ---- input ---- */ 1932 size_t nrow, /* # of rows of matrix */ 1933 size_t ncol, /* # of columns of matrix */ 1934 int xtype, /* CHOLMOD_REAL, _COMPLEX, or _ZOMPLEX */ 1935 /* --------------- */ 1936 cholmod_common *Common 1937 ) ; 1938 1939 cholmod_dense *cholmod_l_zeros (size_t, size_t, int, cholmod_common *) ; 1940 1941 /* -------------------------------------------------------------------------- */ 1942 /* cholmod_ones: allocate a dense matrix and set it to all ones */ 1943 /* -------------------------------------------------------------------------- */ 1944 1945 cholmod_dense *cholmod_ones 1946 ( 1947 /* ---- input ---- */ 1948 size_t nrow, /* # of rows of matrix */ 1949 size_t ncol, /* # of columns of matrix */ 1950 int xtype, /* CHOLMOD_REAL, _COMPLEX, or _ZOMPLEX */ 1951 /* --------------- */ 1952 cholmod_common *Common 1953 ) ; 1954 1955 cholmod_dense *cholmod_l_ones (size_t, size_t, int, cholmod_common *) ; 1956 1957 /* -------------------------------------------------------------------------- */ 1958 /* cholmod_eye: allocate a dense matrix and set it to the identity matrix */ 1959 /* -------------------------------------------------------------------------- */ 1960 1961 cholmod_dense *cholmod_eye 1962 ( 1963 /* ---- input ---- */ 1964 size_t nrow, /* # of rows of matrix */ 1965 size_t ncol, /* # of columns of matrix */ 1966 int xtype, /* CHOLMOD_REAL, _COMPLEX, or _ZOMPLEX */ 1967 /* --------------- */ 1968 cholmod_common *Common 1969 ) ; 1970 1971 cholmod_dense *cholmod_l_eye (size_t, size_t, int, cholmod_common *) ; 1972 1973 /* -------------------------------------------------------------------------- */ 1974 /* cholmod_free_dense: free a dense matrix */ 1975 /* -------------------------------------------------------------------------- */ 1976 1977 int cholmod_free_dense 1978 ( 1979 /* ---- in/out --- */ 1980 cholmod_dense **X, /* dense matrix to deallocate, NULL on output */ 1981 /* --------------- */ 1982 cholmod_common *Common 1983 ) ; 1984 1985 int cholmod_l_free_dense (cholmod_dense **, cholmod_common *) ; 1986 1987 /* -------------------------------------------------------------------------- */ 1988 /* cholmod_ensure_dense: ensure a dense matrix has a given size and type */ 1989 /* -------------------------------------------------------------------------- */ 1990 1991 cholmod_dense *cholmod_ensure_dense 1992 ( 1993 /* ---- input/output ---- */ 1994 cholmod_dense **XHandle, /* matrix handle to check */ 1995 /* ---- input ---- */ 1996 size_t nrow, /* # of rows of matrix */ 1997 size_t ncol, /* # of columns of matrix */ 1998 size_t d, /* leading dimension */ 1999 int xtype, /* CHOLMOD_REAL, _COMPLEX, or _ZOMPLEX */ 2000 /* --------------- */ 2001 cholmod_common *Common 2002 ) ; 2003 2004 cholmod_dense *cholmod_l_ensure_dense (cholmod_dense **, size_t, size_t, size_t, 2005 int, cholmod_common *) ; 2006 2007 /* -------------------------------------------------------------------------- */ 2008 /* cholmod_sparse_to_dense: create a dense matrix copy of a sparse matrix */ 2009 /* -------------------------------------------------------------------------- */ 2010 2011 cholmod_dense *cholmod_sparse_to_dense 2012 ( 2013 /* ---- input ---- */ 2014 cholmod_sparse *A, /* matrix to copy */ 2015 /* --------------- */ 2016 cholmod_common *Common 2017 ) ; 2018 2019 cholmod_dense *cholmod_l_sparse_to_dense (cholmod_sparse *, 2020 cholmod_common *) ; 2021 2022 /* -------------------------------------------------------------------------- */ 2023 /* cholmod_dense_to_sparse: create a sparse matrix copy of a dense matrix */ 2024 /* -------------------------------------------------------------------------- */ 2025 2026 cholmod_sparse *cholmod_dense_to_sparse 2027 ( 2028 /* ---- input ---- */ 2029 cholmod_dense *X, /* matrix to copy */ 2030 int values, /* TRUE if values to be copied, FALSE otherwise */ 2031 /* --------------- */ 2032 cholmod_common *Common 2033 ) ; 2034 2035 cholmod_sparse *cholmod_l_dense_to_sparse (cholmod_dense *, int, 2036 cholmod_common *) ; 2037 2038 /* -------------------------------------------------------------------------- */ 2039 /* cholmod_copy_dense: create a copy of a dense matrix */ 2040 /* -------------------------------------------------------------------------- */ 2041 2042 cholmod_dense *cholmod_copy_dense 2043 ( 2044 /* ---- input ---- */ 2045 cholmod_dense *X, /* matrix to copy */ 2046 /* --------------- */ 2047 cholmod_common *Common 2048 ) ; 2049 2050 cholmod_dense *cholmod_l_copy_dense (cholmod_dense *, cholmod_common *) ; 2051 2052 /* -------------------------------------------------------------------------- */ 2053 /* cholmod_copy_dense2: copy a dense matrix (pre-allocated) */ 2054 /* -------------------------------------------------------------------------- */ 2055 2056 int cholmod_copy_dense2 2057 ( 2058 /* ---- input ---- */ 2059 cholmod_dense *X, /* matrix to copy */ 2060 /* ---- output --- */ 2061 cholmod_dense *Y, /* copy of matrix X */ 2062 /* --------------- */ 2063 cholmod_common *Common 2064 ) ; 2065 2066 int cholmod_l_copy_dense2 (cholmod_dense *, cholmod_dense *, cholmod_common *) ; 2067 2068 /* -------------------------------------------------------------------------- */ 2069 /* cholmod_dense_xtype: change the xtype of a dense matrix */ 2070 /* -------------------------------------------------------------------------- */ 2071 2072 int cholmod_dense_xtype 2073 ( 2074 /* ---- input ---- */ 2075 int to_xtype, /* requested xtype (real, complex,or zomplex) */ 2076 /* ---- in/out --- */ 2077 cholmod_dense *X, /* dense matrix to change */ 2078 /* --------------- */ 2079 cholmod_common *Common 2080 ) ; 2081 2082 int cholmod_l_dense_xtype (int, cholmod_dense *, cholmod_common *) ; 2083 2084 2085 /* ========================================================================== */ 2086 /* === Core/cholmod_triplet ================================================= */ 2087 /* ========================================================================== */ 2088 2089 /* A sparse matrix stored in triplet form. */ 2090 2091 typedef struct cholmod_triplet_struct 2092 { 2093 size_t nrow ; /* the matrix is nrow-by-ncol */ 2094 size_t ncol ; 2095 size_t nzmax ; /* maximum number of entries in the matrix */ 2096 size_t nnz ; /* number of nonzeros in the matrix */ 2097 2098 void *i ; /* i [0..nzmax-1], the row indices */ 2099 void *j ; /* j [0..nzmax-1], the column indices */ 2100 void *x ; /* size nzmax or 2*nzmax, if present */ 2101 void *z ; /* size nzmax, if present */ 2102 2103 int stype ; /* Describes what parts of the matrix are considered: 2104 * 2105 * 0: matrix is "unsymmetric": use both upper and lower triangular parts 2106 * (the matrix may actually be symmetric in pattern and value, but 2107 * both parts are explicitly stored and used). May be square or 2108 * rectangular. 2109 * >0: matrix is square and symmetric. Entries in the lower triangular 2110 * part are transposed and added to the upper triangular part when 2111 * the matrix is converted to cholmod_sparse form. 2112 * <0: matrix is square and symmetric. Entries in the upper triangular 2113 * part are transposed and added to the lower triangular part when 2114 * the matrix is converted to cholmod_sparse form. 2115 * 2116 * Note that stype>0 and stype<0 are different for cholmod_sparse and 2117 * cholmod_triplet. The reason is simple. You can permute a symmetric 2118 * triplet matrix by simply replacing a row and column index with their 2119 * new row and column indices, via an inverse permutation. Suppose 2120 * P = L->Perm is your permutation, and Pinv is an array of size n. 2121 * Suppose a symmetric matrix A is represent by a triplet matrix T, with 2122 * entries only in the upper triangular part. Then the following code: 2123 * 2124 * Ti = T->i ; 2125 * Tj = T->j ; 2126 * for (k = 0 ; k < n ; k++) Pinv [P [k]] = k ; 2127 * for (k = 0 ; k < nz ; k++) Ti [k] = Pinv [Ti [k]] ; 2128 * for (k = 0 ; k < nz ; k++) Tj [k] = Pinv [Tj [k]] ; 2129 * 2130 * creates the triplet form of C=P*A*P'. However, if T initially 2131 * contains just the upper triangular entries (T->stype = 1), after 2132 * permutation it has entries in both the upper and lower triangular 2133 * parts. These entries should be transposed when constructing the 2134 * cholmod_sparse form of A, which is what cholmod_triplet_to_sparse 2135 * does. Thus: 2136 * 2137 * C = cholmod_triplet_to_sparse (T, 0, &Common) ; 2138 * 2139 * will return the matrix C = P*A*P'. 2140 * 2141 * Since the triplet matrix T is so simple to generate, it's quite easy 2142 * to remove entries that you do not want, prior to converting T to the 2143 * cholmod_sparse form. So if you include these entries in T, CHOLMOD 2144 * assumes that there must be a reason (such as the one above). Thus, 2145 * no entry in a triplet matrix is ever ignored. 2146 */ 2147 2148 int itype ; /* CHOLMOD_LONG: i and j are SuiteSparse_long. Otherwise int */ 2149 int xtype ; /* pattern, real, complex, or zomplex */ 2150 int dtype ; /* x and z are double or float */ 2151 2152 } cholmod_triplet ; 2153 2154 /* -------------------------------------------------------------------------- */ 2155 /* cholmod_allocate_triplet: allocate a triplet matrix */ 2156 /* -------------------------------------------------------------------------- */ 2157 2158 cholmod_triplet *cholmod_allocate_triplet 2159 ( 2160 /* ---- input ---- */ 2161 size_t nrow, /* # of rows of T */ 2162 size_t ncol, /* # of columns of T */ 2163 size_t nzmax, /* max # of nonzeros of T */ 2164 int stype, /* stype of T */ 2165 int xtype, /* CHOLMOD_PATTERN, _REAL, _COMPLEX, or _ZOMPLEX */ 2166 /* --------------- */ 2167 cholmod_common *Common 2168 ) ; 2169 2170 cholmod_triplet *cholmod_l_allocate_triplet (size_t, size_t, size_t, int, int, 2171 cholmod_common *) ; 2172 2173 /* -------------------------------------------------------------------------- */ 2174 /* cholmod_free_triplet: free a triplet matrix */ 2175 /* -------------------------------------------------------------------------- */ 2176 2177 int cholmod_free_triplet 2178 ( 2179 /* ---- in/out --- */ 2180 cholmod_triplet **T, /* triplet matrix to deallocate, NULL on output */ 2181 /* --------------- */ 2182 cholmod_common *Common 2183 ) ; 2184 2185 int cholmod_l_free_triplet (cholmod_triplet **, cholmod_common *) ; 2186 2187 /* -------------------------------------------------------------------------- */ 2188 /* cholmod_reallocate_triplet: change the # of entries in a triplet matrix */ 2189 /* -------------------------------------------------------------------------- */ 2190 2191 int cholmod_reallocate_triplet 2192 ( 2193 /* ---- input ---- */ 2194 size_t nznew, /* new # of entries in T */ 2195 /* ---- in/out --- */ 2196 cholmod_triplet *T, /* triplet matrix to modify */ 2197 /* --------------- */ 2198 cholmod_common *Common 2199 ) ; 2200 2201 int cholmod_l_reallocate_triplet (size_t, cholmod_triplet *, cholmod_common *) ; 2202 2203 /* -------------------------------------------------------------------------- */ 2204 /* cholmod_sparse_to_triplet: create a triplet matrix copy of a sparse matrix*/ 2205 /* -------------------------------------------------------------------------- */ 2206 2207 cholmod_triplet *cholmod_sparse_to_triplet 2208 ( 2209 /* ---- input ---- */ 2210 cholmod_sparse *A, /* matrix to copy */ 2211 /* --------------- */ 2212 cholmod_common *Common 2213 ) ; 2214 2215 cholmod_triplet *cholmod_l_sparse_to_triplet (cholmod_sparse *, 2216 cholmod_common *) ; 2217 2218 /* -------------------------------------------------------------------------- */ 2219 /* cholmod_triplet_to_sparse: create a sparse matrix copy of a triplet matrix*/ 2220 /* -------------------------------------------------------------------------- */ 2221 2222 cholmod_sparse *cholmod_triplet_to_sparse 2223 ( 2224 /* ---- input ---- */ 2225 cholmod_triplet *T, /* matrix to copy */ 2226 size_t nzmax, /* allocate at least this much space in output matrix */ 2227 /* --------------- */ 2228 cholmod_common *Common 2229 ) ; 2230 2231 cholmod_sparse *cholmod_l_triplet_to_sparse (cholmod_triplet *, size_t, 2232 cholmod_common *) ; 2233 2234 /* -------------------------------------------------------------------------- */ 2235 /* cholmod_copy_triplet: create a copy of a triplet matrix */ 2236 /* -------------------------------------------------------------------------- */ 2237 2238 cholmod_triplet *cholmod_copy_triplet 2239 ( 2240 /* ---- input ---- */ 2241 cholmod_triplet *T, /* matrix to copy */ 2242 /* --------------- */ 2243 cholmod_common *Common 2244 ) ; 2245 2246 cholmod_triplet *cholmod_l_copy_triplet (cholmod_triplet *, cholmod_common *) ; 2247 2248 /* -------------------------------------------------------------------------- */ 2249 /* cholmod_triplet_xtype: change the xtype of a triplet matrix */ 2250 /* -------------------------------------------------------------------------- */ 2251 2252 int cholmod_triplet_xtype 2253 ( 2254 /* ---- input ---- */ 2255 int to_xtype, /* requested xtype (pattern, real, complex,or zomplex)*/ 2256 /* ---- in/out --- */ 2257 cholmod_triplet *T, /* triplet matrix to change */ 2258 /* --------------- */ 2259 cholmod_common *Common 2260 ) ; 2261 2262 int cholmod_l_triplet_xtype (int, cholmod_triplet *, cholmod_common *) ; 2263 2264 2265 /* ========================================================================== */ 2266 /* === Core/cholmod_memory ================================================== */ 2267 /* ========================================================================== */ 2268 2269 /* The user may make use of these, just like malloc and free. You can even 2270 * malloc an object and safely free it with cholmod_free, and visa versa 2271 * (except that the memory usage statistics will be corrupted). These routines 2272 * do differ from malloc and free. If cholmod_free is given a NULL pointer, 2273 * for example, it does nothing (unlike the ANSI free). cholmod_realloc does 2274 * not return NULL if given a non-NULL pointer and a nonzero size, even if it 2275 * fails (it returns the original pointer and sets an error code in 2276 * Common->status instead). 2277 * 2278 * CHOLMOD keeps track of the amount of memory it has allocated, and so the 2279 * cholmod_free routine also takes the size of the object being freed. This 2280 * is only used for statistics. If you, the user of CHOLMOD, pass the wrong 2281 * size, the only consequence is that the memory usage statistics will be 2282 * corrupted. 2283 */ 2284 2285 void *cholmod_malloc /* returns pointer to the newly malloc'd block */ 2286 ( 2287 /* ---- input ---- */ 2288 size_t n, /* number of items */ 2289 size_t size, /* size of each item */ 2290 /* --------------- */ 2291 cholmod_common *Common 2292 ) ; 2293 2294 void *cholmod_l_malloc (size_t, size_t, cholmod_common *) ; 2295 2296 void *cholmod_calloc /* returns pointer to the newly calloc'd block */ 2297 ( 2298 /* ---- input ---- */ 2299 size_t n, /* number of items */ 2300 size_t size, /* size of each item */ 2301 /* --------------- */ 2302 cholmod_common *Common 2303 ) ; 2304 2305 void *cholmod_l_calloc (size_t, size_t, cholmod_common *) ; 2306 2307 void *cholmod_free /* always returns NULL */ 2308 ( 2309 /* ---- input ---- */ 2310 size_t n, /* number of items */ 2311 size_t size, /* size of each item */ 2312 /* ---- in/out --- */ 2313 void *p, /* block of memory to free */ 2314 /* --------------- */ 2315 cholmod_common *Common 2316 ) ; 2317 2318 void *cholmod_l_free (size_t, size_t, void *, cholmod_common *) ; 2319 2320 void *cholmod_realloc /* returns pointer to reallocated block */ 2321 ( 2322 /* ---- input ---- */ 2323 size_t nnew, /* requested # of items in reallocated block */ 2324 size_t size, /* size of each item */ 2325 /* ---- in/out --- */ 2326 void *p, /* block of memory to realloc */ 2327 size_t *n, /* current size on input, nnew on output if successful*/ 2328 /* --------------- */ 2329 cholmod_common *Common 2330 ) ; 2331 2332 void *cholmod_l_realloc (size_t, size_t, void *, size_t *, cholmod_common *) ; 2333 2334 int cholmod_realloc_multiple 2335 ( 2336 /* ---- input ---- */ 2337 size_t nnew, /* requested # of items in reallocated blocks */ 2338 int nint, /* number of int/SuiteSparse_long blocks */ 2339 int xtype, /* CHOLMOD_PATTERN, _REAL, _COMPLEX, or _ZOMPLEX */ 2340 /* ---- in/out --- */ 2341 void **Iblock, /* int or SuiteSparse_long block */ 2342 void **Jblock, /* int or SuiteSparse_long block */ 2343 void **Xblock, /* complex, double, or float block */ 2344 void **Zblock, /* zomplex case only: double or float block */ 2345 size_t *n, /* current size of the I,J,X,Z blocks on input, 2346 * nnew on output if successful */ 2347 /* --------------- */ 2348 cholmod_common *Common 2349 ) ; 2350 2351 int cholmod_l_realloc_multiple (size_t, int, int, void **, void **, void **, 2352 void **, size_t *, cholmod_common *) ; 2353 2354 /* ========================================================================== */ 2355 /* === version control ====================================================== */ 2356 /* ========================================================================== */ 2357 2358 int cholmod_version /* returns CHOLMOD_VERSION */ 2359 ( 2360 /* output, contents not defined on input. Not used if NULL. 2361 version [0] = CHOLMOD_MAIN_VERSION 2362 version [1] = CHOLMOD_SUB_VERSION 2363 version [2] = CHOLMOD_SUBSUB_VERSION 2364 */ 2365 int version [3] 2366 ) ; 2367 2368 int cholmod_l_version (int version [3]) ; 2369 2370 /* Versions prior to 2.1.1 do not have the above function. The following 2371 code fragment will work with any version of CHOLMOD: 2372 #ifdef CHOLMOD_HAS_VERSION_FUNCTION 2373 v = cholmod_version (NULL) ; 2374 #else 2375 v = CHOLMOD_VERSION ; 2376 #endif 2377 */ 2378 2379 /* ========================================================================== */ 2380 /* === symmetry types ======================================================= */ 2381 /* ========================================================================== */ 2382 2383 #define CHOLMOD_MM_RECTANGULAR 1 2384 #define CHOLMOD_MM_UNSYMMETRIC 2 2385 #define CHOLMOD_MM_SYMMETRIC 3 2386 #define CHOLMOD_MM_HERMITIAN 4 2387 #define CHOLMOD_MM_SKEW_SYMMETRIC 5 2388 #define CHOLMOD_MM_SYMMETRIC_POSDIAG 6 2389 #define CHOLMOD_MM_HERMITIAN_POSDIAG 7 2390 2391 /* ========================================================================== */ 2392 /* === Numerical relop macros =============================================== */ 2393 /* ========================================================================== */ 2394 2395 /* These macros correctly handle the NaN case. 2396 * 2397 * CHOLMOD_IS_NAN(x): 2398 * True if x is NaN. False otherwise. The commonly-existing isnan(x) 2399 * function could be used, but it's not in Kernighan & Ritchie 2nd edition 2400 * (ANSI C89). It may appear in <math.h>, but I'm not certain about 2401 * portability. The expression x != x is true if and only if x is NaN, 2402 * according to the IEEE 754 floating-point standard. 2403 * 2404 * CHOLMOD_IS_ZERO(x): 2405 * True if x is zero. False if x is nonzero, NaN, or +/- Inf. 2406 * This is (x == 0) if the compiler is IEEE 754 compliant. 2407 * 2408 * CHOLMOD_IS_NONZERO(x): 2409 * True if x is nonzero, NaN, or +/- Inf. False if x zero. 2410 * This is (x != 0) if the compiler is IEEE 754 compliant. 2411 * 2412 * CHOLMOD_IS_LT_ZERO(x): 2413 * True if x is < zero or -Inf. False if x is >= 0, NaN, or +Inf. 2414 * This is (x < 0) if the compiler is IEEE 754 compliant. 2415 * 2416 * CHOLMOD_IS_GT_ZERO(x): 2417 * True if x is > zero or +Inf. False if x is <= 0, NaN, or -Inf. 2418 * This is (x > 0) if the compiler is IEEE 754 compliant. 2419 * 2420 * CHOLMOD_IS_LE_ZERO(x): 2421 * True if x is <= zero or -Inf. False if x is > 0, NaN, or +Inf. 2422 * This is (x <= 0) if the compiler is IEEE 754 compliant. 2423 */ 2424 2425 #ifdef CHOLMOD_WINDOWS 2426 2427 /* Yes, this is exceedingly ugly. Blame Microsoft, which hopelessly */ 2428 /* violates the IEEE 754 floating-point standard in a bizarre way. */ 2429 /* If you're using an IEEE 754-compliant compiler, then x != x is true */ 2430 /* iff x is NaN. For Microsoft, (x < x) is true iff x is NaN. */ 2431 /* So either way, this macro safely detects a NaN. */ 2432 #define CHOLMOD_IS_NAN(x) (((x) != (x)) || (((x) < (x)))) 2433 #define CHOLMOD_IS_ZERO(x) (((x) == 0.) && !CHOLMOD_IS_NAN(x)) 2434 #define CHOLMOD_IS_NONZERO(x) (((x) != 0.) || CHOLMOD_IS_NAN(x)) 2435 #define CHOLMOD_IS_LT_ZERO(x) (((x) < 0.) && !CHOLMOD_IS_NAN(x)) 2436 #define CHOLMOD_IS_GT_ZERO(x) (((x) > 0.) && !CHOLMOD_IS_NAN(x)) 2437 #define CHOLMOD_IS_LE_ZERO(x) (((x) <= 0.) && !CHOLMOD_IS_NAN(x)) 2438 2439 #else 2440 2441 /* These all work properly, according to the IEEE 754 standard ... except on */ 2442 /* a PC with windows. Works fine in Linux on the same PC... */ 2443 #define CHOLMOD_IS_NAN(x) ((x) != (x)) 2444 #define CHOLMOD_IS_ZERO(x) ((x) == 0.) 2445 #define CHOLMOD_IS_NONZERO(x) ((x) != 0.) 2446 #define CHOLMOD_IS_LT_ZERO(x) ((x) < 0.) 2447 #define CHOLMOD_IS_GT_ZERO(x) ((x) > 0.) 2448 #define CHOLMOD_IS_LE_ZERO(x) ((x) <= 0.) 2449 2450 #endif 2451 2452 #endif 2453