1 #define PETSCKSP_DLL 2 3 /* 4 Defines a Cholesky factorization preconditioner for any Mat implementation. 5 Presently only provided for MPIRowbs format (i.e. BlockSolve). 6 */ 7 8 #include "src/ksp/pc/impls/factor/icc/icc.h" /*I "petscpc.h" I*/ 9 10 EXTERN_C_BEGIN 11 #undef __FUNCT__ 12 #define __FUNCT__ "PCFactorSetZeroPivot_ICC" 13 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetZeroPivot_ICC(PC pc,PetscReal z) 14 { 15 PC_ICC *icc; 16 17 PetscFunctionBegin; 18 icc = (PC_ICC*)pc->data; 19 icc->info.zeropivot = z; 20 PetscFunctionReturn(0); 21 } 22 EXTERN_C_END 23 24 EXTERN_C_BEGIN 25 #undef __FUNCT__ 26 #define __FUNCT__ "PCFactorSetShiftNonzero_ICC" 27 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftNonzero_ICC(PC pc,PetscReal shift) 28 { 29 PC_ICC *dir; 30 31 PetscFunctionBegin; 32 dir = (PC_ICC*)pc->data; 33 if (shift == (PetscReal) PETSC_DECIDE) { 34 dir->info.shiftnz = 1.e-12; 35 } else { 36 dir->info.shiftnz = shift; 37 } 38 PetscFunctionReturn(0); 39 } 40 EXTERN_C_END 41 42 EXTERN_C_BEGIN 43 #undef __FUNCT__ 44 #define __FUNCT__ "PCFactorSetShiftPd_ICC" 45 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetShiftPd_ICC(PC pc,PetscTruth shift) 46 { 47 PC_ICC *dir; 48 49 PetscFunctionBegin; 50 dir = (PC_ICC*)pc->data; 51 if (shift) { 52 dir->info.shift_fraction = 0.0; 53 dir->info.shiftpd = 1.0; 54 } else { 55 dir->info.shiftpd = 0.0; 56 } 57 PetscFunctionReturn(0); 58 } 59 EXTERN_C_END 60 61 EXTERN_C_BEGIN 62 #undef __FUNCT__ 63 #define __FUNCT__ "PCICCSetMatOrdering_ICC" 64 PetscErrorCode PETSCKSP_DLLEXPORT PCICCSetMatOrdering_ICC(PC pc,MatOrderingType ordering) 65 { 66 PC_ICC *dir = (PC_ICC*)pc->data; 67 PetscErrorCode ierr; 68 69 PetscFunctionBegin; 70 ierr = PetscStrfree(dir->ordering);CHKERRQ(ierr); 71 ierr = PetscStrallocpy(ordering,&dir->ordering);CHKERRQ(ierr); 72 PetscFunctionReturn(0); 73 } 74 EXTERN_C_END 75 76 EXTERN_C_BEGIN 77 #undef __FUNCT__ 78 #define __FUNCT__ "PCFactorSetFill_ICC" 79 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetFill_ICC(PC pc,PetscReal fill) 80 { 81 PC_ICC *dir; 82 83 PetscFunctionBegin; 84 dir = (PC_ICC*)pc->data; 85 dir->info.fill = fill; 86 PetscFunctionReturn(0); 87 } 88 EXTERN_C_END 89 90 EXTERN_C_BEGIN 91 #undef __FUNCT__ 92 #define __FUNCT__ "PCICCSetLevels_ICC" 93 PetscErrorCode PETSCKSP_DLLEXPORT PCICCSetLevels_ICC(PC pc,PetscInt levels) 94 { 95 PC_ICC *icc; 96 97 PetscFunctionBegin; 98 icc = (PC_ICC*)pc->data; 99 icc->info.levels = levels; 100 PetscFunctionReturn(0); 101 } 102 EXTERN_C_END 103 104 #undef __FUNCT__ 105 #define __FUNCT__ "PCICCSetMatOrdering" 106 /*@ 107 PCICCSetMatOrdering - Sets the ordering routine (to reduce fill) to 108 be used it the ICC factorization. 109 110 Collective on PC 111 112 Input Parameters: 113 + pc - the preconditioner context 114 - ordering - the matrix ordering name, for example, MATORDERING_ND or MATORDERING_RCM 115 116 Options Database Key: 117 . -pc_icc_mat_ordering_type <nd,rcm,...> - Sets ordering routine 118 119 Level: intermediate 120 121 .seealso: PCLUSetMatOrdering() 122 123 .keywords: PC, ICC, set, matrix, reordering 124 125 @*/ 126 PetscErrorCode PETSCKSP_DLLEXPORT PCICCSetMatOrdering(PC pc,MatOrderingType ordering) 127 { 128 PetscErrorCode ierr,(*f)(PC,MatOrderingType); 129 130 PetscFunctionBegin; 131 PetscValidHeaderSpecific(pc,PC_COOKIE,1); 132 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCICCSetMatOrdering_C",(void (**)(void))&f);CHKERRQ(ierr); 133 if (f) { 134 ierr = (*f)(pc,ordering);CHKERRQ(ierr); 135 } 136 PetscFunctionReturn(0); 137 } 138 139 #undef __FUNCT__ 140 #define __FUNCT__ "PCICCSetLevels" 141 /*@ 142 PCICCSetLevels - Sets the number of levels of fill to use. 143 144 Collective on PC 145 146 Input Parameters: 147 + pc - the preconditioner context 148 - levels - number of levels of fill 149 150 Options Database Key: 151 . -pc_icc_levels <levels> - Sets fill level 152 153 Level: intermediate 154 155 Concepts: ICC^setting levels of fill 156 157 @*/ 158 PetscErrorCode PETSCKSP_DLLEXPORT PCICCSetLevels(PC pc,PetscInt levels) 159 { 160 PetscErrorCode ierr,(*f)(PC,PetscInt); 161 162 PetscFunctionBegin; 163 PetscValidHeaderSpecific(pc,PC_COOKIE,1); 164 if (levels < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"negative levels"); 165 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCICCSetLevels_C",(void (**)(void))&f);CHKERRQ(ierr); 166 if (f) { 167 ierr = (*f)(pc,levels);CHKERRQ(ierr); 168 } 169 PetscFunctionReturn(0); 170 } 171 172 #undef __FUNCT__ 173 #define __FUNCT__ "PCSetup_ICC" 174 static PetscErrorCode PCSetup_ICC(PC pc) 175 { 176 PC_ICC *icc = (PC_ICC*)pc->data; 177 IS perm,cperm; 178 PetscErrorCode ierr; 179 180 PetscFunctionBegin; 181 ierr = MatGetOrdering(pc->pmat,icc->ordering,&perm,&cperm);CHKERRQ(ierr); 182 183 if (!pc->setupcalled) { 184 ierr = MatICCFactorSymbolic(pc->pmat,perm,&icc->info,&icc->fact);CHKERRQ(ierr); 185 } else if (pc->flag != SAME_NONZERO_PATTERN) { 186 ierr = MatDestroy(icc->fact);CHKERRQ(ierr); 187 ierr = MatICCFactorSymbolic(pc->pmat,perm,&icc->info,&icc->fact);CHKERRQ(ierr); 188 } 189 ierr = ISDestroy(cperm);CHKERRQ(ierr); 190 ierr = ISDestroy(perm);CHKERRQ(ierr); 191 ierr = MatCholeskyFactorNumeric(pc->pmat,&icc->info,&icc->fact);CHKERRQ(ierr); 192 PetscFunctionReturn(0); 193 } 194 195 #undef __FUNCT__ 196 #define __FUNCT__ "PCDestroy_ICC" 197 static PetscErrorCode PCDestroy_ICC(PC pc) 198 { 199 PC_ICC *icc = (PC_ICC*)pc->data; 200 PetscErrorCode ierr; 201 202 PetscFunctionBegin; 203 if (icc->fact) {ierr = MatDestroy(icc->fact);CHKERRQ(ierr);} 204 ierr = PetscStrfree(icc->ordering);CHKERRQ(ierr); 205 ierr = PetscFree(icc);CHKERRQ(ierr); 206 PetscFunctionReturn(0); 207 } 208 209 #undef __FUNCT__ 210 #define __FUNCT__ "PCApply_ICC" 211 static PetscErrorCode PCApply_ICC(PC pc,Vec x,Vec y) 212 { 213 PC_ICC *icc = (PC_ICC*)pc->data; 214 PetscErrorCode ierr; 215 216 PetscFunctionBegin; 217 ierr = MatSolve(icc->fact,x,y);CHKERRQ(ierr); 218 PetscFunctionReturn(0); 219 } 220 221 #undef __FUNCT__ 222 #define __FUNCT__ "PCApplySymmetricLeft_ICC" 223 static PetscErrorCode PCApplySymmetricLeft_ICC(PC pc,Vec x,Vec y) 224 { 225 PetscErrorCode ierr; 226 PC_ICC *icc = (PC_ICC*)pc->data; 227 228 PetscFunctionBegin; 229 ierr = MatForwardSolve(icc->fact,x,y);CHKERRQ(ierr); 230 PetscFunctionReturn(0); 231 } 232 233 #undef __FUNCT__ 234 #define __FUNCT__ "PCApplySymmetricRight_ICC" 235 static PetscErrorCode PCApplySymmetricRight_ICC(PC pc,Vec x,Vec y) 236 { 237 PetscErrorCode ierr; 238 PC_ICC *icc = (PC_ICC*)pc->data; 239 240 PetscFunctionBegin; 241 ierr = MatBackwardSolve(icc->fact,x,y);CHKERRQ(ierr); 242 PetscFunctionReturn(0); 243 } 244 245 #undef __FUNCT__ 246 #define __FUNCT__ "PCGetFactoredMatrix_ICC" 247 static PetscErrorCode PCGetFactoredMatrix_ICC(PC pc,Mat *mat) 248 { 249 PC_ICC *icc = (PC_ICC*)pc->data; 250 251 PetscFunctionBegin; 252 *mat = icc->fact; 253 PetscFunctionReturn(0); 254 } 255 256 #undef __FUNCT__ 257 #define __FUNCT__ "PCSetFromOptions_ICC" 258 static PetscErrorCode PCSetFromOptions_ICC(PC pc) 259 { 260 PC_ICC *icc = (PC_ICC*)pc->data; 261 char tname[256]; 262 PetscTruth flg; 263 PetscErrorCode ierr; 264 PetscFList ordlist; 265 266 PetscFunctionBegin; 267 ierr = MatOrderingRegisterAll(PETSC_NULL);CHKERRQ(ierr); 268 ierr = PetscOptionsHead("ICC Options");CHKERRQ(ierr); 269 ierr = PetscOptionsReal("-pc_icc_levels","levels of fill","PCICCSetLevels",icc->info.levels,&icc->info.levels,&flg);CHKERRQ(ierr); 270 ierr = PetscOptionsReal("-pc_factor_fill","Expected fill in factorization","PCFactorSetFill",icc->info.fill,&icc->info.fill,&flg);CHKERRQ(ierr); 271 ierr = MatGetOrderingList(&ordlist);CHKERRQ(ierr); 272 ierr = PetscOptionsList("-pc_icc_mat_ordering_type","Reorder to reduce nonzeros in ICC","PCICCSetMatOrdering",ordlist,icc->ordering,tname,256,&flg);CHKERRQ(ierr); 273 if (flg) { 274 ierr = PCICCSetMatOrdering(pc,tname);CHKERRQ(ierr); 275 } 276 ierr = PetscOptionsName("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",&flg);CHKERRQ(ierr); 277 if (flg) { 278 ierr = PCFactorSetShiftNonzero(pc,(PetscReal)PETSC_DECIDE);CHKERRQ(ierr); 279 } 280 ierr = PetscOptionsReal("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",icc->info.shiftnz,&icc->info.shiftnz,0);CHKERRQ(ierr); 281 ierr = PetscOptionsName("-pc_factor_shift_positive_definite","Manteuffel shift applied to diagonal","PCICCSetShift",&flg);CHKERRQ(ierr); 282 if (flg) { 283 ierr = PCFactorSetShiftPd(pc,PETSC_TRUE);CHKERRQ(ierr); 284 } else { 285 ierr = PCFactorSetShiftPd(pc,PETSC_FALSE);CHKERRQ(ierr); 286 } 287 ierr = PetscOptionsReal("-pc_factor_zeropivot","Pivot is considered zero if less than","PCFactorSetZeroPivot",icc->info.zeropivot,&icc->info.zeropivot,0);CHKERRQ(ierr); 288 289 ierr = PetscOptionsTail();CHKERRQ(ierr); 290 PetscFunctionReturn(0); 291 } 292 293 #undef __FUNCT__ 294 #define __FUNCT__ "PCView_ICC" 295 static PetscErrorCode PCView_ICC(PC pc,PetscViewer viewer) 296 { 297 PC_ICC *icc = (PC_ICC*)pc->data; 298 PetscErrorCode ierr; 299 PetscTruth isstring,iascii; 300 301 PetscFunctionBegin; 302 ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);CHKERRQ(ierr); 303 ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr); 304 if (iascii) { 305 if (icc->info.levels == 1) { 306 ierr = PetscViewerASCIIPrintf(viewer," ICC: %D level of fill\n",(PetscInt)icc->info.levels);CHKERRQ(ierr); 307 } else { 308 ierr = PetscViewerASCIIPrintf(viewer," ICC: %D levels of fill\n",(PetscInt)icc->info.levels);CHKERRQ(ierr); 309 } 310 ierr = PetscViewerASCIIPrintf(viewer," ICC: max fill ratio allocated %G\n",icc->info.fill);CHKERRQ(ierr); 311 if (icc->info.shiftpd) {ierr = PetscViewerASCIIPrintf(viewer," ICC: using Manteuffel shift\n");CHKERRQ(ierr);} 312 } else if (isstring) { 313 ierr = PetscViewerStringSPrintf(viewer," lvls=%D",(PetscInt)icc->info.levels);CHKERRQ(ierr);CHKERRQ(ierr); 314 } else { 315 SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PCICC",((PetscObject)viewer)->type_name); 316 } 317 PetscFunctionReturn(0); 318 } 319 320 /*MC 321 PCICC - Incomplete Cholesky factorization preconditioners. 322 323 Options Database Keys: 324 + -pc_icc_levels <k> - number of levels of fill for ICC(k) 325 . -pc_icc_in_place - only for ICC(0) with natural ordering, reuses the space of the matrix for 326 its factorization (overwrites original matrix) 327 . -pc_factor_fill <nfill> - expected amount of fill in factored matrix compared to original matrix, nfill > 1 328 . -pc_icc_mat_ordering_type <natural,nd,1wd,rcm,qmd> - set the row/column ordering of the factored matrix 329 . -pc_factor_shift_nonzero <shift> - Sets shift amount or PETSC_DECIDE for the default 330 - -pc_factor_shift_positive_definite [PETSC_TRUE/PETSC_FALSE] - Activate/Deactivate PCFactorSetShiftPd(); the value 331 is optional with PETSC_TRUE being the default 332 333 Level: beginner 334 335 Concepts: incomplete Cholesky factorization 336 337 Notes: Only implemented for some matrix formats. Not implemented in parallel 338 339 For BAIJ matrices this implements a point block ICC. 340 341 The Manteuffel shift is only implemented for matrices with block size 1 342 343 By default, the Manteuffel is applied (for matrices with block size 1). Call PCICCSetShift(pc,PETSC_FALSE); 344 to turn off the shift. 345 346 347 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, PCSOR, MatOrderingType, 348 PCFactorSetZeroPivot(), PCFactorSetShiftNonzero(), PCFactorSetShiftPd(), 349 PCFactorSetFill(), PCICCSetMatOrdering(), PCICCSetReuseOrdering(), 350 PCICCSetLevels(),PCFactorSetShiftNonzero(),PCFactorSetShiftPd(), 351 352 M*/ 353 354 EXTERN_C_BEGIN 355 #undef __FUNCT__ 356 #define __FUNCT__ "PCCreate_ICC" 357 PetscErrorCode PETSCKSP_DLLEXPORT PCCreate_ICC(PC pc) 358 { 359 PetscErrorCode ierr; 360 PC_ICC *icc; 361 362 PetscFunctionBegin; 363 ierr = PetscNew(PC_ICC,&icc);CHKERRQ(ierr); 364 ierr = PetscLogObjectMemory(pc,sizeof(PC_ICC));CHKERRQ(ierr); 365 366 icc->fact = 0; 367 ierr = PetscStrallocpy(MATORDERING_NATURAL,&icc->ordering);CHKERRQ(ierr); 368 ierr = MatFactorInfoInitialize(&icc->info);CHKERRQ(ierr); 369 icc->info.levels = 0; 370 icc->info.fill = 1.0; 371 icc->implctx = 0; 372 373 icc->info.dtcol = PETSC_DEFAULT; 374 icc->info.shiftnz = 0.0; 375 icc->info.shiftpd = 1.0; /* true */ 376 icc->info.shift_fraction = 0.0; 377 icc->info.zeropivot = 1.e-12; 378 pc->data = (void*)icc; 379 380 pc->ops->apply = PCApply_ICC; 381 pc->ops->setup = PCSetup_ICC; 382 pc->ops->destroy = PCDestroy_ICC; 383 pc->ops->setfromoptions = PCSetFromOptions_ICC; 384 pc->ops->view = PCView_ICC; 385 pc->ops->getfactoredmatrix = PCGetFactoredMatrix_ICC; 386 pc->ops->applysymmetricleft = PCApplySymmetricLeft_ICC; 387 pc->ops->applysymmetricright = PCApplySymmetricRight_ICC; 388 389 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetZeroPivot_C","PCFactorSetZeroPivot_ICC", 390 PCFactorSetZeroPivot_ICC);CHKERRQ(ierr); 391 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftNonzero_C","PCFactorSetShiftNonzero_ICC", 392 PCFactorSetShiftNonzero_ICC);CHKERRQ(ierr); 393 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftPd_C","PCFactorSetShiftPd_ICC", 394 PCFactorSetShiftPd_ICC);CHKERRQ(ierr); 395 396 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCICCSetLevels_C","PCICCSetLevels_ICC", 397 PCICCSetLevels_ICC);CHKERRQ(ierr); 398 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetFill_C","PCFactorSetFill_ICC", 399 PCFactorSetFill_ICC);CHKERRQ(ierr); 400 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCICCSetMatOrdering_C","PCICCSetMatOrdering_ICC", 401 PCICCSetMatOrdering_ICC);CHKERRQ(ierr); 402 PetscFunctionReturn(0); 403 } 404 EXTERN_C_END 405 406 407