1 /* ========================================================================== */
2 /* === umf_cholmod ========================================================== */
3 /* ========================================================================== */
4
5 /* umfpack_cholmod: user-defined ordering function to interface UMFPACK
6 * to CHOLMOD.
7 *
8 * This routine is an example of a user-provided ordering function for UMFPACK.
9 *
10 * This function can be passed to umfpack_*_fsymbolic as the
11 * user_ordering function pointer.
12 */
13
14 #include "umf_internal.h"
15 #include "umf_cholmod.h"
16
17 #ifndef NCHOLMOD
18 #include "cholmod.h"
19 #endif
20
21 #if defined (DINT) || defined (ZINT)
22 #define CHOLMOD_start cholmod_start
23 #define CHOLMOD_transpose cholmod_transpose
24 #define CHOLMOD_analyze cholmod_analyze
25 #define CHOLMOD_free_sparse cholmod_free_sparse
26 #define CHOLMOD_free_factor cholmod_free_factor
27 #define CHOLMOD_finish cholmod_finish
28 #define CHOLMOD_print_common cholmod_print_common
29 #else
30 #define CHOLMOD_start cholmod_l_start
31 #define CHOLMOD_transpose cholmod_l_transpose
32 #define CHOLMOD_analyze cholmod_l_analyze
33 #define CHOLMOD_free_sparse cholmod_l_free_sparse
34 #define CHOLMOD_free_factor cholmod_l_free_factor
35 #define CHOLMOD_finish cholmod_l_finish
36 #define CHOLMOD_print_common cholmod_l_print_common
37 #endif
38
UMF_cholmod(Int nrow,Int ncol,Int symmetric,Int Ap[],Int Ai[],Int Perm[],void * user_params,double user_info[3])39 int UMF_cholmod
40 (
41 /* inputs */
42 Int nrow, /* A is nrow-by-ncol */
43 Int ncol, /* A is nrow-by-ncol */
44 Int symmetric, /* if true and nrow=ncol do A+A', else do A'A */
45 Int Ap [ ], /* column pointers, size ncol+1 */
46 Int Ai [ ], /* row indices, size nz = Ap [ncol] */
47 /* output */
48 Int Perm [ ], /* fill-reducing permutation, size ncol */
49 /* user-defined */
50 void *user_params, /* Int array of size 3 */
51 double user_info [3] /* [0]: max col count for L=chol(P(A+A')P')
52 [1]: nnz (L)
53 [2]: flop count for chol, if A real */
54 )
55 {
56 #ifndef NCHOLMOD
57 double dmax, flops, c, lnz ;
58 cholmod_sparse Amatrix, *A, *AT, *S ;
59 cholmod_factor *L ;
60 cholmod_common cm ;
61 Int *P, *ColCount ;
62 Int k, ordering_option, print_level, *params ;
63
64 params = (Int *) user_params ;
65 ordering_option = params [0] ;
66 print_level = params [1] - 1 ;
67 params [2] = -1 ;
68
69 if (Ap == NULL || Ai == NULL || Perm == NULL || nrow < 0 || ncol < 0)
70 {
71 /* invalid inputs */
72 return (FALSE) ;
73 }
74 if (nrow != ncol)
75 {
76 /* force symmetric to be false */
77 symmetric = FALSE ;
78 }
79
80 /* start CHOLMOD */
81 CHOLMOD_start (&cm) ;
82 cm.supernodal = CHOLMOD_SIMPLICIAL ;
83 cm.print = print_level ;
84
85 /* adjust cm based on ordering_option */
86 switch (ordering_option)
87 {
88
89 default:
90 case UMFPACK_ORDERING_AMD:
91 /* AMD on A+A' if symmetric, COLAMD on A otherwise */
92 cm.nmethods = 1 ;
93 cm.method [0].ordering = symmetric ? CHOLMOD_AMD : CHOLMOD_COLAMD ;
94 cm.postorder = TRUE ;
95 break ;
96
97 case UMFPACK_ORDERING_METIS:
98 /* metis on A+A' if symmetric, A'A otherwise */
99 cm.nmethods = 1 ;
100 cm.method [0].ordering = CHOLMOD_METIS ;
101 cm.postorder = TRUE ;
102 break ;
103
104 case UMFPACK_ORDERING_NONE:
105 case UMFPACK_ORDERING_GIVEN:
106 case UMFPACK_ORDERING_USER:
107 /* no ordering. No input permutation here, and no user
108 function, so all these are the same as "none". */
109 cm.nmethods = 1 ;
110 cm.method [0].ordering = CHOLMOD_NATURAL ;
111 cm.postorder = FALSE ;
112 break ;
113
114 case UMFPACK_ORDERING_BEST:
115 /* try AMD, METIS and NESDIS on A+A', or COLAMD(A), METIS(A'A),
116 and NESDIS (A'A) */
117 cm.nmethods = 3 ;
118 cm.method [0].ordering = symmetric ? CHOLMOD_AMD : CHOLMOD_COLAMD ;
119 cm.method [1].ordering = CHOLMOD_METIS ;
120 cm.method [2].ordering = CHOLMOD_NESDIS ;
121 cm.postorder = TRUE ;
122 break ;
123
124 case UMFPACK_ORDERING_CHOLMOD:
125 /* no change to CHOLMOD defaults:
126 Do not use given permutation, since it's not provided.
127 Try AMD. If fill-in and flop count are low, use AMD.
128 Otherwise, try METIS and take the best of AMD and METIS.
129 cm.method [0].ordering = CHOLMOD_GIVEN
130 cm.method [1].ordering = CHOLMOD_AMD
131 cm.method [2].ordering = CHOLMOD_METIS
132 cm.nmethods = 2 if METIS installed, 3 otherwise ('given' is skipped)
133 */
134 break ;
135 }
136
137 /* construct a CHOLMOD version of the input matrix A */
138 A = &Amatrix ;
139 A->nrow = nrow ; /* A is nrow-by-ncol */
140 A->ncol = ncol ;
141 A->nzmax = Ap [ncol] ; /* with nzmax entries */
142 A->packed = TRUE ; /* there is no A->nz array */
143 if (symmetric)
144 {
145 A->stype = 1 ; /* A is symmetric */
146 }
147 else
148 {
149 A->stype = 0 ; /* A is unsymmetric */
150 }
151 A->itype = CHOLMOD_INT ;
152 A->xtype = CHOLMOD_PATTERN ;
153 A->dtype = CHOLMOD_DOUBLE ;
154 A->nz = NULL ;
155 A->p = Ap ; /* column pointers */
156 A->i = Ai ; /* row indices */
157 A->x = NULL ; /* no numerical values */
158 A->z = NULL ;
159 A->sorted = FALSE ; /* columns of A might not be sorted */
160
161 if (symmetric)
162 {
163 /* CHOLMOD with order the symmetric matrix A */
164 AT = NULL ;
165 S = A ;
166 }
167 else
168 {
169 /* S = A'. CHOLMOD will order S*S', which is A'*A */
170 AT = CHOLMOD_transpose (A, 0, &cm) ;
171 S = AT ;
172 }
173
174 /* order and analyze S or S*S' */
175 L = CHOLMOD_analyze (S, &cm) ;
176 CHOLMOD_free_sparse (&AT, &cm) ;
177 if (L == NULL)
178 {
179 return (FALSE) ;
180 }
181
182 /* determine the ordering used */
183 switch (L->ordering)
184 {
185
186 case CHOLMOD_AMD:
187 case CHOLMOD_COLAMD:
188 params [2] = UMFPACK_ORDERING_AMD ;
189 break ;
190
191 case CHOLMOD_METIS:
192 case CHOLMOD_NESDIS:
193 params [2] = UMFPACK_ORDERING_METIS ;
194 break ;
195
196 case CHOLMOD_GIVEN:
197 case CHOLMOD_NATURAL:
198 default:
199 params [2] = UMFPACK_ORDERING_NONE ;
200 break ;
201 }
202
203 /* copy the permutation from L to the output and compute statistics */
204 P = L->Perm ;
205 ColCount = L->ColCount ;
206 dmax = 1 ;
207 lnz = 0 ;
208 flops = 0 ;
209 for (k = 0 ; k < ncol ; k++)
210 {
211 Perm [k] = P [k] ;
212 c = ColCount [k] ;
213 if (c > dmax) dmax = c ;
214 lnz += c ;
215 flops += c*c ;
216 }
217 user_info [0] = dmax ;
218 user_info [1] = lnz ;
219 user_info [2] = flops ;
220
221 CHOLMOD_free_factor (&L, &cm) ;
222 if (print_level > 0)
223 {
224 CHOLMOD_print_common ("for UMFPACK", &cm) ;
225 }
226 CHOLMOD_finish (&cm) ;
227 return (TRUE) ;
228 #else
229 /* CHOLMOD and its supporting packages (CAMD, CCOLAMD, COLAMD, METIS)
230 not installed */
231 return (FALSE) ;
232 #endif
233 }
234