Actual source code: factor.c
1: #define PETSCKSP_DLL
3: #include private/pcimpl.h
7: /*@
8: PCFactorSetZeroPivot - Sets the size at which smaller pivots are declared to be zero
10: Collective on PC
11:
12: Input Parameters:
13: + pc - the preconditioner context
14: - zero - all pivots smaller than this will be considered zero
16: Options Database Key:
17: . -pc_factor_zeropivot <zero> - Sets tolerance for what is considered a zero pivot
19: Level: intermediate
21: .keywords: PC, set, factorization, direct, fill
23: .seealso: PCFactorSetShiftNonzero(), PCFactorSetShiftPd()
24: @*/
25: PetscErrorCode PCFactorSetZeroPivot(PC pc,PetscReal zero)
26: {
27: PetscErrorCode ierr,(*f)(PC,PetscReal);
31: PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetZeroPivot_C",(void (**)(void))&f);
32: if (f) {
33: (*f)(pc,zero);
34: }
35: return(0);
36: }
40: /*@
41: PCFactorSetShiftNonzero - adds this quantity to the diagonal of the matrix during
42: numerical factorization, thus the matrix has nonzero pivots
44: Collective on PC
45:
46: Input Parameters:
47: + pc - the preconditioner context
48: - shift - amount of shift
50: Options Database Key:
51: . -pc_factor_shift_nonzero <shift> - Sets shift amount or PETSC_DECIDE for the default
53: Note: If 0.0 is given, then no shift is used. If a diagonal element is classified as a zero
54: pivot, then the shift is doubled until this is alleviated.
56: Level: intermediate
58: .keywords: PC, set, factorization, direct, fill
60: .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftPd()
61: @*/
62: PetscErrorCode PCFactorSetShiftNonzero(PC pc,PetscReal shift)
63: {
64: PetscErrorCode ierr,(*f)(PC,PetscReal);
68: PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetShiftNonzero_C",(void (**)(void))&f);
69: if (f) {
70: (*f)(pc,shift);
71: }
72: return(0);
73: }
77: /*@
78: PCFactorSetShiftPd - specify whether to use Manteuffel shifting.
79: If a matrix factorisation breaks down because of nonpositive pivots,
80: adding sufficient identity to the diagonal will remedy this.
81: Setting this causes a bisection method to find the minimum shift that
82: will lead to a well-defined matrix factor.
84: Collective on PC
86: Input parameters:
87: + pc - the preconditioner context
88: - shifting - PETSC_TRUE to set shift else PETSC_FALSE
90: Options Database Key:
91: . -pc_factor_shift_positive_definite [PETSC_TRUE/PETSC_FALSE] - Activate/Deactivate PCFactorSetShiftPd(); the value
92: is optional with PETSC_TRUE being the default
94: Level: intermediate
96: .keywords: PC, indefinite, factorization
98: .seealso: PCFactorSetZeroPivot(), PCFactorSetShiftNonzero()
99: @*/
100: PetscErrorCode PCFactorSetShiftPd(PC pc,PetscTruth shift)
101: {
102: PetscErrorCode ierr,(*f)(PC,PetscTruth);
106: PetscObjectQueryFunction((PetscObject)pc,"PCFactorSetShiftPd_C",(void (**)(void))&f);
107: if (f) {
108: (*f)(pc,shift);
109: }
110: return(0);
111: }