1 2 /* 3 Defines a ILU factorization preconditioner for any Mat implementation 4 */ 5 #include <../src/ksp/pc/impls/factor/ilu/ilu.h> /*I "petscpc.h" I*/ 6 7 /* ------------------------------------------------------------------------------------------*/ 8 #undef __FUNCT__ 9 #define __FUNCT__ "PCFactorSetReuseFill_ILU" 10 PetscErrorCode PCFactorSetReuseFill_ILU(PC pc,PetscBool flag) 11 { 12 PC_ILU *lu = (PC_ILU*)pc->data; 13 14 PetscFunctionBegin; 15 lu->reusefill = flag; 16 PetscFunctionReturn(0); 17 } 18 19 #undef __FUNCT__ 20 #define __FUNCT__ "PCFactorReorderForNonzeroDiagonal_ILU" 21 PetscErrorCode PCFactorReorderForNonzeroDiagonal_ILU(PC pc,PetscReal z) 22 { 23 PC_ILU *ilu = (PC_ILU*)pc->data; 24 25 PetscFunctionBegin; 26 ilu->nonzerosalongdiagonal = PETSC_TRUE; 27 if (z == PETSC_DECIDE) ilu->nonzerosalongdiagonaltol = 1.e-10; 28 else ilu->nonzerosalongdiagonaltol = z; 29 PetscFunctionReturn(0); 30 } 31 32 #undef __FUNCT__ 33 #define __FUNCT__ "PCReset_ILU" 34 PetscErrorCode PCReset_ILU(PC pc) 35 { 36 PC_ILU *ilu = (PC_ILU*)pc->data; 37 PetscErrorCode ierr; 38 39 PetscFunctionBegin; 40 if (!ilu->inplace) {ierr = MatDestroy(&((PC_Factor*)ilu)->fact);CHKERRQ(ierr);} 41 if (ilu->row && ilu->col && ilu->row != ilu->col) {ierr = ISDestroy(&ilu->row);CHKERRQ(ierr);} 42 ierr = ISDestroy(&ilu->col);CHKERRQ(ierr); 43 PetscFunctionReturn(0); 44 } 45 46 #undef __FUNCT__ 47 #define __FUNCT__ "PCFactorSetDropTolerance_ILU" 48 PetscErrorCode PCFactorSetDropTolerance_ILU(PC pc,PetscReal dt,PetscReal dtcol,PetscInt dtcount) 49 { 50 PC_ILU *ilu = (PC_ILU*)pc->data; 51 52 PetscFunctionBegin; 53 if (pc->setupcalled && (((PC_Factor*)ilu)->info.dt != dt || ((PC_Factor*)ilu)->info.dtcol != dtcol || ((PC_Factor*)ilu)->info.dtcount != dtcount)) { 54 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Cannot change drop tolerance after using PC"); 55 } 56 ((PC_Factor*)ilu)->info.dt = dt; 57 ((PC_Factor*)ilu)->info.dtcol = dtcol; 58 ((PC_Factor*)ilu)->info.dtcount = dtcount; 59 ((PC_Factor*)ilu)->info.usedt = 1.0; 60 PetscFunctionReturn(0); 61 } 62 63 #undef __FUNCT__ 64 #define __FUNCT__ "PCFactorSetReuseOrdering_ILU" 65 PetscErrorCode PCFactorSetReuseOrdering_ILU(PC pc,PetscBool flag) 66 { 67 PC_ILU *ilu = (PC_ILU*)pc->data; 68 69 PetscFunctionBegin; 70 ilu->reuseordering = flag; 71 PetscFunctionReturn(0); 72 } 73 74 #undef __FUNCT__ 75 #define __FUNCT__ "PCFactorSetUseInPlace_ILU" 76 PetscErrorCode PCFactorSetUseInPlace_ILU(PC pc,PetscBool flg) 77 { 78 PC_ILU *dir = (PC_ILU*)pc->data; 79 80 PetscFunctionBegin; 81 dir->inplace = flg; 82 PetscFunctionReturn(0); 83 } 84 85 #undef __FUNCT__ 86 #define __FUNCT__ "PCFactorGetUseInPlace_ILU" 87 PetscErrorCode PCFactorGetUseInPlace_ILU(PC pc,PetscBool *flg) 88 { 89 PC_ILU *dir = (PC_ILU*)pc->data; 90 91 PetscFunctionBegin; 92 *flg = dir->inplace; 93 PetscFunctionReturn(0); 94 } 95 96 #undef __FUNCT__ 97 #define __FUNCT__ "PCSetFromOptions_ILU" 98 static PetscErrorCode PCSetFromOptions_ILU(PetscOptionItems *PetscOptionsObject,PC pc) 99 { 100 PetscErrorCode ierr; 101 PetscInt itmp; 102 PetscBool flg,set; 103 PC_ILU *ilu = (PC_ILU*)pc->data; 104 PetscReal tol; 105 106 PetscFunctionBegin; 107 ierr = PetscOptionsHead(PetscOptionsObject,"ILU Options");CHKERRQ(ierr); 108 ierr = PCSetFromOptions_Factor(PetscOptionsObject,pc);CHKERRQ(ierr); 109 110 ierr = PetscOptionsInt("-pc_factor_levels","levels of fill","PCFactorSetLevels",(PetscInt)((PC_Factor*)ilu)->info.levels,&itmp,&flg);CHKERRQ(ierr); 111 if (flg) ((PC_Factor*)ilu)->info.levels = itmp; 112 113 ierr = PetscOptionsBool("-pc_factor_diagonal_fill","Allow fill into empty diagonal entry","PCFactorSetAllowDiagonalFill",((PC_Factor*)ilu)->info.diagonal_fill ? PETSC_TRUE : PETSC_FALSE,&flg,&set);CHKERRQ(ierr); 114 if (set) ((PC_Factor*)ilu)->info.diagonal_fill = (PetscReal) flg; 115 ierr = PetscOptionsName("-pc_factor_nonzeros_along_diagonal","Reorder to remove zeros from diagonal","PCFactorReorderForNonzeroDiagonal",&flg);CHKERRQ(ierr); 116 if (flg) { 117 tol = PETSC_DECIDE; 118 ierr = PetscOptionsReal("-pc_factor_nonzeros_along_diagonal","Reorder to remove zeros from diagonal","PCFactorReorderForNonzeroDiagonal",ilu->nonzerosalongdiagonaltol,&tol,0);CHKERRQ(ierr); 119 ierr = PCFactorReorderForNonzeroDiagonal(pc,tol);CHKERRQ(ierr); 120 } 121 122 ierr = PetscOptionsTail();CHKERRQ(ierr); 123 PetscFunctionReturn(0); 124 } 125 126 #undef __FUNCT__ 127 #define __FUNCT__ "PCView_ILU" 128 static PetscErrorCode PCView_ILU(PC pc,PetscViewer viewer) 129 { 130 PC_ILU *ilu = (PC_ILU*)pc->data; 131 PetscErrorCode ierr; 132 PetscBool iascii; 133 134 PetscFunctionBegin; 135 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 136 if (iascii) { 137 if (ilu->inplace) { 138 ierr = PetscViewerASCIIPrintf(viewer," ILU: in-place factorization\n");CHKERRQ(ierr); 139 } else { 140 ierr = PetscViewerASCIIPrintf(viewer," ILU: out-of-place factorization\n");CHKERRQ(ierr); 141 } 142 143 if (ilu->reusefill) {ierr = PetscViewerASCIIPrintf(viewer," ILU: Reusing fill from past factorization\n");CHKERRQ(ierr);} 144 if (ilu->reuseordering) {ierr = PetscViewerASCIIPrintf(viewer," ILU: Reusing reordering from past factorization\n");CHKERRQ(ierr);} 145 } 146 ierr = PCView_Factor(pc,viewer);CHKERRQ(ierr); 147 PetscFunctionReturn(0); 148 } 149 150 #undef __FUNCT__ 151 #define __FUNCT__ "PCSetUp_ILU" 152 static PetscErrorCode PCSetUp_ILU(PC pc) 153 { 154 PetscErrorCode ierr; 155 PC_ILU *ilu = (PC_ILU*)pc->data; 156 MatInfo info; 157 PetscBool flg; 158 const MatSolverPackage stype; 159 160 PetscFunctionBegin; 161 /* ugly hack to change default, since it is not support by some matrix types */ 162 if (((PC_Factor*)ilu)->info.shifttype == (PetscReal)MAT_SHIFT_NONZERO) { 163 ierr = PetscObjectTypeCompare((PetscObject)pc->pmat,MATSEQAIJ,&flg);CHKERRQ(ierr); 164 if (!flg) { 165 ierr = PetscObjectTypeCompare((PetscObject)pc->pmat,MATMPIAIJ,&flg);CHKERRQ(ierr); 166 if (!flg) { 167 ((PC_Factor*)ilu)->info.shifttype = (PetscReal)MAT_SHIFT_INBLOCKS; 168 PetscInfo(pc,"Changing shift type from NONZERO to INBLOCKS because block matrices do not support NONZERO\n");CHKERRQ(ierr); 169 } 170 } 171 } 172 173 ierr = MatSetErrorIfFailure(pc->pmat,pc->erroriffailure);CHKERRQ(ierr); 174 if (ilu->inplace) { 175 if (!pc->setupcalled) { 176 177 /* In-place factorization only makes sense with the natural ordering, 178 so we only need to get the ordering once, even if nonzero structure changes */ 179 ierr = MatGetOrdering(pc->pmat,((PC_Factor*)ilu)->ordering,&ilu->row,&ilu->col);CHKERRQ(ierr); 180 if (ilu->row) {ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)ilu->row);CHKERRQ(ierr);} 181 if (ilu->col) {ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)ilu->col);CHKERRQ(ierr);} 182 } 183 184 /* In place ILU only makes sense with fill factor of 1.0 because 185 cannot have levels of fill */ 186 ((PC_Factor*)ilu)->info.fill = 1.0; 187 ((PC_Factor*)ilu)->info.diagonal_fill = 0.0; 188 189 ierr = MatILUFactor(pc->pmat,ilu->row,ilu->col,&((PC_Factor*)ilu)->info);CHKERRQ(ierr);CHKERRQ(ierr); 190 if (pc->pmat->errortype) { /* Factor() fails */ 191 pc->failedreason = (PCFailedReason)pc->pmat->errortype; 192 PetscFunctionReturn(0); 193 } 194 195 ((PC_Factor*)ilu)->fact = pc->pmat; 196 /* must update the pc record of the matrix state or the PC will attempt to run PCSetUp() yet again */ 197 ierr = PetscObjectStateGet((PetscObject)pc->pmat,&pc->matstate);CHKERRQ(ierr); 198 } else { 199 Mat F; 200 if (!pc->setupcalled) { 201 /* first time in so compute reordering and symbolic factorization */ 202 ierr = MatGetOrdering(pc->pmat,((PC_Factor*)ilu)->ordering,&ilu->row,&ilu->col);CHKERRQ(ierr); 203 if (ilu->row) {ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)ilu->row);CHKERRQ(ierr);} 204 if (ilu->col) {ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)ilu->col);CHKERRQ(ierr);} 205 /* Remove zeros along diagonal? */ 206 if (ilu->nonzerosalongdiagonal) { 207 ierr = MatReorderForNonzeroDiagonal(pc->pmat,ilu->nonzerosalongdiagonaltol,ilu->row,ilu->col);CHKERRQ(ierr); 208 } 209 if (!((PC_Factor*)ilu)->fact) { 210 ierr = MatGetFactor(pc->pmat,((PC_Factor*)ilu)->solvertype,MAT_FACTOR_ILU,&((PC_Factor*)ilu)->fact);CHKERRQ(ierr); 211 } 212 ierr = MatILUFactorSymbolic(((PC_Factor*)ilu)->fact,pc->pmat,ilu->row,ilu->col,&((PC_Factor*)ilu)->info);CHKERRQ(ierr); 213 ierr = MatGetInfo(((PC_Factor*)ilu)->fact,MAT_LOCAL,&info);CHKERRQ(ierr); 214 215 ilu->actualfill = info.fill_ratio_needed; 216 217 ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)((PC_Factor*)ilu)->fact);CHKERRQ(ierr); 218 } else if (pc->flag != SAME_NONZERO_PATTERN) { 219 if (!ilu->reuseordering) { 220 /* compute a new ordering for the ILU */ 221 ierr = ISDestroy(&ilu->row);CHKERRQ(ierr); 222 ierr = ISDestroy(&ilu->col);CHKERRQ(ierr); 223 ierr = MatGetOrdering(pc->pmat,((PC_Factor*)ilu)->ordering,&ilu->row,&ilu->col);CHKERRQ(ierr); 224 if (ilu->row) {ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)ilu->row);CHKERRQ(ierr);} 225 if (ilu->col) {ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)ilu->col);CHKERRQ(ierr);} 226 /* Remove zeros along diagonal? */ 227 if (ilu->nonzerosalongdiagonal) { 228 ierr = MatReorderForNonzeroDiagonal(pc->pmat,ilu->nonzerosalongdiagonaltol,ilu->row,ilu->col);CHKERRQ(ierr); 229 } 230 } 231 ierr = MatDestroy(&((PC_Factor*)ilu)->fact);CHKERRQ(ierr); 232 ierr = MatGetFactor(pc->pmat,((PC_Factor*)ilu)->solvertype,MAT_FACTOR_ILU,&((PC_Factor*)ilu)->fact);CHKERRQ(ierr); 233 ierr = MatILUFactorSymbolic(((PC_Factor*)ilu)->fact,pc->pmat,ilu->row,ilu->col,&((PC_Factor*)ilu)->info);CHKERRQ(ierr); 234 ierr = MatGetInfo(((PC_Factor*)ilu)->fact,MAT_LOCAL,&info);CHKERRQ(ierr); 235 236 ilu->actualfill = info.fill_ratio_needed; 237 238 ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)((PC_Factor*)ilu)->fact);CHKERRQ(ierr); 239 } 240 F = ((PC_Factor*)ilu)->fact; 241 if (F->errortype) { /* FactorSymbolic() fails */ 242 pc->failedreason = (PCFailedReason)F->errortype; 243 PetscFunctionReturn(0); 244 } 245 246 ierr = MatLUFactorNumeric(((PC_Factor*)ilu)->fact,pc->pmat,&((PC_Factor*)ilu)->info);CHKERRQ(ierr); 247 if (F->errortype) { /* FactorNumeric() fails */ 248 pc->failedreason = (PCFailedReason)F->errortype; 249 } 250 } 251 252 ierr = PCFactorGetMatSolverPackage(pc,&stype);CHKERRQ(ierr); 253 if (!stype) { 254 ierr = PCFactorSetMatSolverPackage(pc,((PC_Factor*)ilu)->fact->solvertype);CHKERRQ(ierr); 255 } 256 PetscFunctionReturn(0); 257 } 258 259 #undef __FUNCT__ 260 #define __FUNCT__ "PCDestroy_ILU" 261 static PetscErrorCode PCDestroy_ILU(PC pc) 262 { 263 PC_ILU *ilu = (PC_ILU*)pc->data; 264 PetscErrorCode ierr; 265 266 PetscFunctionBegin; 267 ierr = PCReset_ILU(pc);CHKERRQ(ierr); 268 ierr = PetscFree(((PC_Factor*)ilu)->solvertype);CHKERRQ(ierr); 269 ierr = PetscFree(((PC_Factor*)ilu)->ordering);CHKERRQ(ierr); 270 ierr = PetscFree(pc->data);CHKERRQ(ierr); 271 PetscFunctionReturn(0); 272 } 273 274 #undef __FUNCT__ 275 #define __FUNCT__ "PCApply_ILU" 276 static PetscErrorCode PCApply_ILU(PC pc,Vec x,Vec y) 277 { 278 PC_ILU *ilu = (PC_ILU*)pc->data; 279 PetscErrorCode ierr; 280 281 PetscFunctionBegin; 282 ierr = MatSolve(((PC_Factor*)ilu)->fact,x,y);CHKERRQ(ierr); 283 PetscFunctionReturn(0); 284 } 285 286 #undef __FUNCT__ 287 #define __FUNCT__ "PCApplyTranspose_ILU" 288 static PetscErrorCode PCApplyTranspose_ILU(PC pc,Vec x,Vec y) 289 { 290 PC_ILU *ilu = (PC_ILU*)pc->data; 291 PetscErrorCode ierr; 292 293 PetscFunctionBegin; 294 ierr = MatSolveTranspose(((PC_Factor*)ilu)->fact,x,y);CHKERRQ(ierr); 295 PetscFunctionReturn(0); 296 } 297 298 #undef __FUNCT__ 299 #define __FUNCT__ "PCApplySymmetricLeft_ILU" 300 static PetscErrorCode PCApplySymmetricLeft_ILU(PC pc,Vec x,Vec y) 301 { 302 PetscErrorCode ierr; 303 PC_ILU *icc = (PC_ILU*)pc->data; 304 305 PetscFunctionBegin; 306 ierr = MatForwardSolve(((PC_Factor*)icc)->fact,x,y);CHKERRQ(ierr); 307 PetscFunctionReturn(0); 308 } 309 310 #undef __FUNCT__ 311 #define __FUNCT__ "PCApplySymmetricRight_ILU" 312 static PetscErrorCode PCApplySymmetricRight_ILU(PC pc,Vec x,Vec y) 313 { 314 PetscErrorCode ierr; 315 PC_ILU *icc = (PC_ILU*)pc->data; 316 317 PetscFunctionBegin; 318 ierr = MatBackwardSolve(((PC_Factor*)icc)->fact,x,y);CHKERRQ(ierr); 319 PetscFunctionReturn(0); 320 } 321 322 /*MC 323 PCILU - Incomplete factorization preconditioners. 324 325 Options Database Keys: 326 + -pc_factor_levels <k> - number of levels of fill for ILU(k) 327 . -pc_factor_in_place - only for ILU(0) with natural ordering, reuses the space of the matrix for 328 its factorization (overwrites original matrix) 329 . -pc_factor_diagonal_fill - fill in a zero diagonal even if levels of fill indicate it wouldn't be fill 330 . -pc_factor_reuse_ordering - reuse ordering of factorized matrix from previous factorization 331 . -pc_factor_fill <nfill> - expected amount of fill in factored matrix compared to original matrix, nfill > 1 332 . -pc_factor_nonzeros_along_diagonal - reorder the matrix before factorization to remove zeros from the diagonal, 333 this decreases the chance of getting a zero pivot 334 . -pc_factor_mat_ordering_type <natural,nd,1wd,rcm,qmd> - set the row/column ordering of the factored matrix 335 - -pc_factor_pivot_in_blocks - for block ILU(k) factorization, i.e. with BAIJ matrices with block size larger 336 than 1 the diagonal blocks are factored with partial pivoting (this increases the 337 stability of the ILU factorization 338 339 Level: beginner 340 341 Concepts: incomplete factorization 342 343 Notes: Only implemented for some matrix formats. (for parallel see PCHYPRE for hypre's ILU) 344 345 For BAIJ matrices this implements a point block ILU 346 347 The "symmetric" application of this preconditioner is not actually symmetric since L is not transpose(U) 348 even when the matrix is not symmetric since the U stores the diagonals of the factorization. 349 350 If you are using MATSEQAIJCUSPARSE matrices (or MATMPIAIJCUSPARESE matrices with block Jacobi), factorization 351 is never done on the GPU). 352 353 References: 354 + 1. - T. Dupont, R. Kendall, and H. Rachford. An approximate factorization procedure for solving 355 self adjoint elliptic difference equations. SIAM J. Numer. Anal., 5, 1968. 356 . 2. - T.A. Oliphant. An implicit numerical method for solving two dimensional timedependent diffusion problems. Quart. Appl. Math., 19, 1961. 357 - 3. - TONY F. CHAN AND HENK A. VAN DER VORST, APPROXIMATE AND INCOMPLETE FACTORIZATIONS, 358 Chapter in Parallel Numerical 359 Algorithms, edited by D. Keyes, A. Semah, V. Venkatakrishnan, ICASE/LaRC Interdisciplinary Series in 360 Science and Engineering, Kluwer. 361 362 363 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, PCSOR, MatOrderingType, 364 PCFactorSetZeroPivot(), PCFactorSetShiftSetType(), PCFactorSetAmount(), 365 PCFactorSetDropTolerance(),PCFactorSetFill(), PCFactorSetMatOrderingType(), PCFactorSetReuseOrdering(), 366 PCFactorSetLevels(), PCFactorSetUseInPlace(), PCFactorSetAllowDiagonalFill(), PCFactorSetPivotInBlocks(), 367 PCFactorGetAllowDiagonalFill(), PCFactorGetUseInPlace() 368 369 M*/ 370 371 #undef __FUNCT__ 372 #define __FUNCT__ "PCCreate_ILU" 373 PETSC_EXTERN PetscErrorCode PCCreate_ILU(PC pc) 374 { 375 PetscErrorCode ierr; 376 PC_ILU *ilu; 377 378 PetscFunctionBegin; 379 ierr = PetscNewLog(pc,&ilu);CHKERRQ(ierr); 380 381 ((PC_Factor*)ilu)->fact = 0; 382 ierr = MatFactorInfoInitialize(&((PC_Factor*)ilu)->info);CHKERRQ(ierr); 383 ((PC_Factor*)ilu)->factortype = MAT_FACTOR_ILU; 384 ((PC_Factor*)ilu)->info.levels = 0.; 385 ((PC_Factor*)ilu)->info.fill = 1.0; 386 ilu->col = 0; 387 ilu->row = 0; 388 ilu->inplace = PETSC_FALSE; 389 ierr = PetscStrallocpy(MATORDERINGNATURAL,(char**)&((PC_Factor*)ilu)->ordering);CHKERRQ(ierr); 390 ilu->reuseordering = PETSC_FALSE; 391 ((PC_Factor*)ilu)->info.dt = PETSC_DEFAULT; 392 ((PC_Factor*)ilu)->info.dtcount = PETSC_DEFAULT; 393 ((PC_Factor*)ilu)->info.dtcol = PETSC_DEFAULT; 394 ((PC_Factor*)ilu)->info.shifttype = (PetscReal)MAT_SHIFT_NONE; 395 ((PC_Factor*)ilu)->info.shiftamount = 100.0*PETSC_MACHINE_EPSILON; 396 ((PC_Factor*)ilu)->info.zeropivot = 100.0*PETSC_MACHINE_EPSILON; 397 ((PC_Factor*)ilu)->info.pivotinblocks = 1.0; 398 ilu->reusefill = PETSC_FALSE; 399 ((PC_Factor*)ilu)->info.diagonal_fill = 0.0; 400 pc->data = (void*)ilu; 401 402 pc->ops->reset = PCReset_ILU; 403 pc->ops->destroy = PCDestroy_ILU; 404 pc->ops->apply = PCApply_ILU; 405 pc->ops->applytranspose = PCApplyTranspose_ILU; 406 pc->ops->setup = PCSetUp_ILU; 407 pc->ops->setfromoptions = PCSetFromOptions_ILU; 408 pc->ops->getfactoredmatrix = PCFactorGetMatrix_Factor; 409 pc->ops->view = PCView_ILU; 410 pc->ops->applysymmetricleft = PCApplySymmetricLeft_ILU; 411 pc->ops->applysymmetricright = PCApplySymmetricRight_ILU; 412 pc->ops->applyrichardson = 0; 413 414 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetZeroPivot_C",PCFactorSetZeroPivot_Factor);CHKERRQ(ierr); 415 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetShiftType_C",PCFactorSetShiftType_Factor);CHKERRQ(ierr); 416 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetShiftAmount_C",PCFactorSetShiftAmount_Factor);CHKERRQ(ierr); 417 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetMatSolverPackage_C",PCFactorGetMatSolverPackage_Factor);CHKERRQ(ierr); 418 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetMatSolverPackage_C",PCFactorSetMatSolverPackage_Factor);CHKERRQ(ierr); 419 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetUpMatSolverPackage_C",PCFactorSetUpMatSolverPackage_Factor);CHKERRQ(ierr); 420 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetDropTolerance_C",PCFactorSetDropTolerance_ILU);CHKERRQ(ierr); 421 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetFill_C",PCFactorSetFill_Factor);CHKERRQ(ierr); 422 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetMatOrderingType_C",PCFactorSetMatOrderingType_Factor);CHKERRQ(ierr); 423 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetReuseOrdering_C",PCFactorSetReuseOrdering_ILU);CHKERRQ(ierr); 424 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetReuseFill_C",PCFactorSetReuseFill_ILU);CHKERRQ(ierr); 425 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetLevels_C",PCFactorSetLevels_Factor);CHKERRQ(ierr); 426 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetLevels_C",PCFactorGetLevels_Factor);CHKERRQ(ierr); 427 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetUseInPlace_C",PCFactorSetUseInPlace_ILU);CHKERRQ(ierr); 428 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetUseInPlace_C",PCFactorGetUseInPlace_ILU);CHKERRQ(ierr); 429 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetAllowDiagonalFill_C",PCFactorSetAllowDiagonalFill_Factor);CHKERRQ(ierr); 430 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorGetAllowDiagonalFill_C",PCFactorGetAllowDiagonalFill_Factor);CHKERRQ(ierr); 431 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorSetPivotInBlocks_C",PCFactorSetPivotInBlocks_Factor);CHKERRQ(ierr); 432 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCFactorReorderForNonzeroDiagonal_C",PCFactorReorderForNonzeroDiagonal_ILU);CHKERRQ(ierr); 433 PetscFunctionReturn(0); 434 } 435