1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
10: /*
11: Implements the ST class for preconditioned eigenvalue methods
12: */
14: #include <slepc/private/stimpl.h> /*I "slepcst.h" I*/
16: typedef struct {
17: PetscBool setmat;
18: } ST_PRECOND;
20: static PetscErrorCode STSetDefaultKSP_Precond(ST st) 21: {
23: PC pc;
24: PCType pctype;
25: Mat P;
26: PetscBool t0,t1;
29: KSPGetPC(st->ksp,&pc);
30: PetscObjectGetType((PetscObject)pc,&pctype);
31: STPrecondGetMatForPC(st,&P);
32: if (!pctype && st->A && st->A[0]) {
33: if (P || st->shift_matrix == ST_MATMODE_SHELL) {
34: PCSetType(pc,PCJACOBI);
35: } else {
36: MatHasOperation(st->A[0],MATOP_DUPLICATE,&t0);
37: if (st->nmat>1) {
38: MatHasOperation(st->A[0],MATOP_AXPY,&t1);
39: } else t1 = PETSC_TRUE;
40: PCSetType(pc,(t0 && t1)?PCJACOBI:PCNONE);
41: }
42: }
43: KSPSetErrorIfNotConverged(st->ksp,PETSC_FALSE);
44: return(0);
45: }
47: PetscErrorCode STSetUp_Precond(ST st) 48: {
49: Mat P;
50: PC pc;
51: PetscBool t0,setmat,destroyP=PETSC_FALSE,builtP;
55: /* if the user did not set the shift, use the target value */
56: if (!st->sigma_set) st->sigma = st->defsigma;
58: /* If either pc is none and no matrix has to be set, or pc is shell , exit */
59: if (!st->ksp) { STGetKSP(st,&st->ksp); }
60: KSPGetPC(st->ksp,&pc);
61: PetscObjectTypeCompare((PetscObject)pc,PCSHELL,&t0);
62: if (t0) return(0);
63: PetscObjectTypeCompare((PetscObject)pc,PCNONE,&t0);
64: STPrecondGetKSPHasMat(st,&setmat);
65: if (t0 && !setmat) return(0);
67: /* Check if a user matrix is set */
68: STPrecondGetMatForPC(st,&P);
70: /* If not, create A - shift*B */
71: if (P) {
72: builtP = PETSC_FALSE;
73: destroyP = PETSC_TRUE;
74: PetscObjectReference((PetscObject)P);
75: } else {
76: builtP = PETSC_TRUE;
78: if (!(PetscAbsScalar(st->sigma) < PETSC_MAX_REAL) && st->nmat>1) {
79: P = st->A[1];
80: destroyP = PETSC_FALSE;
81: } else if (st->sigma == 0.0) {
82: P = st->A[0];
83: destroyP = PETSC_FALSE;
84: } else if (PetscAbsScalar(st->sigma) < PETSC_MAX_REAL && st->shift_matrix != ST_MATMODE_SHELL) {
85: if (st->shift_matrix == ST_MATMODE_INPLACE) {
86: P = st->A[0];
87: destroyP = PETSC_FALSE;
88: } else {
89: MatDuplicate(st->A[0],MAT_COPY_VALUES,&P);
90: destroyP = PETSC_TRUE;
91: }
92: if (st->nmat>1) {
93: MatAXPY(P,-st->sigma,st->A[1],st->str);
94: } else {
95: MatShift(P,-st->sigma);
96: }
97: /* TODO: in case of ST_MATMODE_INPLACE should keep the Hermitian flag of st->A and restore at the end */
98: STMatSetHermitian(st,P);
99: } else builtP = PETSC_FALSE;
100: }
102: /* If P was not possible to obtain, set pc to PCNONE */
103: if (!P) {
104: PCSetType(pc,PCNONE);
106: /* If some matrix has to be set to ksp, a shell matrix is created */
107: if (setmat) {
108: STMatShellCreate(st,-st->sigma,0,NULL,NULL,&P);
109: STMatSetHermitian(st,P);
110: destroyP = PETSC_TRUE;
111: }
112: }
114: KSPSetOperators(st->ksp,setmat?P:NULL,P);
116: if (destroyP) {
117: MatDestroy(&P);
118: } else if (st->shift_matrix == ST_MATMODE_INPLACE && builtP) {
119: if (st->sigma != 0.0 && PetscAbsScalar(st->sigma) < PETSC_MAX_REAL) {
120: if (st->nmat>1) {
121: MatAXPY(st->A[0],st->sigma,st->A[1],st->str);
122: } else {
123: MatShift(st->A[0],st->sigma);
124: }
125: }
126: }
127: return(0);
128: }
130: PetscErrorCode STSetShift_Precond(ST st,PetscScalar newshift)131: {
135: /* Nothing to be done if STSetUp has not been called yet */
136: if (!st->state) return(0);
137: st->sigma = newshift;
138: if (st->shift_matrix != ST_MATMODE_SHELL) {
139: STSetUp_Precond(st);
140: }
141: return(0);
142: }
144: static PetscErrorCode STPrecondGetMatForPC_Precond(ST st,Mat *mat)145: {
147: PC pc;
148: PetscBool flag;
151: if (!st->ksp) { STGetKSP(st,&st->ksp); }
152: KSPGetPC(st->ksp,&pc);
153: PCGetOperatorsSet(pc,NULL,&flag);
154: if (flag) {
155: PCGetOperators(pc,NULL,mat);
156: } else *mat = NULL;
157: return(0);
158: }
160: /*@
161: STPrecondGetMatForPC - Returns the matrix previously set by STPrecondSetMatForPC().
163: Not Collective, but the Mat is shared by all processors that share the ST165: Input Parameter:
166: . st - the spectral transformation context
168: Output Parameter:
169: . mat - the matrix that will be used in constructing the preconditioner or
170: NULL if no matrix was set by STPrecondSetMatForPC().
172: Level: advanced
174: .seealso: STPrecondSetMatForPC()
175: @*/
176: PetscErrorCode STPrecondGetMatForPC(ST st,Mat *mat)177: {
183: PetscUseMethod(st,"STPrecondGetMatForPC_C",(ST,Mat*),(st,mat));
184: return(0);
185: }
187: static PetscErrorCode STPrecondSetMatForPC_Precond(ST st,Mat mat)188: {
189: PC pc;
190: Mat A;
191: PetscBool flag;
195: if (!st->ksp) { STGetKSP(st,&st->ksp); }
196: KSPGetPC(st->ksp,&pc);
197: /* Yes, all these lines are needed to safely set mat as the preconditioner
198: matrix in pc */
199: PCGetOperatorsSet(pc,&flag,NULL);
200: if (flag) {
201: PCGetOperators(pc,&A,NULL);
202: PetscObjectReference((PetscObject)A);
203: } else A = NULL;
204: PetscObjectReference((PetscObject)mat);
205: PCSetOperators(pc,A,mat);
206: MatDestroy(&A);
207: MatDestroy(&mat);
208: STPrecondSetKSPHasMat(st,PETSC_TRUE);
209: return(0);
210: }
212: /*@
213: STPrecondSetMatForPC - Sets the matrix that must be used to build the preconditioner.
215: Logically Collective on ST and Mat
217: Input Parameter:
218: + st - the spectral transformation context
219: - mat - the matrix that will be used in constructing the preconditioner
221: Level: advanced
223: Notes:
224: This matrix will be passed to the KSP object (via KSPSetOperators) as
225: the matrix to be used when constructing the preconditioner.
226: If no matrix is set or mat is set to NULL, A - sigma*B will
227: be used to build the preconditioner, being sigma the value set by STSetShift().
229: .seealso: STPrecondSetMatForPC(), STSetShift()
230: @*/
231: PetscErrorCode STPrecondSetMatForPC(ST st,Mat mat)232: {
239: PetscTryMethod(st,"STPrecondSetMatForPC_C",(ST,Mat),(st,mat));
240: return(0);
241: }
243: static PetscErrorCode STPrecondSetKSPHasMat_Precond(ST st,PetscBool setmat)244: {
245: ST_PRECOND *data = (ST_PRECOND*)st->data;
248: data->setmat = setmat;
249: return(0);
250: }
252: /*@
253: STPrecondSetKSPHasMat - Sets a flag indicating that during STSetUp the coefficient
254: matrix of the KSP linear system (A) must be set to be the same matrix as the
255: preconditioner (P).
257: Collective on ST259: Input Parameter:
260: + st - the spectral transformation context
261: - setmat - the flag
263: Notes:
264: In most cases, the preconditioner matrix is used only in the PC object, but
265: in external solvers this matrix must be provided also as the A-matrix in
266: the KSP object.
268: Level: developer
270: .seealso: STPrecondGetKSPHasMat(), STSetShift()
271: @*/
272: PetscErrorCode STPrecondSetKSPHasMat(ST st,PetscBool setmat)273: {
279: PetscTryMethod(st,"STPrecondSetKSPHasMat_C",(ST,PetscBool),(st,setmat));
280: return(0);
281: }
283: static PetscErrorCode STPrecondGetKSPHasMat_Precond(ST st,PetscBool *setmat)284: {
285: ST_PRECOND *data = (ST_PRECOND*)st->data;
288: *setmat = data->setmat;
289: return(0);
290: }
292: /*@
293: STPrecondGetKSPHasMat - Returns the flag indicating if the coefficient
294: matrix of the KSP linear system (A) is set to be the same matrix as the
295: preconditioner (P).
297: Not Collective
299: Input Parameter:
300: . st - the spectral transformation context
302: Output Parameter:
303: . setmat - the flag
305: Level: developer
307: .seealso: STPrecondSetKSPHasMat(), STSetShift()
308: @*/
309: PetscErrorCode STPrecondGetKSPHasMat(ST st,PetscBool *setmat)310: {
316: PetscUseMethod(st,"STPrecondGetKSPHasMat_C",(ST,PetscBool*),(st,setmat));
317: return(0);
318: }
320: PetscErrorCode STDestroy_Precond(ST st)321: {
325: PetscFree(st->data);
326: PetscObjectComposeFunction((PetscObject)st,"STPrecondGetMatForPC_C",NULL);
327: PetscObjectComposeFunction((PetscObject)st,"STPrecondSetMatForPC_C",NULL);
328: PetscObjectComposeFunction((PetscObject)st,"STPrecondGetKSPHasMat_C",NULL);
329: PetscObjectComposeFunction((PetscObject)st,"STPrecondSetKSPHasMat_C",NULL);
330: return(0);
331: }
333: SLEPC_EXTERN PetscErrorCode STCreate_Precond(ST st)334: {
336: ST_PRECOND *ctx;
339: PetscNewLog(st,&ctx);
340: st->data = (void*)ctx;
342: st->ops->getbilinearform = STGetBilinearForm_Default;
343: st->ops->setup = STSetUp_Precond;
344: st->ops->setshift = STSetShift_Precond;
345: st->ops->destroy = STDestroy_Precond;
346: st->ops->setdefaultksp = STSetDefaultKSP_Precond;
348: PetscObjectComposeFunction((PetscObject)st,"STPrecondGetMatForPC_C",STPrecondGetMatForPC_Precond);
349: PetscObjectComposeFunction((PetscObject)st,"STPrecondSetMatForPC_C",STPrecondSetMatForPC_Precond);
350: PetscObjectComposeFunction((PetscObject)st,"STPrecondGetKSPHasMat_C",STPrecondGetKSPHasMat_Precond);
351: PetscObjectComposeFunction((PetscObject)st,"STPrecondSetKSPHasMat_C",STPrecondSetKSPHasMat_Precond);
353: STPrecondSetKSPHasMat_Precond(st,PETSC_TRUE);
354: return(0);
355: }