1 #define PETSCKSP_DLL 2 3 #include "../src/ksp/pc/impls/factor/factor.h" /*I "petscpc.h" I*/ 4 5 #undef __FUNCT__ 6 #define __FUNCT__ "PCFactorSetUpMatSolverPackage" 7 /*@ 8 PCFactorSetUpMatSolverPackage - Can be called after KSPSetOperators() or PCSetOperators(), causes MatGetFactor() to be called so then one may 9 set the options for that particular factorization object. 10 11 Input Parameter: 12 . pc - the preconditioner context 13 14 Notes: After you have called this function (which has to be after the KSPSetOperators() or PCSetOperators()) you can call PCFactorGetMatrix() and then set factor options on that matrix. 15 16 .seealso: PCFactorSetMatSolverPackage(), PCFactorGetMatrix() 17 18 @*/ 19 PetscErrorCode PCFactorSetUpMatSolverPackage(PC pc) 20 { 21 PetscErrorCode ierr; 22 23 PetscFunctionBegin; 24 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 25 ierr = PetscTryMethod(pc,"PCFactorSetUpMatSolverPackage_C",(PC),(pc));CHKERRQ(ierr); 26 PetscFunctionReturn(0); 27 } 28 29 #undef __FUNCT__ 30 #define __FUNCT__ "PCFactorSetZeroPivot" 31 /*@ 32 PCFactorSetZeroPivot - Sets the size at which smaller pivots are declared to be zero 33 34 Logically Collective on PC 35 36 Input Parameters: 37 + pc - the preconditioner context 38 - zero - all pivots smaller than this will be considered zero 39 40 Options Database Key: 41 . -pc_factor_zeropivot <zero> - Sets tolerance for what is considered a zero pivot 42 43 Level: intermediate 44 45 .keywords: PC, set, factorization, direct, fill 46 47 .seealso: PCFactorSetShiftType(), PCFactorSetShiftAmount() 48 @*/ 49 PetscErrorCode PCFactorSetZeroPivot(PC pc,PetscReal zero) 50 { 51 PetscErrorCode ierr; 52 53 PetscFunctionBegin; 54 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 55 PetscValidLogicalCollectiveReal(pc,zero,2); 56 ierr = PetscTryMethod(pc,"PCFactorSetZeroPivot_C",(PC,PetscReal),(pc,zero));CHKERRQ(ierr); 57 PetscFunctionReturn(0); 58 } 59 60 #undef __FUNCT__ 61 #define __FUNCT__ "PCFactorSetShiftType" 62 /*@ 63 PCFactorSetShiftType - adds a particular type of quantity to the diagonal of the matrix during 64 numerical factorization, thus the matrix has nonzero pivots 65 66 Logically Collective on PC 67 68 Input Parameters: 69 + pc - the preconditioner context 70 - shifttype - type of shift; one of MAT_SHIFT_NONE, MAT_SHIFT_NONZERO, MAT_SHIFT_POSITIVE_DEFINITE, MAT_SHIFT_INBLOCKS 71 72 Options Database Key: 73 . -pc_factor_shift_type <shifttype> - Sets shift type or PETSC_DECIDE for the default; use '-help' for a list of available types 74 75 Level: intermediate 76 77 .keywords: PC, set, factorization, 78 79 .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftAmount() 80 @*/ 81 PetscErrorCode PCFactorSetShiftType(PC pc,MatFactorShiftType shifttype) 82 { 83 PetscErrorCode ierr; 84 85 PetscFunctionBegin; 86 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 87 PetscValidLogicalCollectiveEnum(pc,shifttype,2); 88 ierr = PetscTryMethod(pc,"PCFactorSetShiftType_C",(PC,MatFactorShiftType),(pc,shifttype));CHKERRQ(ierr); 89 PetscFunctionReturn(0); 90 } 91 92 #undef __FUNCT__ 93 #define __FUNCT__ "PCFactorSetShiftAmount" 94 /*@ 95 PCFactorSetShiftAmount - adds a quantity to the diagonal of the matrix during 96 numerical factorization, thus the matrix has nonzero pivots 97 98 Logically Collective on PC 99 100 Input Parameters: 101 + pc - the preconditioner context 102 - shiftamount - amount of shift 103 104 Options Database Key: 105 . -pc_factor_shift_amount <shiftamount> - Sets shift amount or PETSC_DECIDE for the default 106 107 Level: intermediate 108 109 .keywords: PC, set, factorization, 110 111 .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftType() 112 @*/ 113 PetscErrorCode PCFactorSetShiftAmount(PC pc,PetscReal shiftamount) 114 { 115 PetscErrorCode ierr; 116 117 PetscFunctionBegin; 118 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 119 PetscValidLogicalCollectiveReal(pc,shiftamount,2); 120 ierr = PetscTryMethod(pc,"PCFactorSetShiftAmount_C",(PC,PetscReal),(pc,shiftamount));CHKERRQ(ierr); 121 PetscFunctionReturn(0); 122 } 123 124 #undef __FUNCT__ 125 #define __FUNCT__ "PCFactorSetDropTolerance" 126 /* 127 PCFactorSetDropTolerance - The preconditioner will use an ILU 128 based on a drop tolerance. (Under development) 129 130 Logically Collective on PC 131 132 Input Parameters: 133 + pc - the preconditioner context 134 . dt - the drop tolerance, try from 1.e-10 to .1 135 . dtcol - tolerance for column pivot, good values [0.1 to 0.01] 136 - maxrowcount - the max number of nonzeros allowed in a row, best value 137 depends on the number of nonzeros in row of original matrix 138 139 Options Database Key: 140 . -pc_factor_drop_tolerance <dt,dtcol,maxrowcount> - Sets drop tolerance 141 142 Level: intermediate 143 144 There are NO default values for the 3 parameters, you must set them with reasonable values for your 145 matrix. We don't know how to compute reasonable values. 146 147 .keywords: PC, levels, reordering, factorization, incomplete, ILU 148 */ 149 PetscErrorCode PCFactorSetDropTolerance(PC pc,PetscReal dt,PetscReal dtcol,PetscInt maxrowcount) 150 { 151 PetscErrorCode ierr; 152 153 PetscFunctionBegin; 154 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 155 PetscValidLogicalCollectiveReal(pc,dtcol,2); 156 PetscValidLogicalCollectiveInt(pc,maxrowcount,3); 157 ierr = PetscTryMethod(pc,"PCFactorSetDropTolerance_C",(PC,PetscReal,PetscReal,PetscInt),(pc,dt,dtcol,maxrowcount));CHKERRQ(ierr); 158 PetscFunctionReturn(0); 159 } 160 161 #undef __FUNCT__ 162 #define __FUNCT__ "PCFactorSetLevels" 163 /*@ 164 PCFactorSetLevels - Sets the number of levels of fill to use. 165 166 Logically Collective on PC 167 168 Input Parameters: 169 + pc - the preconditioner context 170 - levels - number of levels of fill 171 172 Options Database Key: 173 . -pc_factor_levels <levels> - Sets fill level 174 175 Level: intermediate 176 177 .keywords: PC, levels, fill, factorization, incomplete, ILU 178 @*/ 179 PetscErrorCode PCFactorSetLevels(PC pc,PetscInt levels) 180 { 181 PetscErrorCode ierr; 182 183 PetscFunctionBegin; 184 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 185 if (levels < 0) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_OUTOFRANGE,"negative levels"); 186 PetscValidLogicalCollectiveInt(pc,levels,2); 187 ierr = PetscTryMethod(pc,"PCFactorSetLevels_C",(PC,PetscInt),(pc,levels));CHKERRQ(ierr); 188 PetscFunctionReturn(0); 189 } 190 191 #undef __FUNCT__ 192 #define __FUNCT__ "PCFactorSetAllowDiagonalFill" 193 /*@ 194 PCFactorSetAllowDiagonalFill - Causes all diagonal matrix entries to be 195 treated as level 0 fill even if there is no non-zero location. 196 197 Logically Collective on PC 198 199 Input Parameters: 200 + pc - the preconditioner context 201 202 Options Database Key: 203 . -pc_factor_diagonal_fill 204 205 Notes: 206 Does not apply with 0 fill. 207 208 Level: intermediate 209 210 .keywords: PC, levels, fill, factorization, incomplete, ILU 211 @*/ 212 PetscErrorCode PCFactorSetAllowDiagonalFill(PC pc) 213 { 214 PetscErrorCode ierr; 215 216 PetscFunctionBegin; 217 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 218 ierr = PetscTryMethod(pc,"PCFactorSetAllowDiagonalFill_C",(PC),(pc));CHKERRQ(ierr); 219 PetscFunctionReturn(0); 220 } 221 222 #undef __FUNCT__ 223 #define __FUNCT__ "PCFactorReorderForNonzeroDiagonal" 224 /*@ 225 PCFactorReorderForNonzeroDiagonal - reorders rows/columns of matrix to remove zeros from diagonal 226 227 Logically Collective on PC 228 229 Input Parameters: 230 + pc - the preconditioner context 231 - tol - diagonal entries smaller than this in absolute value are considered zero 232 233 Options Database Key: 234 . -pc_factor_nonzeros_along_diagonal 235 236 Level: intermediate 237 238 .keywords: PC, set, factorization, direct, fill 239 240 .seealso: PCFactorSetFill(), PCFactorSetShiftNonzero(), PCFactorSetZeroPivot(), MatReorderForNonzeroDiagonal() 241 @*/ 242 PetscErrorCode PCFactorReorderForNonzeroDiagonal(PC pc,PetscReal rtol) 243 { 244 PetscErrorCode ierr; 245 246 PetscFunctionBegin; 247 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 248 PetscValidLogicalCollectiveReal(pc,rtol,2); 249 ierr = PetscTryMethod(pc,"PCFactorReorderForNonzeroDiagonal_C",(PC,PetscReal),(pc,rtol));CHKERRQ(ierr); 250 PetscFunctionReturn(0); 251 } 252 253 #undef __FUNCT__ 254 #define __FUNCT__ "PCFactorSetMatSolverPackage" 255 /*@C 256 PCFactorSetMatSolverPackage - sets the software that is used to perform the factorization 257 258 Logically Collective on PC 259 260 Input Parameters: 261 + pc - the preconditioner context 262 - stype - for example, spooles, superlu, superlu_dist 263 264 Options Database Key: 265 . -pc_factor_mat_solver_package <stype> - spooles, petsc, superlu, superlu_dist, mumps 266 267 Level: intermediate 268 269 Note: 270 By default this will use the PETSc factorization if it exists 271 272 273 .keywords: PC, set, factorization, direct, fill 274 275 .seealso: MatGetFactor(), MatSolverPackage, PCFactorGetMatSolverPackage() 276 277 @*/ 278 PetscErrorCode PCFactorSetMatSolverPackage(PC pc,const MatSolverPackage stype) 279 { 280 PetscErrorCode ierr; 281 282 PetscFunctionBegin; 283 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 284 ierr = PetscTryMethod(pc,"PCFactorSetMatSolverPackage_C",(PC,const MatSolverPackage),(pc,stype));CHKERRQ(ierr); 285 PetscFunctionReturn(0); 286 } 287 288 #undef __FUNCT__ 289 #define __FUNCT__ "PCFactorGetMatSolverPackage" 290 /*@C 291 PCFactorGetMatSolverPackage - gets the software that is used to perform the factorization 292 293 Not Collective 294 295 Input Parameter: 296 . pc - the preconditioner context 297 298 Output Parameter: 299 . stype - for example, spooles, superlu, superlu_dist 300 301 Level: intermediate 302 303 304 .keywords: PC, set, factorization, direct, fill 305 306 .seealso: MatGetFactor(), MatSolverPackage, PCFactorGetMatSolverPackage() 307 308 @*/ 309 PetscErrorCode PCFactorGetMatSolverPackage(PC pc,const MatSolverPackage *stype) 310 { 311 PetscErrorCode ierr; 312 313 PetscFunctionBegin; 314 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 315 ierr = PetscTryMethod(pc,"PCFactorGetMatSolverPackage_C",(PC,const MatSolverPackage*),(pc,stype));CHKERRQ(ierr); 316 PetscFunctionReturn(0); 317 } 318 319 #undef __FUNCT__ 320 #define __FUNCT__ "PCFactorSetFill" 321 /*@ 322 PCFactorSetFill - Indicate the amount of fill you expect in the factored matrix, 323 fill = number nonzeros in factor/number nonzeros in original matrix. 324 325 Not Collective, each process can expect a different amount of fill 326 327 Input Parameters: 328 + pc - the preconditioner context 329 - fill - amount of expected fill 330 331 Options Database Key: 332 . -pc_factor_fill <fill> - Sets fill amount 333 334 Level: intermediate 335 336 Note: 337 For sparse matrix factorizations it is difficult to predict how much 338 fill to expect. By running with the option -info PETSc will print the 339 actual amount of fill used; allowing you to set the value accurately for 340 future runs. Default PETSc uses a value of 5.0 341 342 .keywords: PC, set, factorization, direct, fill 343 344 @*/ 345 PetscErrorCode PCFactorSetFill(PC pc,PetscReal fill) 346 { 347 PetscErrorCode ierr; 348 349 PetscFunctionBegin; 350 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 351 if (fill < 1.0) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Fill factor cannot be less then 1.0"); 352 ierr = PetscTryMethod(pc,"PCFactorSetFill_C",(PC,PetscReal),(pc,fill));CHKERRQ(ierr); 353 PetscFunctionReturn(0); 354 } 355 356 #undef __FUNCT__ 357 #define __FUNCT__ "PCFactorSetUseInPlace" 358 /*@ 359 PCFactorSetUseInPlace - Tells the system to do an in-place factorization. 360 For dense matrices, this enables the solution of much larger problems. 361 For sparse matrices the factorization cannot be done truly in-place 362 so this does not save memory during the factorization, but after the matrix 363 is factored, the original unfactored matrix is freed, thus recovering that 364 space. 365 366 Logically Collective on PC 367 368 Input Parameters: 369 . pc - the preconditioner context 370 371 Options Database Key: 372 . -pc_factor_in_place - Activates in-place factorization 373 374 Notes: 375 PCFactorSetUseInplace() can only be used with the KSP method KSPPREONLY or when 376 a different matrix is provided for the multiply and the preconditioner in 377 a call to KSPSetOperators(). 378 This is because the Krylov space methods require an application of the 379 matrix multiplication, which is not possible here because the matrix has 380 been factored in-place, replacing the original matrix. 381 382 Level: intermediate 383 384 .keywords: PC, set, factorization, direct, inplace, in-place, LU 385 386 .seealso: PCILUSetUseInPlace() 387 @*/ 388 PetscErrorCode PCFactorSetUseInPlace(PC pc) 389 { 390 PetscErrorCode ierr; 391 392 PetscFunctionBegin; 393 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 394 ierr = PetscTryMethod(pc,"PCFactorSetUseInPlace_C",(PC),(pc));CHKERRQ(ierr); 395 PetscFunctionReturn(0); 396 } 397 398 #undef __FUNCT__ 399 #define __FUNCT__ "PCFactorSetMatOrderingType" 400 /*@C 401 PCFactorSetMatOrderingType - Sets the ordering routine (to reduce fill) to 402 be used in the LU factorization. 403 404 Logically Collective on PC 405 406 Input Parameters: 407 + pc - the preconditioner context 408 - ordering - the matrix ordering name, for example, MATORDERINGND or MATORDERINGRCM 409 410 Options Database Key: 411 . -pc_factor_mat_ordering_type <nd,rcm,...> - Sets ordering routine 412 413 Level: intermediate 414 415 Notes: nested dissection is used by default 416 417 For Cholesky and ICC and the SBAIJ format reorderings are not available, 418 since only the upper triangular part of the matrix is stored. You can use the 419 SeqAIJ format in this case to get reorderings. 420 421 @*/ 422 PetscErrorCode PCFactorSetMatOrderingType(PC pc,const MatOrderingType ordering) 423 { 424 PetscErrorCode ierr; 425 426 PetscFunctionBegin; 427 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 428 ierr = PetscTryMethod(pc,"PCFactorSetMatOrderingType_C",(PC,const MatOrderingType),(pc,ordering));CHKERRQ(ierr); 429 PetscFunctionReturn(0); 430 } 431 432 #undef __FUNCT__ 433 #define __FUNCT__ "PCFactorSetColumnPivot" 434 /*@ 435 PCFactorSetColumnPivot - Determines when column pivoting is done during matrix factorization. 436 For PETSc dense matrices column pivoting is always done, for PETSc sparse matrices 437 it is never done. For the MATLAB and SuperLU factorization this is used. 438 439 Logically Collective on PC 440 441 Input Parameters: 442 + pc - the preconditioner context 443 - dtcol - 0.0 implies no pivoting, 1.0 complete pivoting (slower, requires more memory but more stable) 444 445 Options Database Key: 446 . -pc_factor_pivoting <dtcol> 447 448 Level: intermediate 449 450 .seealso: PCILUSetMatOrdering(), PCFactorSetPivotInBlocks() 451 @*/ 452 PetscErrorCode PCFactorSetColumnPivot(PC pc,PetscReal dtcol) 453 { 454 PetscErrorCode ierr; 455 456 PetscFunctionBegin; 457 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 458 PetscValidLogicalCollectiveReal(pc,dtcol,2); 459 ierr = PetscTryMethod(pc,"PCFactorSetColumnPivot_C",(PC,PetscReal),(pc,dtcol));CHKERRQ(ierr); 460 PetscFunctionReturn(0); 461 } 462 463 #undef __FUNCT__ 464 #define __FUNCT__ "PCFactorSetPivotInBlocks" 465 /*@ 466 PCFactorSetPivotInBlocks - Determines if pivoting is done while factoring each block 467 with BAIJ or SBAIJ matrices 468 469 Logically Collective on PC 470 471 Input Parameters: 472 + pc - the preconditioner context 473 - pivot - PETSC_TRUE or PETSC_FALSE 474 475 Options Database Key: 476 . -pc_factor_pivot_in_blocks <true,false> 477 478 Level: intermediate 479 480 .seealso: PCILUSetMatOrdering(), PCFactorSetColumnPivot() 481 @*/ 482 PetscErrorCode PCFactorSetPivotInBlocks(PC pc,PetscBool pivot) 483 { 484 PetscErrorCode ierr; 485 486 PetscFunctionBegin; 487 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 488 PetscValidLogicalCollectiveBool(pc,pivot,2); 489 ierr = PetscTryMethod(pc,"PCFactorSetPivotInBlocks_C",(PC,PetscBool),(pc,pivot));CHKERRQ(ierr); 490 PetscFunctionReturn(0); 491 } 492 493 #undef __FUNCT__ 494 #define __FUNCT__ "PCFactorSetReuseFill" 495 /*@ 496 PCFactorSetReuseFill - When matrices with same different nonzero structure are factored, 497 this causes later ones to use the fill ratio computed in the initial factorization. 498 499 Logically Collective on PC 500 501 Input Parameters: 502 + pc - the preconditioner context 503 - flag - PETSC_TRUE to reuse else PETSC_FALSE 504 505 Options Database Key: 506 . -pc_factor_reuse_fill - Activates PCFactorSetReuseFill() 507 508 Level: intermediate 509 510 .keywords: PC, levels, reordering, factorization, incomplete, Cholesky 511 512 .seealso: PCFactorSetReuseOrdering() 513 @*/ 514 PetscErrorCode PCFactorSetReuseFill(PC pc,PetscBool flag) 515 { 516 PetscErrorCode ierr; 517 518 PetscFunctionBegin; 519 PetscValidHeaderSpecific(pc,PC_CLASSID,2); 520 PetscValidLogicalCollectiveBool(pc,flag,2); 521 ierr = PetscTryMethod(pc,"PCFactorSetReuseFill_C",(PC,PetscBool),(pc,flag));CHKERRQ(ierr); 522 PetscFunctionReturn(0); 523 } 524