1 2 #include "src/ksp/pc/pcimpl.h" /*I "petscpc.h" I*/ 3 4 /* Options Database Keys: ??? 5 . -pc_ilu_damping - add damping to diagonal to prevent zero (or very small) pivots 6 . -pc_ilu_shift - apply Manteuffel shift to diagonal to force positive definite preconditioner 7 . -pc_ilu_zeropivot <tol> - set tolerance for what is considered a zero pivot 8 */ 9 10 #undef __FUNCT__ 11 #define __FUNCT__ "PCFactorSetZeroPivot" 12 /*@ 13 PCFactorSetZeroPivot - Sets the size at which smaller pivots are declared to be zero 14 15 Collective on PC 16 17 Input Parameters: 18 + pc - the preconditioner context 19 - zero - all pivots smaller than this will be considered zero 20 21 Options Database Key: 22 . -pc_factor_zeropivot <zero> - Sets tolerance for what is considered a zero pivot 23 24 Level: intermediate 25 26 .keywords: PC, set, factorization, direct, fill 27 28 .seealso: PCFactorSetShiftNonzero(), PCFactorSetShiftPd() 29 @*/ 30 PetscErrorCode PCFactorSetZeroPivot(PC pc,PetscReal zero) 31 { 32 PetscErrorCode ierr,(*f)(PC,PetscReal); 33 34 PetscFunctionBegin; 35 PetscValidHeaderSpecific(pc,PC_COOKIE,1); 36 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetZeroPivot_C",(void (**)(void))&f);CHKERRQ(ierr); 37 if (f) { 38 ierr = (*f)(pc,zero);CHKERRQ(ierr); 39 } 40 PetscFunctionReturn(0); 41 } 42 43 #undef __FUNCT__ 44 #define __FUNCT__ "PCFactorSetShiftNonzero" 45 /*@ 46 PCFactorSetShiftNonzero - adds this quantity to the diagonal of the matrix during 47 numerical factorization, thus the matrix has nonzero pivots 48 49 Collective on PC 50 51 Input Parameters: 52 + pc - the preconditioner context 53 - shift - amount of shift 54 55 Options Database Key: 56 . -pc_factor_shift_nonzero <shift> - Sets shift amount or PETSC_DECIDE for the default 57 58 Note: If 0.0 is given, then no shift is used. If a diagonal element is classified as a zero 59 pivot, then the shift is doubled until this is alleviated. 60 61 Level: intermediate 62 63 .keywords: PC, set, factorization, direct, fill 64 65 .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftPd() 66 @*/ 67 PetscErrorCode PCFactorSetShiftNonzero(PC pc,PetscReal shift) 68 { 69 PetscErrorCode ierr,(*f)(PC,PetscReal); 70 71 PetscFunctionBegin; 72 PetscValidHeaderSpecific(pc,PC_COOKIE,1); 73 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetShiftNonzero_C",(void (**)(void))&f);CHKERRQ(ierr); 74 if (f) { 75 ierr = (*f)(pc,shift);CHKERRQ(ierr); 76 } 77 PetscFunctionReturn(0); 78 } 79 80 #undef __FUNCT__ 81 #define __FUNCT__ "PCFactorSetShiftPd" 82 /*@ 83 PCFactorSetShiftPd - specify whether to use Manteuffel shifting. 84 If a matrix factorisation breaks down because of nonpositive pivots, 85 adding sufficient identity to the diagonal will remedy this. 86 Setting this causes a bisection method to find the minimum shift that 87 will lead to a well-defined matrix factor. 88 89 Collective on PC 90 91 Input parameters: 92 + pc - the preconditioner context 93 - shifting - PETSC_TRUE to set shift else PETSC_FALSE 94 95 Options Database Key: 96 . -pc_factor_shift_positive_definite [PETSC_TRUE/PETSC_FALSE] - Activate/Deactivate PCFactorSetShiftPd(); the value 97 is optional with PETSC_TRUE being the default 98 99 Level: intermediate 100 101 .keywords: PC, indefinite, factorization 102 103 .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftNonzero() 104 @*/ 105 PetscErrorCode PCFactorSetShiftPd(PC pc,PetscTruth shift) 106 { 107 PetscErrorCode ierr,(*f)(PC,PetscTruth); 108 109 PetscFunctionBegin; 110 PetscValidHeaderSpecific(pc,PC_COOKIE,1); 111 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetShiftPd_C",(void (**)(void))&f);CHKERRQ(ierr); 112 if (f) { 113 ierr = (*f)(pc,shift);CHKERRQ(ierr); 114 } 115 PetscFunctionReturn(0); 116 } 117