1 /* ========================================================================== */ 2 /* === Include/cholmod_supernodal.h ========================================= */ 3 /* ========================================================================== */ 4 5 /* ----------------------------------------------------------------------------- 6 * CHOLMOD/Include/cholmod_supernodal.h. 7 * Copyright (C) 2005-2006, Timothy A. Davis 8 * http://www.suitesparse.com 9 * -------------------------------------------------------------------------- */ 10 11 /* CHOLMOD Supernodal module. 12 * 13 * Supernodal analysis, factorization, and solve. The simplest way to use 14 * these routines is via the Cholesky module. It does not provide any 15 * fill-reducing orderings, but does accept the orderings computed by the 16 * Cholesky module. It does not require the Cholesky module itself, however. 17 * 18 * Primary routines: 19 * ----------------- 20 * cholmod_super_symbolic supernodal symbolic analysis 21 * cholmod_super_numeric supernodal numeric factorization 22 * cholmod_super_lsolve supernodal Lx=b solve 23 * cholmod_super_ltsolve supernodal L'x=b solve 24 * 25 * Prototypes for the BLAS and LAPACK routines that CHOLMOD uses are listed 26 * below, including how they are used in CHOLMOD. 27 * 28 * BLAS routines: 29 * -------------- 30 * dtrsv solve Lx=b or L'x=b, L non-unit diagonal, x and b stride-1 31 * dtrsm solve LX=B or L'X=b, L non-unit diagonal 32 * dgemv y=y-A*x or y=y-A'*x (x and y stride-1) 33 * dgemm C=A*B', C=C-A*B, or C=C-A'*B 34 * dsyrk C=tril(A*A') 35 * 36 * LAPACK routines: 37 * ---------------- 38 * dpotrf LAPACK: A=chol(tril(A)) 39 * 40 * Requires the Core module, and two external packages: LAPACK and the BLAS. 41 * Optionally used by the Cholesky module. 42 */ 43 44 #ifndef CHOLMOD_SUPERNODAL_H 45 #define CHOLMOD_SUPERNODAL_H 46 47 #include "cholmod_core.h" 48 49 /* -------------------------------------------------------------------------- */ 50 /* cholmod_super_symbolic */ 51 /* -------------------------------------------------------------------------- */ 52 53 /* Analyzes A, AA', or A(:,f)*A(:,f)' in preparation for a supernodal numeric 54 * factorization. The user need not call this directly; cholmod_analyze is 55 * a "simple" wrapper for this routine. 56 */ 57 58 int cholmod_super_symbolic 59 ( 60 /* ---- input ---- */ 61 cholmod_sparse *A, /* matrix to analyze */ 62 cholmod_sparse *F, /* F = A' or A(:,f)' */ 63 int *Parent, /* elimination tree */ 64 /* ---- in/out --- */ 65 cholmod_factor *L, /* simplicial symbolic on input, 66 * supernodal symbolic on output */ 67 /* --------------- */ 68 cholmod_common *Common 69 ) ; 70 71 int cholmod_l_super_symbolic (cholmod_sparse *, cholmod_sparse *, 72 SuiteSparse_long *, cholmod_factor *, cholmod_common *) ; 73 74 /* -------------------------------------------------------------------------- */ 75 /* cholmod_super_symbolic2 */ 76 /* -------------------------------------------------------------------------- */ 77 78 /* Analyze for supernodal Cholesky or multifrontal QR */ 79 80 int cholmod_super_symbolic2 81 ( 82 /* ---- input ---- */ 83 int for_whom, /* FOR_SPQR (0): for SPQR but not GPU-accelerated 84 FOR_CHOLESKY (1): for Cholesky (GPU or not) 85 FOR_SPQRGPU (2): for SPQR with GPU acceleration */ 86 cholmod_sparse *A, /* matrix to analyze */ 87 cholmod_sparse *F, /* F = A' or A(:,f)' */ 88 int *Parent, /* elimination tree */ 89 /* ---- in/out --- */ 90 cholmod_factor *L, /* simplicial symbolic on input, 91 * supernodal symbolic on output */ 92 /* --------------- */ 93 cholmod_common *Common 94 ) ; 95 96 int cholmod_l_super_symbolic2 (int, cholmod_sparse *, cholmod_sparse *, 97 SuiteSparse_long *, cholmod_factor *, cholmod_common *) ; 98 99 /* -------------------------------------------------------------------------- */ 100 /* cholmod_super_numeric */ 101 /* -------------------------------------------------------------------------- */ 102 103 /* Computes the numeric LL' factorization of A, AA', or A(:,f)*A(:,f)' using 104 * a BLAS-based supernodal method. The user need not call this directly; 105 * cholmod_factorize is a "simple" wrapper for this routine. 106 */ 107 108 int cholmod_super_numeric 109 ( 110 /* ---- input ---- */ 111 cholmod_sparse *A, /* matrix to factorize */ 112 cholmod_sparse *F, /* F = A' or A(:,f)' */ 113 double beta [2], /* beta*I is added to diagonal of matrix to factorize */ 114 /* ---- in/out --- */ 115 cholmod_factor *L, /* factorization */ 116 /* --------------- */ 117 cholmod_common *Common 118 ) ; 119 120 int cholmod_l_super_numeric (cholmod_sparse *, cholmod_sparse *, double *, 121 cholmod_factor *, cholmod_common *) ; 122 123 /* -------------------------------------------------------------------------- */ 124 /* cholmod_super_lsolve */ 125 /* -------------------------------------------------------------------------- */ 126 127 /* Solve Lx=b where L is from a supernodal numeric factorization. The user 128 * need not call this routine directly. cholmod_solve is a "simple" wrapper 129 * for this routine. */ 130 131 int cholmod_super_lsolve 132 ( 133 /* ---- input ---- */ 134 cholmod_factor *L, /* factor to use for the forward solve */ 135 /* ---- output ---- */ 136 cholmod_dense *X, /* b on input, solution to Lx=b on output */ 137 /* ---- workspace */ 138 cholmod_dense *E, /* workspace of size nrhs*(L->maxesize) */ 139 /* --------------- */ 140 cholmod_common *Common 141 ) ; 142 143 int cholmod_l_super_lsolve (cholmod_factor *, cholmod_dense *, cholmod_dense *, 144 cholmod_common *) ; 145 146 /* -------------------------------------------------------------------------- */ 147 /* cholmod_super_ltsolve */ 148 /* -------------------------------------------------------------------------- */ 149 150 /* Solve L'x=b where L is from a supernodal numeric factorization. The user 151 * need not call this routine directly. cholmod_solve is a "simple" wrapper 152 * for this routine. */ 153 154 int cholmod_super_ltsolve 155 ( 156 /* ---- input ---- */ 157 cholmod_factor *L, /* factor to use for the backsolve */ 158 /* ---- output ---- */ 159 cholmod_dense *X, /* b on input, solution to L'x=b on output */ 160 /* ---- workspace */ 161 cholmod_dense *E, /* workspace of size nrhs*(L->maxesize) */ 162 /* --------------- */ 163 cholmod_common *Common 164 ) ; 165 166 int cholmod_l_super_ltsolve (cholmod_factor *, cholmod_dense *, cholmod_dense *, 167 cholmod_common *) ; 168 169 #endif 170