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