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__ "PCFactorSetMatOrdering_ICC" 64 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetMatOrdering_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__ "PCFactorSetLevels_ICC" 93 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetLevels_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__ "PCSetup_ICC" 106 static PetscErrorCode PCSetup_ICC(PC pc) 107 { 108 PC_ICC *icc = (PC_ICC*)pc->data; 109 IS perm,cperm; 110 PetscErrorCode ierr; 111 112 PetscFunctionBegin; 113 ierr = MatGetOrdering(pc->pmat,icc->ordering,&perm,&cperm);CHKERRQ(ierr); 114 115 if (!pc->setupcalled) { 116 ierr = MatICCFactorSymbolic(pc->pmat,perm,&icc->info,&icc->fact);CHKERRQ(ierr); 117 } else if (pc->flag != SAME_NONZERO_PATTERN) { 118 ierr = MatDestroy(icc->fact);CHKERRQ(ierr); 119 ierr = MatICCFactorSymbolic(pc->pmat,perm,&icc->info,&icc->fact);CHKERRQ(ierr); 120 } 121 ierr = ISDestroy(cperm);CHKERRQ(ierr); 122 ierr = ISDestroy(perm);CHKERRQ(ierr); 123 ierr = MatCholeskyFactorNumeric(pc->pmat,&icc->info,&icc->fact);CHKERRQ(ierr); 124 PetscFunctionReturn(0); 125 } 126 127 #undef __FUNCT__ 128 #define __FUNCT__ "PCDestroy_ICC" 129 static PetscErrorCode PCDestroy_ICC(PC pc) 130 { 131 PC_ICC *icc = (PC_ICC*)pc->data; 132 PetscErrorCode ierr; 133 134 PetscFunctionBegin; 135 if (icc->fact) {ierr = MatDestroy(icc->fact);CHKERRQ(ierr);} 136 ierr = PetscStrfree(icc->ordering);CHKERRQ(ierr); 137 ierr = PetscFree(icc);CHKERRQ(ierr); 138 PetscFunctionReturn(0); 139 } 140 141 #undef __FUNCT__ 142 #define __FUNCT__ "PCApply_ICC" 143 static PetscErrorCode PCApply_ICC(PC pc,Vec x,Vec y) 144 { 145 PC_ICC *icc = (PC_ICC*)pc->data; 146 PetscErrorCode ierr; 147 148 PetscFunctionBegin; 149 ierr = MatSolve(icc->fact,x,y);CHKERRQ(ierr); 150 PetscFunctionReturn(0); 151 } 152 153 #undef __FUNCT__ 154 #define __FUNCT__ "PCApplySymmetricLeft_ICC" 155 static PetscErrorCode PCApplySymmetricLeft_ICC(PC pc,Vec x,Vec y) 156 { 157 PetscErrorCode ierr; 158 PC_ICC *icc = (PC_ICC*)pc->data; 159 160 PetscFunctionBegin; 161 ierr = MatForwardSolve(icc->fact,x,y);CHKERRQ(ierr); 162 PetscFunctionReturn(0); 163 } 164 165 #undef __FUNCT__ 166 #define __FUNCT__ "PCApplySymmetricRight_ICC" 167 static PetscErrorCode PCApplySymmetricRight_ICC(PC pc,Vec x,Vec y) 168 { 169 PetscErrorCode ierr; 170 PC_ICC *icc = (PC_ICC*)pc->data; 171 172 PetscFunctionBegin; 173 ierr = MatBackwardSolve(icc->fact,x,y);CHKERRQ(ierr); 174 PetscFunctionReturn(0); 175 } 176 177 #undef __FUNCT__ 178 #define __FUNCT__ "PCGetFactoredMatrix_ICC" 179 static PetscErrorCode PCGetFactoredMatrix_ICC(PC pc,Mat *mat) 180 { 181 PC_ICC *icc = (PC_ICC*)pc->data; 182 183 PetscFunctionBegin; 184 *mat = icc->fact; 185 PetscFunctionReturn(0); 186 } 187 188 #undef __FUNCT__ 189 #define __FUNCT__ "PCSetFromOptions_ICC" 190 static PetscErrorCode PCSetFromOptions_ICC(PC pc) 191 { 192 PC_ICC *icc = (PC_ICC*)pc->data; 193 char tname[256]; 194 PetscTruth flg; 195 PetscErrorCode ierr; 196 PetscFList ordlist; 197 198 PetscFunctionBegin; 199 ierr = MatOrderingRegisterAll(PETSC_NULL);CHKERRQ(ierr); 200 ierr = PetscOptionsHead("ICC Options");CHKERRQ(ierr); 201 ierr = PetscOptionsReal("-pc_factor_levels","levels of fill","PCFactorSetLevels",icc->info.levels,&icc->info.levels,&flg);CHKERRQ(ierr); 202 ierr = PetscOptionsReal("-pc_factor_fill","Expected fill in factorization","PCFactorSetFill",icc->info.fill,&icc->info.fill,&flg);CHKERRQ(ierr); 203 ierr = MatGetOrderingList(&ordlist);CHKERRQ(ierr); 204 ierr = PetscOptionsList("-pc_factor_mat_ordering_type","Reorder to reduce nonzeros in ICC","PCFactorSetMatOrdering",ordlist,icc->ordering,tname,256,&flg);CHKERRQ(ierr); 205 if (flg) { 206 ierr = PCFactorSetMatOrdering(pc,tname);CHKERRQ(ierr); 207 } 208 ierr = PetscOptionsName("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",&flg);CHKERRQ(ierr); 209 if (flg) { 210 ierr = PCFactorSetShiftNonzero(pc,(PetscReal)PETSC_DECIDE);CHKERRQ(ierr); 211 } 212 ierr = PetscOptionsReal("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",icc->info.shiftnz,&icc->info.shiftnz,0);CHKERRQ(ierr); 213 ierr = PetscOptionsName("-pc_factor_shift_positive_definite","Manteuffel shift applied to diagonal","PCFactorSetShift",&flg);CHKERRQ(ierr); 214 if (flg) { 215 ierr = PCFactorSetShiftPd(pc,PETSC_TRUE);CHKERRQ(ierr); 216 } else { 217 ierr = PCFactorSetShiftPd(pc,PETSC_FALSE);CHKERRQ(ierr); 218 } 219 ierr = PetscOptionsReal("-pc_factor_zeropivot","Pivot is considered zero if less than","PCFactorSetZeroPivot",icc->info.zeropivot,&icc->info.zeropivot,0);CHKERRQ(ierr); 220 221 ierr = PetscOptionsTail();CHKERRQ(ierr); 222 PetscFunctionReturn(0); 223 } 224 225 #undef __FUNCT__ 226 #define __FUNCT__ "PCView_ICC" 227 static PetscErrorCode PCView_ICC(PC pc,PetscViewer viewer) 228 { 229 PC_ICC *icc = (PC_ICC*)pc->data; 230 PetscErrorCode ierr; 231 PetscTruth isstring,iascii; 232 233 PetscFunctionBegin; 234 ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);CHKERRQ(ierr); 235 ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr); 236 if (iascii) { 237 if (icc->info.levels == 1) { 238 ierr = PetscViewerASCIIPrintf(viewer," ICC: %D level of fill\n",(PetscInt)icc->info.levels);CHKERRQ(ierr); 239 } else { 240 ierr = PetscViewerASCIIPrintf(viewer," ICC: %D levels of fill\n",(PetscInt)icc->info.levels);CHKERRQ(ierr); 241 } 242 ierr = PetscViewerASCIIPrintf(viewer," ICC: max fill ratio allocated %G\n",icc->info.fill);CHKERRQ(ierr); 243 if (icc->info.shiftpd) {ierr = PetscViewerASCIIPrintf(viewer," ICC: using Manteuffel shift\n");CHKERRQ(ierr);} 244 } else if (isstring) { 245 ierr = PetscViewerStringSPrintf(viewer," lvls=%D",(PetscInt)icc->info.levels);CHKERRQ(ierr);CHKERRQ(ierr); 246 } else { 247 SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PCICC",((PetscObject)viewer)->type_name); 248 } 249 PetscFunctionReturn(0); 250 } 251 252 /*MC 253 PCICC - Incomplete Cholesky factorization preconditioners. 254 255 Options Database Keys: 256 + -pc_factor_levels <k> - number of levels of fill for ICC(k) 257 . -pc_factor_in_place - only for ICC(0) with natural ordering, reuses the space of the matrix for 258 its factorization (overwrites original matrix) 259 . -pc_factor_fill <nfill> - expected amount of fill in factored matrix compared to original matrix, nfill > 1 260 . -pc_factor_mat_ordering_type <natural,nd,1wd,rcm,qmd> - set the row/column ordering of the factored matrix 261 . -pc_factor_shift_nonzero <shift> - Sets shift amount or PETSC_DECIDE for the default 262 - -pc_factor_shift_positive_definite [PETSC_TRUE/PETSC_FALSE] - Activate/Deactivate PCFactorSetShiftPd(); the value 263 is optional with PETSC_TRUE being the default 264 265 Level: beginner 266 267 Concepts: incomplete Cholesky factorization 268 269 Notes: Only implemented for some matrix formats. Not implemented in parallel 270 271 For BAIJ matrices this implements a point block ICC. 272 273 The Manteuffel shift is only implemented for matrices with block size 1 274 275 By default, the Manteuffel is applied (for matrices with block size 1). Call PCFactorSetShiftPd(pc,PETSC_FALSE); 276 to turn off the shift. 277 278 279 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, PCSOR, MatOrderingType, 280 PCFactorSetZeroPivot(), PCFactorSetShiftNonzero(), PCFactorSetShiftPd(), 281 PCFactorSetFill(), PCFactorSetMatOrdering(), PCFactorSetReuseOrdering(), 282 PCFactorSetLevels(),PCFactorSetShiftNonzero(),PCFactorSetShiftPd(), 283 284 M*/ 285 286 EXTERN_C_BEGIN 287 #undef __FUNCT__ 288 #define __FUNCT__ "PCCreate_ICC" 289 PetscErrorCode PETSCKSP_DLLEXPORT PCCreate_ICC(PC pc) 290 { 291 PetscErrorCode ierr; 292 PC_ICC *icc; 293 294 PetscFunctionBegin; 295 ierr = PetscNew(PC_ICC,&icc);CHKERRQ(ierr); 296 ierr = PetscLogObjectMemory(pc,sizeof(PC_ICC));CHKERRQ(ierr); 297 298 icc->fact = 0; 299 ierr = PetscStrallocpy(MATORDERING_NATURAL,&icc->ordering);CHKERRQ(ierr); 300 ierr = MatFactorInfoInitialize(&icc->info);CHKERRQ(ierr); 301 icc->info.levels = 0; 302 icc->info.fill = 1.0; 303 icc->implctx = 0; 304 305 icc->info.dtcol = PETSC_DEFAULT; 306 icc->info.shiftnz = 0.0; 307 icc->info.shiftpd = 1.0; /* true */ 308 icc->info.shift_fraction = 0.0; 309 icc->info.zeropivot = 1.e-12; 310 pc->data = (void*)icc; 311 312 pc->ops->apply = PCApply_ICC; 313 pc->ops->setup = PCSetup_ICC; 314 pc->ops->destroy = PCDestroy_ICC; 315 pc->ops->setfromoptions = PCSetFromOptions_ICC; 316 pc->ops->view = PCView_ICC; 317 pc->ops->getfactoredmatrix = PCGetFactoredMatrix_ICC; 318 pc->ops->applysymmetricleft = PCApplySymmetricLeft_ICC; 319 pc->ops->applysymmetricright = PCApplySymmetricRight_ICC; 320 321 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetZeroPivot_C","PCFactorSetZeroPivot_ICC", 322 PCFactorSetZeroPivot_ICC);CHKERRQ(ierr); 323 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftNonzero_C","PCFactorSetShiftNonzero_ICC", 324 PCFactorSetShiftNonzero_ICC);CHKERRQ(ierr); 325 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftPd_C","PCFactorSetShiftPd_ICC", 326 PCFactorSetShiftPd_ICC);CHKERRQ(ierr); 327 328 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetLevels_C","PCFactorSetLevels_ICC", 329 PCFactorSetLevels_ICC);CHKERRQ(ierr); 330 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetFill_C","PCFactorSetFill_ICC", 331 PCFactorSetFill_ICC);CHKERRQ(ierr); 332 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetMatOrdering_C","PCFactorSetMatOrdering_ICC", 333 PCFactorSetMatOrdering_ICC);CHKERRQ(ierr); 334 PetscFunctionReturn(0); 335 } 336 EXTERN_C_END 337 338 339