1 #define PETSCKSP_DLL 2 3 /* 4 Defines a direct factorization preconditioner for any Mat implementation 5 Note: this need not be consided a preconditioner since it supplies 6 a direct solver. 7 */ 8 #include "private/pcimpl.h" /*I "petscpc.h" I*/ 9 10 typedef struct { 11 Mat fact; /* factored matrix */ 12 PetscReal actualfill; /* actual fill in factor */ 13 PetscTruth inplace; /* flag indicating in-place factorization */ 14 IS row,col; /* index sets used for reordering */ 15 MatOrderingType ordering; /* matrix ordering */ 16 PetscTruth reuseordering; /* reuses previous reordering computed */ 17 PetscTruth reusefill; /* reuse fill from previous Cholesky */ 18 MatFactorInfo info; 19 MatSolverType solvertype; 20 } PC_Cholesky; 21 22 EXTERN_C_BEGIN 23 #undef __FUNCT__ 24 #define __FUNCT__ "PCFactorSetZeroPivot_Cholesky" 25 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetZeroPivot_Cholesky(PC pc,PetscReal z) 26 { 27 PC_Cholesky *ch; 28 29 PetscFunctionBegin; 30 ch = (PC_Cholesky*)pc->data; 31 ch->info.zeropivot = z; 32 PetscFunctionReturn(0); 33 } 34 EXTERN_C_END 35 36 EXTERN_C_BEGIN 37 #undef __FUNCT__ 38 #define __FUNCT__ "PCFactorSetShiftNonzero_Cholesky" 39 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftNonzero_Cholesky(PC pc,PetscReal shift) 40 { 41 PC_Cholesky *dir; 42 43 PetscFunctionBegin; 44 dir = (PC_Cholesky*)pc->data; 45 if (shift == (PetscReal) PETSC_DECIDE) { 46 dir->info.shiftnz = 1.e-12; 47 } else { 48 dir->info.shiftnz = shift; 49 } 50 PetscFunctionReturn(0); 51 } 52 EXTERN_C_END 53 54 EXTERN_C_BEGIN 55 #undef __FUNCT__ 56 #define __FUNCT__ "PCFactorSetShiftPd_Cholesky" 57 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftPd_Cholesky(PC pc,PetscTruth shift) 58 { 59 PC_Cholesky *dir; 60 61 PetscFunctionBegin; 62 dir = (PC_Cholesky*)pc->data; 63 if (shift) { 64 dir->info.shift_fraction = 0.0; 65 dir->info.shiftpd = 1.0; 66 } else { 67 dir->info.shiftpd = 0.0; 68 } 69 PetscFunctionReturn(0); 70 } 71 EXTERN_C_END 72 73 EXTERN_C_BEGIN 74 #undef __FUNCT__ 75 #define __FUNCT__ "PCFactorSetReuseOrdering_Cholesky" 76 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseOrdering_Cholesky(PC pc,PetscTruth flag) 77 { 78 PC_Cholesky *lu; 79 80 PetscFunctionBegin; 81 lu = (PC_Cholesky*)pc->data; 82 lu->reuseordering = flag; 83 PetscFunctionReturn(0); 84 } 85 EXTERN_C_END 86 87 EXTERN_C_BEGIN 88 #undef __FUNCT__ 89 #define __FUNCT__ "PCFactorSetReuseFill_Cholesky" 90 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseFill_Cholesky(PC pc,PetscTruth flag) 91 { 92 PC_Cholesky *lu; 93 94 PetscFunctionBegin; 95 lu = (PC_Cholesky*)pc->data; 96 lu->reusefill = flag; 97 PetscFunctionReturn(0); 98 } 99 EXTERN_C_END 100 101 #undef __FUNCT__ 102 #define __FUNCT__ "PCSetFromOptions_Cholesky" 103 static PetscErrorCode PCSetFromOptions_Cholesky(PC pc) 104 { 105 PC_Cholesky *lu = (PC_Cholesky*)pc->data; 106 PetscErrorCode ierr; 107 PetscTruth flg; 108 char tname[256], solvertype[64]; 109 PetscFList ordlist; 110 111 PetscFunctionBegin; 112 ierr = MatOrderingRegisterAll(PETSC_NULL);CHKERRQ(ierr); 113 ierr = PetscOptionsHead("Cholesky options");CHKERRQ(ierr); 114 ierr = PetscOptionsName("-pc_factor_in_place","Form Cholesky in the same memory as the matrix","PCFactorSetUseInPlace",&flg);CHKERRQ(ierr); 115 if (flg) { 116 ierr = PCFactorSetUseInPlace(pc);CHKERRQ(ierr); 117 } 118 ierr = PetscOptionsReal("-pc_factor_fill","Expected non-zeros in Cholesky/non-zeros in matrix","PCFactorSetFill",lu->info.fill,&lu->info.fill,0);CHKERRQ(ierr); 119 120 ierr = PetscOptionsName("-pc_factor_reuse_fill","Use fill from previous factorization","PCFactorSetReuseFill",&flg);CHKERRQ(ierr); 121 if (flg) { 122 ierr = PCFactorSetReuseFill(pc,PETSC_TRUE);CHKERRQ(ierr); 123 } 124 ierr = PetscOptionsName("-pc_factor_reuse_ordering","Reuse ordering from previous factorization","PCFactorSetReuseOrdering",&flg);CHKERRQ(ierr); 125 if (flg) { 126 ierr = PCFactorSetReuseOrdering(pc,PETSC_TRUE);CHKERRQ(ierr); 127 } 128 129 ierr = MatGetOrderingList(&ordlist);CHKERRQ(ierr); 130 ierr = PetscOptionsList("-pc_factor_mat_ordering_type","Reordering to reduce nonzeros in Cholesky","PCFactorSetMatOrderingType",ordlist,lu->ordering,tname,256,&flg);CHKERRQ(ierr); 131 if (flg) { 132 ierr = PCFactorSetMatOrderingType(pc,tname);CHKERRQ(ierr); 133 } 134 135 /* maybe should have MatGetSolverTypes(Mat,&list) like the ordering list */ 136 ierr = PetscOptionsString("-pc_mat_solver_type","Specific Cholesky solver to use","MatGetFactor",lu->solvertype,solvertype,64,&flg);CHKERRQ(ierr); 137 if (flg) { 138 ierr = PetscStrallocpy(solvertype,&lu->solvertype);CHKERRQ(ierr); 139 } 140 141 ierr = PetscOptionsName("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",&flg);CHKERRQ(ierr); 142 if (flg) { 143 ierr = PCFactorSetShiftNonzero(pc,(PetscReal) PETSC_DECIDE);CHKERRQ(ierr); 144 } 145 ierr = PetscOptionsReal("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",lu->info.shiftnz,&lu->info.shiftnz,0);CHKERRQ(ierr); 146 ierr = PetscOptionsName("-pc_factor_shift_positive_definite","Manteuffel shift applied to diagonal","PCFactorSetShiftPd",&flg);CHKERRQ(ierr); 147 if (flg) { 148 ierr = PCFactorSetShiftPd(pc,PETSC_TRUE);CHKERRQ(ierr); 149 } 150 ierr = PetscOptionsReal("-pc_factor_zeropivot","Pivot is considered zero if less than","PCFactorSetZeroPivot",lu->info.zeropivot,&lu->info.zeropivot,0);CHKERRQ(ierr); 151 152 ierr = PetscOptionsTail();CHKERRQ(ierr); 153 PetscFunctionReturn(0); 154 } 155 156 #undef __FUNCT__ 157 #define __FUNCT__ "PCView_Cholesky" 158 static PetscErrorCode PCView_Cholesky(PC pc,PetscViewer viewer) 159 { 160 PC_Cholesky *lu = (PC_Cholesky*)pc->data; 161 PetscErrorCode ierr; 162 PetscTruth iascii,isstring; 163 164 PetscFunctionBegin; 165 ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr); 166 ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);CHKERRQ(ierr); 167 if (iascii) { 168 169 if (lu->inplace) {ierr = PetscViewerASCIIPrintf(viewer," Cholesky: in-place factorization\n");CHKERRQ(ierr);} 170 else {ierr = PetscViewerASCIIPrintf(viewer," Cholesky: out-of-place factorization\n");CHKERRQ(ierr);} 171 ierr = PetscViewerASCIIPrintf(viewer," matrix ordering: %s\n",lu->ordering);CHKERRQ(ierr); 172 if (lu->fact) { 173 ierr = PetscViewerASCIIPrintf(viewer," Cholesky: factor fill ratio needed %G\n",lu->actualfill);CHKERRQ(ierr); 174 ierr = PetscViewerASCIIPrintf(viewer," Factored matrix follows\n");CHKERRQ(ierr); 175 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 176 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 177 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 178 ierr = PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);CHKERRQ(ierr); 179 ierr = MatView(lu->fact,viewer);CHKERRQ(ierr); 180 ierr = PetscViewerPopFormat(viewer);CHKERRQ(ierr); 181 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 182 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 183 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 184 } 185 if (lu->reusefill) {ierr = PetscViewerASCIIPrintf(viewer," Reusing fill from past factorization\n");CHKERRQ(ierr);} 186 if (lu->reuseordering) {ierr = PetscViewerASCIIPrintf(viewer," Reusing reordering from past factorization\n");CHKERRQ(ierr);} 187 } else if (isstring) { 188 ierr = PetscViewerStringSPrintf(viewer," order=%s",lu->ordering);CHKERRQ(ierr);CHKERRQ(ierr); 189 } else { 190 SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PCCHOLESKY",((PetscObject)viewer)->type_name); 191 } 192 PetscFunctionReturn(0); 193 } 194 195 #undef __FUNCT__ 196 #define __FUNCT__ "PCFactorGetMatrix_Cholesky" 197 static PetscErrorCode PCFactorGetMatrix_Cholesky(PC pc,Mat *mat) 198 { 199 PC_Cholesky *dir = (PC_Cholesky*)pc->data; 200 201 PetscFunctionBegin; 202 if (!dir->fact) SETERRQ(PETSC_ERR_ORDER,"Matrix not yet factored; call after KSPSetUp() or PCSetUp()"); 203 *mat = dir->fact; 204 PetscFunctionReturn(0); 205 } 206 207 #undef __FUNCT__ 208 #define __FUNCT__ "PCSetUp_Cholesky" 209 static PetscErrorCode PCSetUp_Cholesky(PC pc) 210 { 211 PetscErrorCode ierr; 212 PetscTruth flg; 213 PC_Cholesky *dir = (PC_Cholesky*)pc->data; 214 215 PetscFunctionBegin; 216 if (dir->reusefill && pc->setupcalled) dir->info.fill = dir->actualfill; 217 218 if (dir->inplace) { 219 if (dir->row && dir->col && (dir->row != dir->col)) { 220 ierr = ISDestroy(dir->row);CHKERRQ(ierr); 221 dir->row = 0; 222 } 223 if (dir->col) { 224 ierr = ISDestroy(dir->col);CHKERRQ(ierr); 225 dir->col = 0; 226 } 227 ierr = MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);CHKERRQ(ierr); 228 if (dir->col && (dir->row != dir->col)) { /* only use row ordering for SBAIJ */ 229 ierr = ISDestroy(dir->col);CHKERRQ(ierr); 230 dir->col=0; 231 } 232 if (dir->row) {ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr);} 233 ierr = MatCholeskyFactor(pc->pmat,dir->row,&dir->info);CHKERRQ(ierr); 234 dir->fact = pc->pmat; 235 } else { 236 MatInfo info; 237 if (!pc->setupcalled) { 238 ierr = MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);CHKERRQ(ierr); 239 /* check if dir->row == dir->col */ 240 ierr = ISEqual(dir->row,dir->col,&flg);CHKERRQ(ierr); 241 if (!flg) SETERRQ(PETSC_ERR_ARG_INCOMP,"row and column permutations must equal"); 242 ierr = ISDestroy(dir->col);CHKERRQ(ierr); /* only pass one ordering into CholeskyFactor */ 243 dir->col=0; 244 245 ierr = PetscOptionsHasName(((PetscObject)pc)->prefix,"-pc_factor_nonzeros_along_diagonal",&flg);CHKERRQ(ierr); 246 if (flg) { 247 PetscReal tol = 1.e-10; 248 ierr = PetscOptionsGetReal(((PetscObject)pc)->prefix,"-pc_factor_nonzeros_along_diagonal",&tol,PETSC_NULL);CHKERRQ(ierr); 249 ierr = MatReorderForNonzeroDiagonal(pc->pmat,tol,dir->row,dir->row);CHKERRQ(ierr); 250 } 251 if (dir->row) {ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr);} 252 ierr = MatGetFactor(pc->pmat,dir->solvertype,MAT_FACTOR_CHOLESKY,&dir->fact);CHKERRQ(ierr); 253 ierr = MatCholeskyFactorSymbolic(pc->pmat,dir->row,&dir->info,&dir->fact);CHKERRQ(ierr); 254 ierr = MatGetInfo(dir->fact,MAT_LOCAL,&info);CHKERRQ(ierr); 255 dir->actualfill = info.fill_ratio_needed; 256 ierr = PetscLogObjectParent(pc,dir->fact);CHKERRQ(ierr); 257 } else if (pc->flag != SAME_NONZERO_PATTERN) { 258 if (!dir->reuseordering) { 259 if (dir->row && dir->col && (dir->row != dir->col)) { 260 ierr = ISDestroy(dir->row);CHKERRQ(ierr); 261 dir->row = 0; 262 } 263 if (dir->col) { 264 ierr = ISDestroy(dir->col);CHKERRQ(ierr); 265 dir->col =0; 266 } 267 ierr = MatGetOrdering(pc->pmat,dir->ordering,&dir->row,&dir->col);CHKERRQ(ierr); 268 if (dir->col && (dir->row != dir->col)) { /* only use row ordering for SBAIJ */ 269 ierr = ISDestroy(dir->col);CHKERRQ(ierr); 270 dir->col=0; 271 } 272 ierr = PetscOptionsHasName(((PetscObject)pc)->prefix,"-pc_factor_nonzeros_along_diagonal",&flg);CHKERRQ(ierr); 273 if (flg) { 274 PetscReal tol = 1.e-10; 275 ierr = PetscOptionsGetReal(((PetscObject)pc)->prefix,"-pc_factor_nonzeros_along_diagonal",&tol,PETSC_NULL);CHKERRQ(ierr); 276 ierr = MatReorderForNonzeroDiagonal(pc->pmat,tol,dir->row,dir->row);CHKERRQ(ierr); 277 } 278 if (dir->row) {ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr);} 279 } 280 ierr = MatDestroy(dir->fact);CHKERRQ(ierr); 281 ierr = MatGetFactor(pc->pmat,dir->solvertype,MAT_FACTOR_CHOLESKY,&dir->fact);CHKERRQ(ierr); 282 ierr = MatCholeskyFactorSymbolic(pc->pmat,dir->row,&dir->info,&dir->fact);CHKERRQ(ierr); 283 ierr = MatGetInfo(dir->fact,MAT_LOCAL,&info);CHKERRQ(ierr); 284 dir->actualfill = info.fill_ratio_needed; 285 ierr = PetscLogObjectParent(pc,dir->fact);CHKERRQ(ierr); 286 } 287 ierr = MatCholeskyFactorNumeric(pc->pmat,&dir->info,&dir->fact);CHKERRQ(ierr); 288 } 289 PetscFunctionReturn(0); 290 } 291 292 #undef __FUNCT__ 293 #define __FUNCT__ "PCDestroy_Cholesky" 294 static PetscErrorCode PCDestroy_Cholesky(PC pc) 295 { 296 PC_Cholesky *dir = (PC_Cholesky*)pc->data; 297 PetscErrorCode ierr; 298 299 PetscFunctionBegin; 300 if (!dir->inplace && dir->fact) {ierr = MatDestroy(dir->fact);CHKERRQ(ierr);} 301 if (dir->row) {ierr = ISDestroy(dir->row);CHKERRQ(ierr);} 302 if (dir->col) {ierr = ISDestroy(dir->col);CHKERRQ(ierr);} 303 ierr = PetscStrfree(dir->ordering);CHKERRQ(ierr); 304 ierr = PetscStrfree(dir->solvertype);CHKERRQ(ierr); 305 ierr = PetscFree(dir);CHKERRQ(ierr); 306 PetscFunctionReturn(0); 307 } 308 309 #undef __FUNCT__ 310 #define __FUNCT__ "PCApply_Cholesky" 311 static PetscErrorCode PCApply_Cholesky(PC pc,Vec x,Vec y) 312 { 313 PC_Cholesky *dir = (PC_Cholesky*)pc->data; 314 PetscErrorCode ierr; 315 316 PetscFunctionBegin; 317 if (dir->inplace) {ierr = MatSolve(pc->pmat,x,y);CHKERRQ(ierr);} 318 else {ierr = MatSolve(dir->fact,x,y);CHKERRQ(ierr);} 319 PetscFunctionReturn(0); 320 } 321 322 #undef __FUNCT__ 323 #define __FUNCT__ "PCApplyTranspose_Cholesky" 324 static PetscErrorCode PCApplyTranspose_Cholesky(PC pc,Vec x,Vec y) 325 { 326 PC_Cholesky *dir = (PC_Cholesky*)pc->data; 327 PetscErrorCode ierr; 328 329 PetscFunctionBegin; 330 if (dir->inplace) {ierr = MatSolveTranspose(pc->pmat,x,y);CHKERRQ(ierr);} 331 else {ierr = MatSolveTranspose(dir->fact,x,y);CHKERRQ(ierr);} 332 PetscFunctionReturn(0); 333 } 334 335 /* -----------------------------------------------------------------------------------*/ 336 337 EXTERN_C_BEGIN 338 #undef __FUNCT__ 339 #define __FUNCT__ "PCFactorSetFill_Cholesky" 340 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetFill_Cholesky(PC pc,PetscReal fill) 341 { 342 PC_Cholesky *dir; 343 344 PetscFunctionBegin; 345 dir = (PC_Cholesky*)pc->data; 346 dir->info.fill = fill; 347 PetscFunctionReturn(0); 348 } 349 EXTERN_C_END 350 351 EXTERN_C_BEGIN 352 #undef __FUNCT__ 353 #define __FUNCT__ "PCFactorSetUseInPlace_Cholesky" 354 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetUseInPlace_Cholesky(PC pc) 355 { 356 PC_Cholesky *dir; 357 358 PetscFunctionBegin; 359 dir = (PC_Cholesky*)pc->data; 360 dir->inplace = PETSC_TRUE; 361 PetscFunctionReturn(0); 362 } 363 EXTERN_C_END 364 365 EXTERN_C_BEGIN 366 #undef __FUNCT__ 367 #define __FUNCT__ "PCFactorSetMatOrderingType_Cholesky" 368 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetMatOrderingType_Cholesky(PC pc,MatOrderingType ordering) 369 { 370 PC_Cholesky *dir = (PC_Cholesky*)pc->data; 371 PetscErrorCode ierr; 372 373 PetscFunctionBegin; 374 ierr = PetscStrfree(dir->ordering);CHKERRQ(ierr); 375 ierr = PetscStrallocpy(ordering,&dir->ordering);CHKERRQ(ierr); 376 PetscFunctionReturn(0); 377 } 378 EXTERN_C_END 379 380 /* -----------------------------------------------------------------------------------*/ 381 382 #undef __FUNCT__ 383 #define __FUNCT__ "PCFactorSetReuseOrdering" 384 /*@ 385 PCFactorSetReuseOrdering - When similar matrices are factored, this 386 causes the ordering computed in the first factor to be used for all 387 following factors. 388 389 Collective on PC 390 391 Input Parameters: 392 + pc - the preconditioner context 393 - flag - PETSC_TRUE to reuse else PETSC_FALSE 394 395 Options Database Key: 396 . -pc_factor_reuse_ordering - Activate PCFactorSetReuseOrdering() 397 398 Level: intermediate 399 400 .keywords: PC, levels, reordering, factorization, incomplete, LU 401 402 .seealso: PCFactorSetReuseFill() 403 @*/ 404 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseOrdering(PC pc,PetscTruth flag) 405 { 406 PetscErrorCode ierr,(*f)(PC,PetscTruth); 407 408 PetscFunctionBegin; 409 PetscValidHeaderSpecific(pc,PC_COOKIE,1); 410 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetReuseOrdering_C",(void (**)(void))&f);CHKERRQ(ierr); 411 if (f) { 412 ierr = (*f)(pc,flag);CHKERRQ(ierr); 413 } 414 PetscFunctionReturn(0); 415 } 416 417 #undef __FUNCT__ 418 #define __FUNCT__ "PCFactorSetReuseFill" 419 /*@ 420 PCFactorSetReuseFill - When matrices with same different nonzero structure are factored, 421 this causes later ones to use the fill ratio computed in the initial factorization. 422 423 Collective on PC 424 425 Input Parameters: 426 + pc - the preconditioner context 427 - flag - PETSC_TRUE to reuse else PETSC_FALSE 428 429 Options Database Key: 430 . -pc_factor_reuse_fill - Activates PCFactorSetReuseFill() 431 432 Level: intermediate 433 434 .keywords: PC, levels, reordering, factorization, incomplete, Cholesky 435 436 .seealso: PCFactorSetReuseOrdering() 437 @*/ 438 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseFill(PC pc,PetscTruth flag) 439 { 440 PetscErrorCode ierr,(*f)(PC,PetscTruth); 441 442 PetscFunctionBegin; 443 PetscValidHeaderSpecific(pc,PC_COOKIE,2); 444 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetReuseFill_C",(void (**)(void))&f);CHKERRQ(ierr); 445 if (f) { 446 ierr = (*f)(pc,flag);CHKERRQ(ierr); 447 } 448 PetscFunctionReturn(0); 449 } 450 451 /*MC 452 PCCHOLESKY - Uses a direct solver, based on Cholesky factorization, as a preconditioner 453 454 Options Database Keys: 455 + -pc_factor_reuse_ordering - Activate PCFactorSetReuseOrdering() 456 . -pc_factor_reuse_fill - Activates PCFactorSetReuseFill() 457 . -pc_factor_fill <fill> - Sets fill amount 458 . -pc_factor_in_place - Activates in-place factorization 459 . -pc_factor_mat_ordering_type <nd,rcm,...> - Sets ordering routine 460 . -pc_factor_shift_nonzero <shift> - Sets shift amount or PETSC_DECIDE for the default 461 - -pc_factor_shift_positive_definite [PETSC_TRUE/PETSC_FALSE] - Activate/Deactivate PCFactorSetShiftPd(); the value 462 is optional with PETSC_TRUE being the default 463 464 Notes: Not all options work for all matrix formats 465 466 Level: beginner 467 468 Concepts: Cholesky factorization, direct solver 469 470 Notes: Usually this will compute an "exact" solution in one iteration and does 471 not need a Krylov method (i.e. you can use -ksp_type preonly, or 472 KSPSetType(ksp,KSPPREONLY) for the Krylov method 473 474 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 475 PCILU, PCLU, PCICC, PCFactorSetReuseOrdering(), PCFactorSetReuseFill(), PCFactorGetMatrix(), 476 PCFactorSetFill(), PCFactorSetShiftNonzero(), PCFactorSetShiftPd(), 477 PCFactorSetUseInPlace(), PCFactorSetMatOrderingType(),PCFactorSetShiftNonzero(),PCFactorSetShiftPd() 478 479 M*/ 480 481 EXTERN_C_BEGIN 482 #undef __FUNCT__ 483 #define __FUNCT__ "PCCreate_Cholesky" 484 PetscErrorCode PETSCKSP_DLLEXPORT PCCreate_Cholesky(PC pc) 485 { 486 PetscErrorCode ierr; 487 PC_Cholesky *dir; 488 489 PetscFunctionBegin; 490 ierr = PetscNewLog(pc,PC_Cholesky,&dir);CHKERRQ(ierr); 491 492 dir->fact = 0; 493 dir->inplace = PETSC_FALSE; 494 ierr = MatFactorInfoInitialize(&dir->info);CHKERRQ(ierr); 495 dir->info.fill = 5.0; 496 dir->info.shiftnz = 0.0; 497 dir->info.shiftpd = 0.0; /* false */ 498 dir->info.shift_fraction = 0.0; 499 dir->info.pivotinblocks = 1.0; 500 dir->col = 0; 501 dir->row = 0; 502 ierr = PetscStrallocpy(MATORDERING_NATURAL,&dir->ordering);CHKERRQ(ierr); 503 ierr = PetscStrallocpy("petsc",&dir->solvertype);CHKERRQ(ierr); 504 dir->reusefill = PETSC_FALSE; 505 dir->reuseordering = PETSC_FALSE; 506 pc->data = (void*)dir; 507 508 pc->ops->destroy = PCDestroy_Cholesky; 509 pc->ops->apply = PCApply_Cholesky; 510 pc->ops->applytranspose = PCApplyTranspose_Cholesky; 511 pc->ops->setup = PCSetUp_Cholesky; 512 pc->ops->setfromoptions = PCSetFromOptions_Cholesky; 513 pc->ops->view = PCView_Cholesky; 514 pc->ops->applyrichardson = 0; 515 pc->ops->getfactoredmatrix = PCFactorGetMatrix_Cholesky; 516 517 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetZeroPivot_C","PCFactorSetZeroPivot_Cholesky", 518 PCFactorSetZeroPivot_Cholesky);CHKERRQ(ierr); 519 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftNonzero_C","PCFactorSetShiftNonzero_Cholesky", 520 PCFactorSetShiftNonzero_Cholesky);CHKERRQ(ierr); 521 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftPd_C","PCFactorSetShiftPd_Cholesky", 522 PCFactorSetShiftPd_Cholesky);CHKERRQ(ierr); 523 524 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetFill_C","PCFactorSetFill_Cholesky", 525 PCFactorSetFill_Cholesky);CHKERRQ(ierr); 526 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetUseInPlace_C","PCFactorSetUseInPlace_Cholesky", 527 PCFactorSetUseInPlace_Cholesky);CHKERRQ(ierr); 528 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetMatOrderingType_C","PCFactorSetMatOrderingType_Cholesky", 529 PCFactorSetMatOrderingType_Cholesky);CHKERRQ(ierr); 530 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseOrdering_C","PCFactorSetReuseOrdering_Cholesky", 531 PCFactorSetReuseOrdering_Cholesky);CHKERRQ(ierr); 532 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseFill_C","PCFactorSetReuseFill_Cholesky", 533 PCFactorSetReuseFill_Cholesky);CHKERRQ(ierr); 534 PetscFunctionReturn(0); 535 } 536 EXTERN_C_END 537