1 2 /* 3 Provides an interface to the CHOLMOD 1.7.1 sparse solver 4 5 When build with PETSC_USE_64BIT_INDICES this will use UF_Long as the 6 integer type in UMFPACK, otherwise it will use int. This means 7 all integers in this file as simply declared as PetscInt. Also it means 8 that UMFPACK UL_Long version MUST be built with 64 bit integers when used. 9 10 */ 11 12 #include "../src/mat/impls/sbaij/seq/sbaij.h" 13 #include "../src/mat/impls/sbaij/seq/cholmod/cholmodimpl.h" 14 15 /* 16 This is a terrible hack, but it allows the error handler to retain a context. 17 Note that this hack really cannot be made both reentrant and concurrent. 18 */ 19 static Mat static_F; 20 21 #undef __FUNCT__ 22 #define __FUNCT__ "CholmodErrorHandler" 23 static void CholmodErrorHandler(int status,const char *file,int line,const char *message) 24 { 25 26 PetscFunctionBegin; 27 if (status > CHOLMOD_OK) { 28 PetscInfo4(static_F,"CHOLMOD warning %d at %s:%d: %s",status,file,line,message); 29 } else if (status == CHOLMOD_OK) { /* Documentation says this can happen, but why? */ 30 PetscInfo3(static_F,"CHOLMOD OK at %s:%d: %s",file,line,message); 31 } else { 32 PetscErrorPrintf("CHOLMOD error %d at %s:%d: %s\n",status,file,line,message); 33 } 34 PetscFunctionReturnVoid(); 35 } 36 37 #undef __FUNCT__ 38 #define __FUNCT__ "CholmodStart" 39 PetscErrorCode CholmodStart(Mat F) 40 { 41 PetscErrorCode ierr; 42 Mat_CHOLMOD *chol=(Mat_CHOLMOD*)F->spptr; 43 cholmod_common *c; 44 PetscBool flg; 45 46 PetscFunctionBegin; 47 if (chol->common) PetscFunctionReturn(0); 48 ierr = PetscMalloc(sizeof(*chol->common),&chol->common);CHKERRQ(ierr); 49 ierr = !cholmod_X_start(chol->common);CHKERRQ(ierr); 50 c = chol->common; 51 c->error_handler = CholmodErrorHandler; 52 53 #define CHOLMOD_OPTION_DOUBLE(name,help) do { \ 54 PetscReal tmp = (PetscReal)c->name; \ 55 ierr = PetscOptionsReal("-mat_cholmod_" #name,help,"None",tmp,&tmp,0);CHKERRQ(ierr); \ 56 c->name = (double)tmp; \ 57 } while (0) 58 #define CHOLMOD_OPTION_INT(name,help) do { \ 59 PetscInt tmp = (PetscInt)c->name; \ 60 ierr = PetscOptionsInt("-mat_cholmod_" #name,help,"None",tmp,&tmp,0);CHKERRQ(ierr); \ 61 c->name = (int)tmp; \ 62 } while (0) 63 #define CHOLMOD_OPTION_SIZE_T(name,help) do { \ 64 PetscInt tmp = (PetscInt)c->name; \ 65 ierr = PetscOptionsInt("-mat_cholmod_" #name,help,"None",tmp,&tmp,0);CHKERRQ(ierr); \ 66 if (tmp < 0) SETERRQ(((PetscObject)F)->comm,PETSC_ERR_ARG_OUTOFRANGE,"value must be positive"); \ 67 c->name = (size_t)tmp; \ 68 } while (0) 69 #define CHOLMOD_OPTION_TRUTH(name,help) do { \ 70 PetscBool tmp = (PetscBool)!!c->name; \ 71 ierr = PetscOptionsBool("-mat_cholmod_" #name,help,"None",tmp,&tmp,0);CHKERRQ(ierr); \ 72 c->name = (int)tmp; \ 73 } while (0) 74 75 ierr = PetscOptionsBegin(((PetscObject)F)->comm,((PetscObject)F)->prefix,"CHOLMOD Options","Mat");CHKERRQ(ierr); 76 /* CHOLMOD handles first-time packing and refactor-packing separately, but we usually want them to be the same. */ 77 chol->pack = (PetscBool)c->final_pack; 78 ierr = PetscOptionsBool("-mat_cholmod_pack","Pack factors after factorization [disable for frequent repeat factorization]","None",chol->pack,&chol->pack,0);CHKERRQ(ierr); 79 c->final_pack = (int)chol->pack; 80 81 CHOLMOD_OPTION_DOUBLE(dbound,"Minimum absolute value of diagonal entries of D"); 82 CHOLMOD_OPTION_DOUBLE(grow0,"Global growth ratio when factors are modified"); 83 CHOLMOD_OPTION_DOUBLE(grow1,"Column growth ratio when factors are modified"); 84 CHOLMOD_OPTION_SIZE_T(grow2,"Affine column growth constant when factors are modified"); 85 CHOLMOD_OPTION_SIZE_T(maxrank,"Max rank of update, larger values are faster but use more memory [2,4,8]"); 86 { 87 static const char *const list[] = {"SIMPLICIAL","AUTO","SUPERNODAL","MatCholmodFactorType","MAT_CHOLMOD_FACTOR_",0}; 88 PetscEnum choice = (PetscEnum)c->supernodal; 89 ierr = PetscOptionsEnum("-mat_cholmod_factor","Factorization method","None",list,(PetscEnum)c->supernodal,&choice,0);CHKERRQ(ierr); 90 c->supernodal = (int)choice; 91 } 92 if (c->supernodal) CHOLMOD_OPTION_DOUBLE(supernodal_switch,"flop/nnz_L threshold for switching to supernodal factorization"); 93 CHOLMOD_OPTION_TRUTH(final_asis,"Leave factors \"as is\""); 94 CHOLMOD_OPTION_TRUTH(final_pack,"Pack the columns when finished (use FALSE if the factors will be updated later)"); 95 if (!c->final_asis) { 96 CHOLMOD_OPTION_TRUTH(final_super,"Leave supernodal factors instead of converting to simplicial"); 97 CHOLMOD_OPTION_TRUTH(final_ll,"Turn LDL' factorization into LL'"); 98 CHOLMOD_OPTION_TRUTH(final_monotonic,"Ensure columns are monotonic when done"); 99 CHOLMOD_OPTION_TRUTH(final_resymbol,"Remove numerically zero values resulting from relaxed supernodal amalgamation"); 100 } 101 { 102 PetscReal tmp[] = {(PetscReal)c->zrelax[0],(PetscReal)c->zrelax[1],(PetscReal)c->zrelax[2]}; 103 PetscInt n = 3; 104 ierr = PetscOptionsRealArray("-mat_cholmod_zrelax","3 real supernodal relaxed amalgamation parameters","None",tmp,&n,&flg);CHKERRQ(ierr); 105 if (flg && n != 3) SETERRQ(((PetscObject)F)->comm,PETSC_ERR_ARG_OUTOFRANGE,"must provide exactly 3 parameters to -mat_cholmod_zrelax"); 106 if (flg) while (n--) c->zrelax[n] = (double)tmp[n]; 107 } 108 { 109 PetscInt n,tmp[] = {(PetscInt)c->nrelax[0],(PetscInt)c->nrelax[1],(PetscInt)c->nrelax[2]}; 110 ierr = PetscOptionsIntArray("-mat_cholmod_nrelax","3 size_t supernodal relaxed amalgamation parameters","None",tmp,&n,&flg);CHKERRQ(ierr); 111 if (flg && n != 3) SETERRQ(((PetscObject)F)->comm,PETSC_ERR_ARG_OUTOFRANGE,"must provide exactly 3 parameters to -mat_cholmod_nrelax"); 112 if (flg) while (n--) c->nrelax[n] = (size_t)tmp[n]; 113 } 114 CHOLMOD_OPTION_TRUTH(prefer_upper,"Work with upper triangular form [faster when using fill-reducing ordering, slower in natural ordering]"); 115 CHOLMOD_OPTION_TRUTH(default_nesdis,"Use NESDIS instead of METIS for nested dissection"); 116 CHOLMOD_OPTION_INT(print,"Verbosity level"); 117 ierr = PetscOptionsEnd();CHKERRQ(ierr); 118 PetscFunctionReturn(0); 119 } 120 121 #undef __FUNCT__ 122 #define __FUNCT__ "MatWrapCholmod_seqsbaij" 123 static PetscErrorCode MatWrapCholmod_seqsbaij(Mat A,PetscBool values,cholmod_sparse *C,PetscBool *aijalloc) 124 { 125 Mat_SeqSBAIJ *sbaij = (Mat_SeqSBAIJ*)A->data; 126 PetscErrorCode ierr; 127 128 PetscFunctionBegin; 129 ierr = PetscMemzero(C,sizeof(*C));CHKERRQ(ierr); 130 /* CHOLMOD uses column alignment, SBAIJ stores the upper factor, so we pass it on as a lower factor, swapping the meaning of row and column */ 131 C->nrow = (size_t)A->cmap->n; 132 C->ncol = (size_t)A->rmap->n; 133 C->nzmax = (size_t)sbaij->maxnz; 134 C->p = sbaij->i; 135 C->i = sbaij->j; 136 C->x = sbaij->a; 137 C->stype = -1; 138 C->itype = CHOLMOD_INT_TYPE; 139 C->xtype = CHOLMOD_SCALAR_TYPE; 140 C->dtype = CHOLMOD_DOUBLE; 141 C->sorted = 1; 142 C->packed = 1; 143 *aijalloc = PETSC_FALSE; 144 PetscFunctionReturn(0); 145 } 146 147 #undef __FUNCT__ 148 #define __FUNCT__ "VecWrapCholmod" 149 static PetscErrorCode VecWrapCholmod(Vec X,cholmod_dense *Y) 150 { 151 PetscErrorCode ierr; 152 PetscScalar *x; 153 PetscInt n; 154 155 PetscFunctionBegin; 156 ierr = PetscMemzero(Y,sizeof(*Y));CHKERRQ(ierr); 157 ierr = VecGetArray(X,&x);CHKERRQ(ierr); 158 ierr = VecGetSize(X,&n);CHKERRQ(ierr); 159 Y->x = (double*)x; 160 Y->nrow = n; 161 Y->ncol = 1; 162 Y->nzmax = n; 163 Y->d = n; 164 Y->x = (double*)x; 165 Y->xtype = CHOLMOD_SCALAR_TYPE; 166 Y->dtype = CHOLMOD_DOUBLE; 167 PetscFunctionReturn(0); 168 } 169 170 #undef __FUNCT__ 171 #define __FUNCT__ "MatDestroy_CHOLMOD" 172 PetscErrorCode MatDestroy_CHOLMOD(Mat F) 173 { 174 PetscErrorCode ierr; 175 Mat_CHOLMOD *chol=(Mat_CHOLMOD*)F->spptr; 176 177 PetscFunctionBegin; 178 ierr = !cholmod_X_free_factor(&chol->factor,chol->common);CHKERRQ(ierr); 179 ierr = !cholmod_X_finish(chol->common);CHKERRQ(ierr); 180 ierr = PetscFree(chol->common);CHKERRQ(ierr); 181 ierr = PetscFree(chol->matrix);CHKERRQ(ierr); 182 ierr = (*chol->Destroy)(F);CHKERRQ(ierr); 183 PetscFunctionReturn(0); 184 } 185 186 static PetscErrorCode MatSolve_CHOLMOD(Mat,Vec,Vec); 187 188 static const char *const CholmodOrderingMethods[] = {"User","AMD","METIS","NESDIS(default)","Natural","NESDIS(small=20000)","NESDIS(small=4,no constrained)","NESDIS()"}; 189 190 #undef __FUNCT__ 191 #define __FUNCT__ "MatFactorInfo_CHOLMOD" 192 static PetscErrorCode MatFactorInfo_CHOLMOD(Mat F,PetscViewer viewer) 193 { 194 Mat_CHOLMOD *chol = (Mat_CHOLMOD*)F->spptr; 195 const cholmod_common *c = chol->common; 196 PetscErrorCode ierr; 197 PetscInt i; 198 199 PetscFunctionBegin; 200 if (F->ops->solve != MatSolve_CHOLMOD) PetscFunctionReturn(0); 201 ierr = PetscViewerASCIIPrintf(viewer,"CHOLMOD run parameters:\n");CHKERRQ(ierr); 202 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 203 ierr = PetscViewerASCIIPrintf(viewer,"Pack factors after symbolic factorization: %s\n",chol->pack?"TRUE":"FALSE");CHKERRQ(ierr); 204 ierr = PetscViewerASCIIPrintf(viewer,"Common.dbound %g (Smallest absolute value of diagonal entries of D)\n",c->dbound);CHKERRQ(ierr); 205 ierr = PetscViewerASCIIPrintf(viewer,"Common.grow0 %g\n",c->grow0);CHKERRQ(ierr); 206 ierr = PetscViewerASCIIPrintf(viewer,"Common.grow1 %g\n",c->grow1);CHKERRQ(ierr); 207 ierr = PetscViewerASCIIPrintf(viewer,"Common.grow2 %u\n",(unsigned)c->grow2);CHKERRQ(ierr); 208 ierr = PetscViewerASCIIPrintf(viewer,"Common.maxrank %u\n",(unsigned)c->maxrank);CHKERRQ(ierr); 209 ierr = PetscViewerASCIIPrintf(viewer,"Common.supernodal_switch %g\n",c->supernodal_switch);CHKERRQ(ierr); 210 ierr = PetscViewerASCIIPrintf(viewer,"Common.supernodal %d\n",c->supernodal);CHKERRQ(ierr); 211 ierr = PetscViewerASCIIPrintf(viewer,"Common.final_asis %d\n",c->final_asis);CHKERRQ(ierr); 212 ierr = PetscViewerASCIIPrintf(viewer,"Common.final_super %d\n",c->final_super);CHKERRQ(ierr); 213 ierr = PetscViewerASCIIPrintf(viewer,"Common.final_ll %d\n",c->final_ll);CHKERRQ(ierr); 214 ierr = PetscViewerASCIIPrintf(viewer,"Common.final_pack %d\n",c->final_pack);CHKERRQ(ierr); 215 ierr = PetscViewerASCIIPrintf(viewer,"Common.final_monotonic %d\n",c->final_monotonic);CHKERRQ(ierr); 216 ierr = PetscViewerASCIIPrintf(viewer,"Common.final_resymbol %d\n",c->final_resymbol);CHKERRQ(ierr); 217 ierr = PetscViewerASCIIPrintf(viewer,"Common.zrelax [%g,%g,%g]\n",c->zrelax[0],c->zrelax[1],c->zrelax[2]);CHKERRQ(ierr); 218 ierr = PetscViewerASCIIPrintf(viewer,"Common.nrelax [%u,%u,%u]\n",(unsigned)c->nrelax[0],(unsigned)c->nrelax[1],(unsigned)c->nrelax[2]);CHKERRQ(ierr); 219 ierr = PetscViewerASCIIPrintf(viewer,"Common.prefer_upper %d\n",c->prefer_upper);CHKERRQ(ierr); 220 ierr = PetscViewerASCIIPrintf(viewer,"Common.print %d\n",c->print);CHKERRQ(ierr); 221 for (i=0; i<c->nmethods; i++) { 222 ierr = PetscViewerASCIIPrintf(viewer,"Ordering method %D%s:\n",i,i==c->selected?" [SELECTED]":"");CHKERRQ(ierr); 223 ierr = PetscViewerASCIIPrintf(viewer," lnz %g, fl %g, prune_dense %g, prune_dense2 %g\n", 224 c->method[i].lnz,c->method[i].fl,c->method[i].prune_dense,c->method[i].prune_dense2);CHKERRQ(ierr); 225 } 226 ierr = PetscViewerASCIIPrintf(viewer,"Common.postorder %d\n",c->postorder);CHKERRQ(ierr); 227 ierr = PetscViewerASCIIPrintf(viewer,"Common.default_nesdis %d (use NESDIS instead of METIS for nested dissection)\n",c->default_nesdis);CHKERRQ(ierr); 228 /* Statistics */ 229 ierr = PetscViewerASCIIPrintf(viewer,"Common.fl %g (flop count from most recent analysis)\n",c->fl);CHKERRQ(ierr); 230 ierr = PetscViewerASCIIPrintf(viewer,"Common.lnz %g (fundamental nz in L)\n",c->lnz);CHKERRQ(ierr); 231 ierr = PetscViewerASCIIPrintf(viewer,"Common.anz %g\n",c->anz);CHKERRQ(ierr); 232 ierr = PetscViewerASCIIPrintf(viewer,"Common.modfl %g (flop count from most recent update)\n",c->modfl);CHKERRQ(ierr); 233 ierr = PetscViewerASCIIPrintf(viewer,"Common.malloc_count %g (number of live objects)\n",(double)c->malloc_count);CHKERRQ(ierr);CHKERRQ(ierr); 234 ierr = PetscViewerASCIIPrintf(viewer,"Common.memory_usage %g (peak memory usage in bytes)\n",(double)c->memory_usage);CHKERRQ(ierr);CHKERRQ(ierr); 235 ierr = PetscViewerASCIIPrintf(viewer,"Common.memory_inuse %g (current memory usage in bytes)\n",(double)c->memory_inuse);CHKERRQ(ierr);CHKERRQ(ierr); 236 ierr = PetscViewerASCIIPrintf(viewer,"Common.nrealloc_col %g (number of column reallocations)\n",c->nrealloc_col);CHKERRQ(ierr);CHKERRQ(ierr); 237 ierr = PetscViewerASCIIPrintf(viewer,"Common.nrealloc_factor %g (number of factor reallocations due to column reallocations)\n",c->nrealloc_factor);CHKERRQ(ierr);CHKERRQ(ierr); 238 ierr = PetscViewerASCIIPrintf(viewer,"Common.ndbounds_hit %g (number of times diagonal was modified by dbound)\n",c->ndbounds_hit);CHKERRQ(ierr);CHKERRQ(ierr); 239 ierr = PetscViewerASCIIPrintf(viewer,"Common.rowfacfl %g (number of flops in last call to cholmod_rowfac)\n",c->rowfacfl);CHKERRQ(ierr);CHKERRQ(ierr); 240 ierr = PetscViewerASCIIPrintf(viewer,"Common.aatfl %g (number of flops to compute A(:,f)*A(:,f)')\n",c->aatfl);CHKERRQ(ierr);CHKERRQ(ierr); 241 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 242 PetscFunctionReturn(0); 243 } 244 245 #undef __FUNCT__ 246 #define __FUNCT__ "MatView_CHOLMOD" 247 PetscErrorCode MatView_CHOLMOD(Mat F,PetscViewer viewer) 248 { 249 PetscErrorCode ierr; 250 PetscBool iascii; 251 PetscViewerFormat format; 252 253 PetscFunctionBegin; 254 ierr = MatView_SeqSBAIJ(F,viewer);CHKERRQ(ierr); 255 ierr = PetscTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 256 if (iascii) { 257 ierr = PetscViewerGetFormat(viewer,&format);CHKERRQ(ierr); 258 if (format == PETSC_VIEWER_ASCII_INFO) { 259 ierr = MatFactorInfo_CHOLMOD(F,viewer);CHKERRQ(ierr); 260 } 261 } 262 PetscFunctionReturn(0); 263 } 264 265 #undef __FUNCT__ 266 #define __FUNCT__ "MatSolve_CHOLMOD" 267 static PetscErrorCode MatSolve_CHOLMOD(Mat F,Vec B,Vec X) 268 { 269 Mat_CHOLMOD *chol = (Mat_CHOLMOD*)F->spptr; 270 cholmod_dense cholB,*cholX; 271 PetscScalar *x; 272 PetscErrorCode ierr; 273 274 PetscFunctionBegin; 275 ierr = VecWrapCholmod(B,&cholB);CHKERRQ(ierr); 276 static_F = F; 277 cholX = cholmod_X_solve(CHOLMOD_A,chol->factor,&cholB,chol->common); 278 if (!cholX) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"CHOLMOD failed"); 279 ierr = VecGetArray(X,&x);CHKERRQ(ierr); 280 ierr = PetscMemcpy(x,cholX->x,cholX->nrow*sizeof(*x));CHKERRQ(ierr); 281 ierr = !cholmod_X_free_dense(&cholX,chol->common);CHKERRQ(ierr); 282 ierr = VecRestoreArray(X,&x);CHKERRQ(ierr); 283 PetscFunctionReturn(0); 284 } 285 286 #undef __FUNCT__ 287 #define __FUNCT__ "MatCholeskyFactorNumeric_CHOLMOD" 288 static PetscErrorCode MatCholeskyFactorNumeric_CHOLMOD(Mat F,Mat A,const MatFactorInfo *info) 289 { 290 Mat_CHOLMOD *chol = (Mat_CHOLMOD*)F->spptr; 291 cholmod_sparse cholA; 292 PetscBool aijalloc; 293 PetscErrorCode ierr; 294 295 PetscFunctionBegin; 296 ierr = (*chol->Wrap)(A,PETSC_TRUE,&cholA,&aijalloc);CHKERRQ(ierr); 297 static_F = F; 298 ierr = !cholmod_X_factorize(&cholA,chol->factor,chol->common); 299 if (ierr) SETERRQ1(((PetscObject)F)->comm,PETSC_ERR_LIB,"CHOLMOD factorization failed with status %d",chol->common->status); 300 if (chol->common->status == CHOLMOD_NOT_POSDEF) 301 SETERRQ1(((PetscObject)F)->comm,PETSC_ERR_MAT_CH_ZRPVT,"CHOLMOD detected that the matrix is not positive definite, failure at column %u",(unsigned)chol->factor->minor); 302 303 if (aijalloc) {ierr = PetscFree3(cholA.p,cholA.i,cholA.x);CHKERRQ(ierr);} 304 305 F->ops->solve = MatSolve_CHOLMOD; 306 F->ops->solvetranspose = MatSolve_CHOLMOD; 307 PetscFunctionReturn(0); 308 } 309 310 #undef __FUNCT__ 311 #define __FUNCT__ "MatCholeskyFactorSymbolic_CHOLMOD" 312 PetscErrorCode MatCholeskyFactorSymbolic_CHOLMOD(Mat F,Mat A,IS perm,const MatFactorInfo *info) 313 { 314 Mat_CHOLMOD *chol = (Mat_CHOLMOD*)F->spptr; 315 PetscErrorCode ierr; 316 cholmod_sparse cholA; 317 PetscBool aijalloc; 318 PetscInt *fset = 0; 319 size_t fsize = 0; 320 321 PetscFunctionBegin; 322 ierr = (*chol->Wrap)(A,PETSC_FALSE,&cholA,&aijalloc);CHKERRQ(ierr); 323 static_F = F; 324 if (chol->factor) { 325 ierr = !cholmod_X_resymbol(&cholA,fset,fsize,(int)chol->pack,chol->factor,chol->common); 326 if (ierr) SETERRQ1(((PetscObject)F)->comm,PETSC_ERR_LIB,"CHOLMOD analysis failed with status %d",chol->common->status); 327 } else if (perm) { 328 const PetscInt *ip; 329 ierr = ISGetIndices(perm,&ip);CHKERRQ(ierr); 330 chol->factor = cholmod_X_analyze_p(&cholA,(PetscInt*)ip,fset,fsize,chol->common); 331 if (!chol->factor) SETERRQ1(((PetscObject)F)->comm,PETSC_ERR_LIB,"CHOLMOD analysis failed with status %d",chol->common->status); 332 ierr = ISRestoreIndices(perm,&ip);CHKERRQ(ierr); 333 } else { 334 chol->factor = cholmod_X_analyze(&cholA,chol->common); 335 if (!chol->factor) SETERRQ1(((PetscObject)F)->comm,PETSC_ERR_LIB,"CHOLMOD analysis failed with status %d",chol->common->status); 336 } 337 338 if (aijalloc) {ierr = PetscFree3(cholA.p,cholA.i,cholA.x);CHKERRQ(ierr);} 339 340 F->ops->choleskyfactornumeric = MatCholeskyFactorNumeric_CHOLMOD; 341 PetscFunctionReturn(0); 342 } 343 344 EXTERN_C_BEGIN 345 #undef __FUNCT__ 346 #define __FUNCT__ "MatFactorGetSolverPackage_seqsbaij_cholmod" 347 PetscErrorCode MatFactorGetSolverPackage_seqsbaij_cholmod(Mat A,const MatSolverPackage *type) 348 { 349 PetscFunctionBegin; 350 *type = MATSOLVERCHOLMOD; 351 PetscFunctionReturn(0); 352 } 353 EXTERN_C_END 354 355 /*MC 356 MATSOLVERCHOLMOD = "cholmod" - A matrix type providing direct solvers (Cholesky) for sequential matrices 357 via the external package CHOLMOD. 358 359 ./configure --download-cholmod to install PETSc to use CHOLMOD 360 361 Consult CHOLMOD documentation for more information about the Common parameters 362 which correspond to the options database keys below. 363 364 Options Database Keys: 365 -mat_cholmod_dbound <0>: Minimum absolute value of diagonal entries of D (None) 366 -mat_cholmod_grow0 <1.2>: Global growth ratio when factors are modified (None) 367 -mat_cholmod_grow1 <1.2>: Column growth ratio when factors are modified (None) 368 -mat_cholmod_grow2 <5>: Affine column growth constant when factors are modified (None) 369 -mat_cholmod_maxrank <8>: Max rank of update, larger values are faster but use more memory [2,4,8] (None) 370 -mat_cholmod_factor <AUTO> (choose one of) SIMPLICIAL AUTO SUPERNODAL 371 -mat_cholmod_supernodal_switch <40>: flop/nnz_L threshold for switching to supernodal factorization (None) 372 -mat_cholmod_final_asis: <TRUE> Leave factors "as is" (None) 373 -mat_cholmod_final_pack: <TRUE> Pack the columns when finished (use FALSE if the factors will be updated later) (None) 374 -mat_cholmod_zrelax <0.8>: 3 real supernodal relaxed amalgamation parameters (None) 375 -mat_cholmod_nrelax <4>: 3 size_t supernodal relaxed amalgamation parameters (None) 376 -mat_cholmod_prefer_upper: <TRUE> Work with upper triangular form (faster when using fill-reducing ordering, slower in natural ordering) (None) 377 -mat_cholmod_print <3>: Verbosity level (None) 378 379 Level: beginner 380 381 .seealso: PCCHOLESKY, PCFactorSetMatSolverPackage(), MatSolverPackage 382 M*/ 383 EXTERN_C_BEGIN 384 #undef __FUNCT__ 385 #define __FUNCT__ "MatGetFactor_seqsbaij_cholmod" 386 PetscErrorCode MatGetFactor_seqsbaij_cholmod(Mat A,MatFactorType ftype,Mat *F) 387 { 388 Mat B; 389 Mat_CHOLMOD *chol; 390 PetscErrorCode ierr; 391 PetscInt m=A->rmap->n,n=A->cmap->n,bs; 392 393 PetscFunctionBegin; 394 if (ftype != MAT_FACTOR_CHOLESKY) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_SUP,"CHOLMOD cannot do %s factorization with SBAIJ, only %s", 395 MatFactorTypes[ftype],MatFactorTypes[MAT_FACTOR_CHOLESKY]); 396 ierr = MatGetBlockSize(A,&bs);CHKERRQ(ierr); 397 if (bs != 1) SETERRQ1(((PetscObject)A)->comm,PETSC_ERR_SUP,"CHOLMOD only supports block size=1, given %D",bs); 398 /* Create the factorization matrix F */ 399 ierr = MatCreate(((PetscObject)A)->comm,&B);CHKERRQ(ierr); 400 ierr = MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,m,n);CHKERRQ(ierr); 401 ierr = MatSetType(B,((PetscObject)A)->type_name);CHKERRQ(ierr); 402 ierr = MatSeqSBAIJSetPreallocation(B,1,0,PETSC_NULL);CHKERRQ(ierr); 403 ierr = PetscNewLog(B,Mat_CHOLMOD,&chol);CHKERRQ(ierr); 404 chol->Wrap = MatWrapCholmod_seqsbaij; 405 chol->Destroy = MatDestroy_SeqSBAIJ; 406 B->spptr = chol; 407 408 B->ops->view = MatView_CHOLMOD; 409 B->ops->choleskyfactorsymbolic = MatCholeskyFactorSymbolic_CHOLMOD; 410 B->ops->destroy = MatDestroy_CHOLMOD; 411 ierr = PetscObjectComposeFunctionDynamic((PetscObject)B,"MatFactorGetSolverPackage_C","MatFactorGetSolverPackage_seqsbaij_cholmod",MatFactorGetSolverPackage_seqsbaij_cholmod);CHKERRQ(ierr); 412 B->factortype = MAT_FACTOR_CHOLESKY; 413 B->assembled = PETSC_TRUE; /* required by -ksp_view */ 414 B->preallocated = PETSC_TRUE; 415 416 ierr = CholmodStart(B);CHKERRQ(ierr); 417 *F = B; 418 PetscFunctionReturn(0); 419 } 420 EXTERN_C_END 421