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 9 #include "../src/ksp/pc/impls/factor/lu/lu.h" /*I "petscpc.h" I*/ 10 11 EXTERN_C_BEGIN 12 #undef __FUNCT__ 13 #define __FUNCT__ "PCFactorReorderForNonzeroDiagonal_LU" 14 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorReorderForNonzeroDiagonal_LU(PC pc,PetscReal z) 15 { 16 PC_LU *lu = (PC_LU*)pc->data; 17 18 PetscFunctionBegin; 19 lu->nonzerosalongdiagonal = PETSC_TRUE; 20 if (z == PETSC_DECIDE) { 21 lu->nonzerosalongdiagonaltol = 1.e-10; 22 } else { 23 lu->nonzerosalongdiagonaltol = z; 24 } 25 PetscFunctionReturn(0); 26 } 27 EXTERN_C_END 28 29 EXTERN_C_BEGIN 30 #undef __FUNCT__ 31 #define __FUNCT__ "PCFactorSetReuseOrdering_LU" 32 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseOrdering_LU(PC pc,PetscTruth flag) 33 { 34 PC_LU *lu = (PC_LU*)pc->data; 35 36 PetscFunctionBegin; 37 lu->reuseordering = flag; 38 PetscFunctionReturn(0); 39 } 40 EXTERN_C_END 41 42 EXTERN_C_BEGIN 43 #undef __FUNCT__ 44 #define __FUNCT__ "PCFactorSetReuseFill_LU" 45 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetReuseFill_LU(PC pc,PetscTruth flag) 46 { 47 PC_LU *lu = (PC_LU*)pc->data; 48 49 PetscFunctionBegin; 50 lu->reusefill = flag; 51 PetscFunctionReturn(0); 52 } 53 EXTERN_C_END 54 55 #undef __FUNCT__ 56 #define __FUNCT__ "PCSetFromOptions_LU" 57 static PetscErrorCode PCSetFromOptions_LU(PC pc) 58 { 59 PC_LU *lu = (PC_LU*)pc->data; 60 PetscErrorCode ierr; 61 PetscTruth flg = PETSC_FALSE,set; 62 char tname[256], solvertype[64]; 63 PetscFList ordlist; 64 PetscReal tol; 65 66 PetscFunctionBegin; 67 if (!MatOrderingRegisterAllCalled) {ierr = MatOrderingRegisterAll(PETSC_NULL);CHKERRQ(ierr);} 68 ierr = PetscOptionsHead("LU options");CHKERRQ(ierr); 69 ierr = PetscOptionsTruth("-pc_factor_in_place","Form LU in the same memory as the matrix","PCFactorSetUseInPlace",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 70 if (flg) { 71 ierr = PCFactorSetUseInPlace(pc);CHKERRQ(ierr); 72 } 73 ierr = PetscOptionsReal("-pc_factor_fill","Expected non-zeros in LU/non-zeros in matrix","PCFactorSetFill",((PC_Factor*)lu)->info.fill,&((PC_Factor*)lu)->info.fill,0);CHKERRQ(ierr); 74 75 ierr = PetscOptionsName("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",&flg);CHKERRQ(ierr); 76 if (flg) { 77 ierr = PCFactorSetShiftNonzero(pc,(PetscReal) PETSC_DECIDE);CHKERRQ(ierr); 78 } 79 ierr = PetscOptionsReal("-pc_factor_shift_nonzero","Shift added to diagonal","PCFactorSetShiftNonzero",((PC_Factor*)lu)->info.shiftnz,&((PC_Factor*)lu)->info.shiftnz,0);CHKERRQ(ierr); 80 flg = PETSC_FALSE; 81 ierr = PetscOptionsTruth("-pc_factor_shift_positive_definite","Manteuffel shift applied to diagonal","PCFactorSetShiftPd",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 82 if (flg) { 83 ierr = PCFactorSetShiftPd(pc,PETSC_TRUE);CHKERRQ(ierr); 84 } 85 ierr = PetscOptionsReal("-pc_factor_zeropivot","Pivot is considered zero if less than","PCFactorSetZeroPivot",((PC_Factor*)lu)->info.zeropivot,&((PC_Factor*)lu)->info.zeropivot,0);CHKERRQ(ierr); 86 87 flg = PETSC_FALSE; 88 ierr = PetscOptionsTruth("-pc_factor_reuse_fill","Use fill from previous factorization","PCFactorSetReuseFill",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 89 if (flg) { 90 ierr = PCFactorSetReuseFill(pc,PETSC_TRUE);CHKERRQ(ierr); 91 } 92 flg = PETSC_FALSE; 93 ierr = PetscOptionsTruth("-pc_factor_reuse_ordering","Reuse ordering from previous factorization","PCFactorSetReuseOrdering",flg,&flg,PETSC_NULL);CHKERRQ(ierr); 94 if (flg) { 95 ierr = PCFactorSetReuseOrdering(pc,PETSC_TRUE);CHKERRQ(ierr); 96 } 97 98 ierr = MatGetOrderingList(&ordlist);CHKERRQ(ierr); 99 ierr = PetscOptionsList("-pc_factor_mat_ordering_type","Reordering to reduce nonzeros in LU","PCFactorSetMatOrderingType",ordlist,((PC_Factor*)lu)->ordering,tname,256,&flg);CHKERRQ(ierr); 100 if (flg) { 101 ierr = PCFactorSetMatOrderingType(pc,tname);CHKERRQ(ierr); 102 } 103 104 /* maybe should have MatGetSolverTypes(Mat,&list) like the ordering list */ 105 ierr = PetscOptionsString("-pc_factor_mat_solver_package","Specific LU solver to use","MatGetFactor",((PC_Factor*)lu)->solvertype,solvertype,64,&flg);CHKERRQ(ierr); 106 if (flg) { 107 ierr = PCFactorSetMatSolverPackage(pc,solvertype);CHKERRQ(ierr); 108 } 109 110 ierr = PetscOptionsName("-pc_factor_nonzeros_along_diagonal","Reorder to remove zeros from diagonal","PCFactorReorderForNonzeroDiagonal",&flg);CHKERRQ(ierr); 111 if (flg) { 112 tol = PETSC_DECIDE; 113 ierr = PetscOptionsReal("-pc_factor_nonzeros_along_diagonal","Reorder to remove zeros from diagonal","PCFactorReorderForNonzeroDiagonal",lu->nonzerosalongdiagonaltol,&tol,0);CHKERRQ(ierr); 114 ierr = PCFactorReorderForNonzeroDiagonal(pc,tol);CHKERRQ(ierr); 115 } 116 117 ierr = PetscOptionsReal("-pc_factor_pivoting","Pivoting tolerance (used only for some factorization)","PCFactorSetPivoting",((PC_Factor*)lu)->info.dtcol,&((PC_Factor*)lu)->info.dtcol,&flg);CHKERRQ(ierr); 118 119 flg = ((PC_Factor*)lu)->info.pivotinblocks ? PETSC_TRUE : PETSC_FALSE; 120 ierr = PetscOptionsTruth("-pc_factor_pivot_in_blocks","Pivot inside matrix blocks for BAIJ and SBAIJ","PCFactorSetPivotInBlocks",flg,&flg,&set);CHKERRQ(ierr); 121 if (set) { 122 ierr = PCFactorSetPivotInBlocks(pc,flg);CHKERRQ(ierr); 123 } 124 ierr = PetscOptionsTail();CHKERRQ(ierr); 125 PetscFunctionReturn(0); 126 } 127 128 #undef __FUNCT__ 129 #define __FUNCT__ "PCView_LU" 130 static PetscErrorCode PCView_LU(PC pc,PetscViewer viewer) 131 { 132 PC_LU *lu = (PC_LU*)pc->data; 133 PetscErrorCode ierr; 134 PetscTruth iascii,isstring; 135 136 PetscFunctionBegin; 137 ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr); 138 ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);CHKERRQ(ierr); 139 if (iascii) { 140 141 if (lu->inplace) {ierr = PetscViewerASCIIPrintf(viewer," LU: in-place factorization\n");CHKERRQ(ierr);} 142 else {ierr = PetscViewerASCIIPrintf(viewer," LU: out-of-place factorization\n");CHKERRQ(ierr);} 143 ierr = PetscViewerASCIIPrintf(viewer," matrix ordering: %s\n",((PC_Factor*)lu)->ordering);CHKERRQ(ierr); 144 ierr = PetscViewerASCIIPrintf(viewer," LU: tolerance for zero pivot %G\n",((PC_Factor*)lu)->info.zeropivot);CHKERRQ(ierr); 145 if (((PC_Factor*)lu)->info.shiftpd) {ierr = PetscViewerASCIIPrintf(viewer," LU: using Manteuffel shift\n");CHKERRQ(ierr);} 146 if (((PC_Factor*)lu)->fact) { 147 ierr = PetscViewerASCIIPrintf(viewer," LU: factor fill ratio needed %G\n",lu->actualfill);CHKERRQ(ierr); 148 ierr = PetscViewerASCIIPrintf(viewer," Factored matrix follows\n");CHKERRQ(ierr); 149 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 150 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 151 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 152 ierr = PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);CHKERRQ(ierr); 153 ierr = MatView(((PC_Factor*)lu)->fact,viewer);CHKERRQ(ierr); 154 ierr = PetscViewerPopFormat(viewer);CHKERRQ(ierr); 155 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 156 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 157 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 158 } 159 if (lu->reusefill) {ierr = PetscViewerASCIIPrintf(viewer," Reusing fill from past factorization\n");CHKERRQ(ierr);} 160 if (lu->reuseordering) {ierr = PetscViewerASCIIPrintf(viewer," Reusing reordering from past factorization\n");CHKERRQ(ierr);} 161 } else if (isstring) { 162 ierr = PetscViewerStringSPrintf(viewer," order=%s",((PC_Factor*)lu)->ordering);CHKERRQ(ierr);CHKERRQ(ierr); 163 } else { 164 SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PCLU",((PetscObject)viewer)->type_name); 165 } 166 PetscFunctionReturn(0); 167 } 168 169 #undef __FUNCT__ 170 #define __FUNCT__ "PCSetUp_LU" 171 static PetscErrorCode PCSetUp_LU(PC pc) 172 { 173 PetscErrorCode ierr; 174 PC_LU *dir = (PC_LU*)pc->data; 175 176 PetscFunctionBegin; 177 if (dir->reusefill && pc->setupcalled) ((PC_Factor*)dir)->info.fill = dir->actualfill; 178 179 if (dir->inplace) { 180 if (dir->row && dir->col && dir->row != dir->col) {ierr = ISDestroy(dir->row);CHKERRQ(ierr);} 181 if (dir->col) {ierr = ISDestroy(dir->col);CHKERRQ(ierr);} 182 ierr = MatGetOrdering(pc->pmat,((PC_Factor*)dir)->ordering,&dir->row,&dir->col);CHKERRQ(ierr); 183 if (dir->row) { 184 ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr); 185 ierr = PetscLogObjectParent(pc,dir->col);CHKERRQ(ierr); 186 } 187 ierr = MatLUFactor(pc->pmat,dir->row,dir->col,&((PC_Factor*)dir)->info);CHKERRQ(ierr); 188 ((PC_Factor*)dir)->fact = pc->pmat; 189 } else { 190 MatInfo info; 191 if (!pc->setupcalled) { 192 ierr = MatGetOrdering(pc->pmat,((PC_Factor*)dir)->ordering,&dir->row,&dir->col);CHKERRQ(ierr); 193 if (dir->nonzerosalongdiagonal) { 194 ierr = MatReorderForNonzeroDiagonal(pc->pmat,dir->nonzerosalongdiagonaltol,dir->row,dir->col);CHKERRQ(ierr); 195 } 196 if (dir->row) { 197 ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr); 198 ierr = PetscLogObjectParent(pc,dir->col);CHKERRQ(ierr); 199 } 200 ierr = MatGetFactor(pc->pmat,((PC_Factor*)dir)->solvertype,MAT_FACTOR_LU,&((PC_Factor*)dir)->fact);CHKERRQ(ierr); 201 ierr = MatLUFactorSymbolic(((PC_Factor*)dir)->fact,pc->pmat,dir->row,dir->col,&((PC_Factor*)dir)->info);CHKERRQ(ierr); 202 ierr = MatGetInfo(((PC_Factor*)dir)->fact,MAT_LOCAL,&info);CHKERRQ(ierr); 203 dir->actualfill = info.fill_ratio_needed; 204 ierr = PetscLogObjectParent(pc,((PC_Factor*)dir)->fact);CHKERRQ(ierr); 205 } else if (pc->flag != SAME_NONZERO_PATTERN) { 206 if (!dir->reuseordering) { 207 if (dir->row && dir->col && dir->row != dir->col) {ierr = ISDestroy(dir->row);CHKERRQ(ierr);} 208 if (dir->col) {ierr = ISDestroy(dir->col);CHKERRQ(ierr);} 209 ierr = MatGetOrdering(pc->pmat,((PC_Factor*)dir)->ordering,&dir->row,&dir->col);CHKERRQ(ierr); 210 if (dir->nonzerosalongdiagonal) { 211 ierr = MatReorderForNonzeroDiagonal(pc->pmat,dir->nonzerosalongdiagonaltol,dir->row,dir->col);CHKERRQ(ierr); 212 } 213 if (dir->row) { 214 ierr = PetscLogObjectParent(pc,dir->row);CHKERRQ(ierr); 215 ierr = PetscLogObjectParent(pc,dir->col);CHKERRQ(ierr); 216 } 217 } 218 ierr = MatDestroy(((PC_Factor*)dir)->fact);CHKERRQ(ierr); 219 ierr = MatGetFactor(pc->pmat,((PC_Factor*)dir)->solvertype,MAT_FACTOR_LU,&((PC_Factor*)dir)->fact);CHKERRQ(ierr); 220 ierr = MatLUFactorSymbolic(((PC_Factor*)dir)->fact,pc->pmat,dir->row,dir->col,&((PC_Factor*)dir)->info);CHKERRQ(ierr); 221 ierr = MatGetInfo(((PC_Factor*)dir)->fact,MAT_LOCAL,&info);CHKERRQ(ierr); 222 dir->actualfill = info.fill_ratio_needed; 223 ierr = PetscLogObjectParent(pc,((PC_Factor*)dir)->fact);CHKERRQ(ierr); 224 } 225 ierr = MatLUFactorNumeric(((PC_Factor*)dir)->fact,pc->pmat,&((PC_Factor*)dir)->info);CHKERRQ(ierr); 226 } 227 PetscFunctionReturn(0); 228 } 229 230 #undef __FUNCT__ 231 #define __FUNCT__ "PCDestroy_LU" 232 static PetscErrorCode PCDestroy_LU(PC pc) 233 { 234 PC_LU *dir = (PC_LU*)pc->data; 235 PetscErrorCode ierr; 236 237 PetscFunctionBegin; 238 if (!dir->inplace && ((PC_Factor*)dir)->fact) {ierr = MatDestroy(((PC_Factor*)dir)->fact);CHKERRQ(ierr);} 239 if (dir->row && dir->col && dir->row != dir->col) {ierr = ISDestroy(dir->row);CHKERRQ(ierr);} 240 if (dir->col) {ierr = ISDestroy(dir->col);CHKERRQ(ierr);} 241 ierr = PetscStrfree(((PC_Factor*)dir)->ordering);CHKERRQ(ierr); 242 ierr = PetscStrfree(((PC_Factor*)dir)->solvertype);CHKERRQ(ierr); 243 ierr = PetscFree(dir);CHKERRQ(ierr); 244 PetscFunctionReturn(0); 245 } 246 247 #undef __FUNCT__ 248 #define __FUNCT__ "PCApply_LU" 249 static PetscErrorCode PCApply_LU(PC pc,Vec x,Vec y) 250 { 251 PC_LU *dir = (PC_LU*)pc->data; 252 PetscErrorCode ierr; 253 254 PetscFunctionBegin; 255 if (dir->inplace) {ierr = MatSolve(pc->pmat,x,y);CHKERRQ(ierr);} 256 else {ierr = MatSolve(((PC_Factor*)dir)->fact,x,y);CHKERRQ(ierr);} 257 PetscFunctionReturn(0); 258 } 259 260 #undef __FUNCT__ 261 #define __FUNCT__ "PCApplyTranspose_LU" 262 static PetscErrorCode PCApplyTranspose_LU(PC pc,Vec x,Vec y) 263 { 264 PC_LU *dir = (PC_LU*)pc->data; 265 PetscErrorCode ierr; 266 267 PetscFunctionBegin; 268 if (dir->inplace) {ierr = MatSolveTranspose(pc->pmat,x,y);CHKERRQ(ierr);} 269 else {ierr = MatSolveTranspose(((PC_Factor*)dir)->fact,x,y);CHKERRQ(ierr);} 270 PetscFunctionReturn(0); 271 } 272 273 /* -----------------------------------------------------------------------------------*/ 274 275 EXTERN_C_BEGIN 276 #undef __FUNCT__ 277 #define __FUNCT__ "PCFactorSetUseInPlace_LU" 278 PetscErrorCode PETSCKSP_DLLEXPORT PCFactorSetUseInPlace_LU(PC pc) 279 { 280 PC_LU *dir = (PC_LU*)pc->data; 281 282 PetscFunctionBegin; 283 dir->inplace = PETSC_TRUE; 284 PetscFunctionReturn(0); 285 } 286 EXTERN_C_END 287 288 /* ------------------------------------------------------------------------ */ 289 290 /*MC 291 PCLU - Uses a direct solver, based on LU factorization, as a preconditioner 292 293 Options Database Keys: 294 + -pc_factor_reuse_ordering - Activate PCFactorSetReuseOrdering() 295 . -pc_factor_mat_solver_package - Actives PCFactorSetMatSolverPackage() to choose the direct solver, like spooles 296 . -pc_factor_reuse_fill - Activates PCFactorSetReuseFill() 297 . -pc_factor_fill <fill> - Sets fill amount 298 . -pc_factor_in_place - Activates in-place factorization 299 . -pc_factor_mat_ordering_type <nd,rcm,...> - Sets ordering routine 300 . -pc_factor_pivot_in_blocks <true,false> - allow pivoting within the small blocks during factorization (may increase 301 stability of factorization. 302 . -pc_factor_shift_nonzero <shift> - Sets shift amount or PETSC_DECIDE for the default 303 . -pc_factor_shift_positive_definite [PETSC_TRUE/PETSC_FALSE] - Activate/Deactivate PCFactorSetShiftPd(); the value 304 is optional with PETSC_TRUE being the default 305 - -pc_factor_nonzeros_along_diagonal - permutes the rows and columns to try to put nonzero value along the 306 diagonal. 307 308 Notes: Not all options work for all matrix formats 309 Run with -help to see additional options for particular matrix formats or factorization 310 algorithms 311 312 Level: beginner 313 314 Concepts: LU factorization, direct solver 315 316 Notes: Usually this will compute an "exact" solution in one iteration and does 317 not need a Krylov method (i.e. you can use -ksp_type preonly, or 318 KSPSetType(ksp,KSPPREONLY) for the Krylov method 319 320 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 321 PCILU, PCCHOLESKY, PCICC, PCFactorSetReuseOrdering(), PCFactorSetReuseFill(), PCFactorGetMatrix(), 322 PCFactorSetFill(), PCFactorSetUseInPlace(), PCFactorSetMatOrderingType(), PCFactorSetPivoting(), 323 PCFactorSetPivotingInBlocks(),PCFactorSetShiftNonzero(),PCFactorSetShiftPd(), PCFactorReorderForNonzeroDiagonal() 324 M*/ 325 326 EXTERN_C_BEGIN 327 #undef __FUNCT__ 328 #define __FUNCT__ "PCCreate_LU" 329 PetscErrorCode PETSCKSP_DLLEXPORT PCCreate_LU(PC pc) 330 { 331 PetscErrorCode ierr; 332 PetscMPIInt size; 333 PC_LU *dir; 334 335 PetscFunctionBegin; 336 ierr = PetscNewLog(pc,PC_LU,&dir);CHKERRQ(ierr); 337 338 ierr = MatFactorInfoInitialize(&((PC_Factor*)dir)->info);CHKERRQ(ierr); 339 ((PC_Factor*)dir)->fact = 0; 340 dir->inplace = PETSC_FALSE; 341 dir->nonzerosalongdiagonal = PETSC_FALSE; 342 343 ((PC_Factor*)dir)->info.fill = 5.0; 344 ((PC_Factor*)dir)->info.dtcol = 1.e-6; /* default to pivoting; this is only thing PETSc LU supports */ 345 ((PC_Factor*)dir)->info.shiftnz = 0.0; 346 ((PC_Factor*)dir)->info.zeropivot = 1.e-12; 347 ((PC_Factor*)dir)->info.pivotinblocks = 1.0; 348 ((PC_Factor*)dir)->info.shiftpd = 0.0; /* false */ 349 dir->col = 0; 350 dir->row = 0; 351 352 ierr = PetscStrallocpy(MAT_SOLVER_PETSC,&((PC_Factor*)dir)->solvertype);CHKERRQ(ierr); 353 ierr = MPI_Comm_size(((PetscObject)pc)->comm,&size);CHKERRQ(ierr); 354 if (size == 1) { 355 ierr = PetscStrallocpy(MATORDERING_ND,&((PC_Factor*)dir)->ordering);CHKERRQ(ierr); 356 } else { 357 ierr = PetscStrallocpy(MATORDERING_NATURAL,&((PC_Factor*)dir)->ordering);CHKERRQ(ierr); 358 } 359 dir->reusefill = PETSC_FALSE; 360 dir->reuseordering = PETSC_FALSE; 361 pc->data = (void*)dir; 362 363 pc->ops->destroy = PCDestroy_LU; 364 pc->ops->apply = PCApply_LU; 365 pc->ops->applytranspose = PCApplyTranspose_LU; 366 pc->ops->setup = PCSetUp_LU; 367 pc->ops->setfromoptions = PCSetFromOptions_LU; 368 pc->ops->view = PCView_LU; 369 pc->ops->applyrichardson = 0; 370 pc->ops->getfactoredmatrix = PCFactorGetMatrix_Factor; 371 372 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorGetMatSolverPackage_C","PCFactorGetMatSolverPackage_Factor", 373 PCFactorGetMatSolverPackage_Factor);CHKERRQ(ierr); 374 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetMatSolverPackage_C","PCFactorSetMatSolverPackage_Factor", 375 PCFactorSetMatSolverPackage_Factor);CHKERRQ(ierr); 376 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetZeroPivot_C","PCFactorSetZeroPivot_Factor", 377 PCFactorSetZeroPivot_Factor);CHKERRQ(ierr); 378 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftNonzero_C","PCFactorSetShiftNonzero_Factor", 379 PCFactorSetShiftNonzero_Factor);CHKERRQ(ierr); 380 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetShiftPd_C","PCFactorSetShiftPd_Factor", 381 PCFactorSetShiftPd_Factor);CHKERRQ(ierr); 382 383 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetFill_C","PCFactorSetFill_Factor", 384 PCFactorSetFill_Factor);CHKERRQ(ierr); 385 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetUseInPlace_C","PCFactorSetUseInPlace_LU", 386 PCFactorSetUseInPlace_LU);CHKERRQ(ierr); 387 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetMatOrderingType_C","PCFactorSetMatOrderingType_Factor", 388 PCFactorSetMatOrderingType_Factor);CHKERRQ(ierr); 389 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseOrdering_C","PCFactorSetReuseOrdering_LU", 390 PCFactorSetReuseOrdering_LU);CHKERRQ(ierr); 391 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetReuseFill_C","PCFactorSetReuseFill_LU", 392 PCFactorSetReuseFill_LU);CHKERRQ(ierr); 393 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetPivoting_C","PCFactorSetPivoting_Factor", 394 PCFactorSetPivoting_Factor);CHKERRQ(ierr); 395 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorSetPivotInBlocks_C","PCFactorSetPivotInBlocks_Factor", 396 PCFactorSetPivotInBlocks_Factor);CHKERRQ(ierr); 397 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFactorReorderForNonzeroDiagonal_C","PCFactorReorderForNonzeroDiagonal_LU", 398 PCFactorReorderForNonzeroDiagonal_LU);CHKERRQ(ierr); 399 PetscFunctionReturn(0); 400 } 401 EXTERN_C_END 402