Actual source code: sor.c
1: /*
2: Defines a (S)SOR preconditioner for any Mat implementation
3: */
4: #include <petsc/private/pcimpl.h>
6: typedef struct {
7: PetscInt its; /* inner iterations, number of sweeps */
8: PetscInt lits; /* local inner iterations, number of sweeps applied by the local matrix mat->A */
9: MatSORType sym; /* forward, reverse, symmetric etc. */
10: PetscReal omega;
11: PetscReal fshift;
12: } PC_SOR;
14: static PetscErrorCode PCDestroy_SOR(PC pc)
15: {
19: PetscFree(pc->data);
20: return(0);
21: }
23: static PetscErrorCode PCApply_SOR(PC pc,Vec x,Vec y)
24: {
25: PC_SOR *jac = (PC_SOR*)pc->data;
27: PetscInt flag = jac->sym | SOR_ZERO_INITIAL_GUESS;
30: MatSOR(pc->pmat,x,jac->omega,(MatSORType)flag,jac->fshift,jac->its,jac->lits,y);
31: MatFactorGetError(pc->pmat,(MatFactorError*)&pc->failedreason);
32: return(0);
33: }
35: static PetscErrorCode PCApplyTranspose_SOR(PC pc,Vec x,Vec y)
36: {
37: PC_SOR *jac = (PC_SOR*)pc->data;
39: PetscInt flag = jac->sym | SOR_ZERO_INITIAL_GUESS;
40: PetscBool set,sym;
43: MatIsSymmetricKnown(pc->pmat,&set,&sym);
44: if (!set || !sym || (jac->sym != SOR_SYMMETRIC_SWEEP && jac->sym != SOR_LOCAL_SYMMETRIC_SWEEP)) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Can only apply transpose of SOR if matrix is symmetric and sweep is symmetric");
45: MatSOR(pc->pmat,x,jac->omega,(MatSORType)flag,jac->fshift,jac->its,jac->lits,y);
46: MatFactorGetError(pc->pmat,(MatFactorError*)&pc->failedreason);
47: return(0);
48: }
50: static PetscErrorCode PCApplyRichardson_SOR(PC pc,Vec b,Vec y,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt its,PetscBool guesszero,PetscInt *outits,PCRichardsonConvergedReason *reason)
51: {
52: PC_SOR *jac = (PC_SOR*)pc->data;
54: MatSORType stype = jac->sym;
57: PetscInfo1(pc,"Warning, convergence critera ignored, using %D iterations\n",its);
58: if (guesszero) stype = (MatSORType) (stype | SOR_ZERO_INITIAL_GUESS);
59: MatSOR(pc->pmat,b,jac->omega,stype,jac->fshift,its*jac->its,jac->lits,y);
60: MatFactorGetError(pc->pmat,(MatFactorError*)&pc->failedreason);
61: *outits = its;
62: *reason = PCRICHARDSON_CONVERGED_ITS;
63: return(0);
64: }
66: PetscErrorCode PCSetFromOptions_SOR(PetscOptionItems *PetscOptionsObject,PC pc)
67: {
68: PC_SOR *jac = (PC_SOR*)pc->data;
70: PetscBool flg;
73: PetscOptionsHead(PetscOptionsObject,"(S)SOR options");
74: PetscOptionsReal("-pc_sor_omega","relaxation factor (0 < omega < 2)","PCSORSetOmega",jac->omega,&jac->omega,NULL);
75: PetscOptionsReal("-pc_sor_diagonal_shift","Add to the diagonal entries","",jac->fshift,&jac->fshift,NULL);
76: PetscOptionsInt("-pc_sor_its","number of inner SOR iterations","PCSORSetIterations",jac->its,&jac->its,NULL);
77: PetscOptionsInt("-pc_sor_lits","number of local inner SOR iterations","PCSORSetIterations",jac->lits,&jac->lits,NULL);
78: PetscOptionsBoolGroupBegin("-pc_sor_symmetric","SSOR, not SOR","PCSORSetSymmetric",&flg);
79: if (flg) {PCSORSetSymmetric(pc,SOR_SYMMETRIC_SWEEP);}
80: PetscOptionsBoolGroup("-pc_sor_backward","use backward sweep instead of forward","PCSORSetSymmetric",&flg);
81: if (flg) {PCSORSetSymmetric(pc,SOR_BACKWARD_SWEEP);}
82: PetscOptionsBoolGroup("-pc_sor_forward","use forward sweep","PCSORSetSymmetric",&flg);
83: if (flg) {PCSORSetSymmetric(pc,SOR_FORWARD_SWEEP);}
84: PetscOptionsBoolGroup("-pc_sor_local_symmetric","use SSOR separately on each processor","PCSORSetSymmetric",&flg);
85: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_SYMMETRIC_SWEEP);}
86: PetscOptionsBoolGroup("-pc_sor_local_backward","use backward sweep locally","PCSORSetSymmetric",&flg);
87: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_BACKWARD_SWEEP);}
88: PetscOptionsBoolGroupEnd("-pc_sor_local_forward","use forward sweep locally","PCSORSetSymmetric",&flg);
89: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_FORWARD_SWEEP);}
90: PetscOptionsTail();
91: return(0);
92: }
94: PetscErrorCode PCView_SOR(PC pc,PetscViewer viewer)
95: {
96: PC_SOR *jac = (PC_SOR*)pc->data;
97: MatSORType sym = jac->sym;
98: const char *sortype;
100: PetscBool iascii;
103: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
104: if (iascii) {
105: if (sym & SOR_ZERO_INITIAL_GUESS) {PetscViewerASCIIPrintf(viewer," zero initial guess\n");}
106: if (sym == SOR_APPLY_UPPER) sortype = "apply_upper";
107: else if (sym == SOR_APPLY_LOWER) sortype = "apply_lower";
108: else if (sym & SOR_EISENSTAT) sortype = "Eisenstat";
109: else if ((sym & SOR_SYMMETRIC_SWEEP) == SOR_SYMMETRIC_SWEEP) sortype = "symmetric";
110: else if (sym & SOR_BACKWARD_SWEEP) sortype = "backward";
111: else if (sym & SOR_FORWARD_SWEEP) sortype = "forward";
112: else if ((sym & SOR_LOCAL_SYMMETRIC_SWEEP) == SOR_LOCAL_SYMMETRIC_SWEEP) sortype = "local_symmetric";
113: else if (sym & SOR_LOCAL_FORWARD_SWEEP) sortype = "local_forward";
114: else if (sym & SOR_LOCAL_BACKWARD_SWEEP) sortype = "local_backward";
115: else sortype = "unknown";
116: PetscViewerASCIIPrintf(viewer," type = %s, iterations = %D, local iterations = %D, omega = %g\n",sortype,jac->its,jac->lits,(double)jac->omega);
117: }
118: return(0);
119: }
121: /* ------------------------------------------------------------------------------*/
122: static PetscErrorCode PCSORSetSymmetric_SOR(PC pc,MatSORType flag)
123: {
124: PC_SOR *jac = (PC_SOR*)pc->data;
127: jac->sym = flag;
128: return(0);
129: }
131: static PetscErrorCode PCSORSetOmega_SOR(PC pc,PetscReal omega)
132: {
133: PC_SOR *jac = (PC_SOR*)pc->data;
136: if (omega >= 2.0 || omega <= 0.0) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Relaxation out of range");
137: jac->omega = omega;
138: return(0);
139: }
141: static PetscErrorCode PCSORSetIterations_SOR(PC pc,PetscInt its,PetscInt lits)
142: {
143: PC_SOR *jac = (PC_SOR*)pc->data;
146: jac->its = its;
147: jac->lits = lits;
148: return(0);
149: }
151: static PetscErrorCode PCSORGetSymmetric_SOR(PC pc,MatSORType *flag)
152: {
153: PC_SOR *jac = (PC_SOR*)pc->data;
156: *flag = jac->sym;
157: return(0);
158: }
160: static PetscErrorCode PCSORGetOmega_SOR(PC pc,PetscReal *omega)
161: {
162: PC_SOR *jac = (PC_SOR*)pc->data;
165: *omega = jac->omega;
166: return(0);
167: }
169: static PetscErrorCode PCSORGetIterations_SOR(PC pc,PetscInt *its,PetscInt *lits)
170: {
171: PC_SOR *jac = (PC_SOR*)pc->data;
174: if (its) *its = jac->its;
175: if (lits) *lits = jac->lits;
176: return(0);
177: }
179: /* ------------------------------------------------------------------------------*/
180: /*@
181: PCSORGetSymmetric - Gets the form the SOR preconditioner is using; backward, or forward relaxation. The local variants perform SOR on
182: each processor. By default forward relaxation is used.
184: Logically Collective on PC
186: Input Parameter:
187: . pc - the preconditioner context
189: Output Parameter:
190: . flag - one of the following
191: .vb
192: SOR_FORWARD_SWEEP
193: SOR_BACKWARD_SWEEP
194: SOR_SYMMETRIC_SWEEP
195: SOR_LOCAL_FORWARD_SWEEP
196: SOR_LOCAL_BACKWARD_SWEEP
197: SOR_LOCAL_SYMMETRIC_SWEEP
198: .ve
200: Options Database Keys:
201: + -pc_sor_symmetric - Activates symmetric version
202: . -pc_sor_backward - Activates backward version
203: . -pc_sor_local_forward - Activates local forward version
204: . -pc_sor_local_symmetric - Activates local symmetric version
205: - -pc_sor_local_backward - Activates local backward version
207: Notes:
208: To use the Eisenstat trick with SSOR, employ the PCEISENSTAT preconditioner,
209: which can be chosen with the option
210: . -pc_type eisenstat - Activates Eisenstat trick
212: Level: intermediate
214: .seealso: PCEisenstatSetOmega(), PCSORSetIterations(), PCSORSetOmega(), PCSORSetSymmetric()
215: @*/
216: PetscErrorCode PCSORGetSymmetric(PC pc,MatSORType *flag)
217: {
222: PetscUseMethod(pc,"PCSORGetSymmetric_C",(PC,MatSORType*),(pc,flag));
223: return(0);
224: }
226: /*@
227: PCSORGetOmega - Gets the SOR relaxation coefficient, omega
228: (where omega = 1.0 by default).
230: Logically Collective on PC
232: Input Parameter:
233: . pc - the preconditioner context
235: Output Parameter:
236: . omega - relaxation coefficient (0 < omega < 2).
238: Options Database Key:
239: . -pc_sor_omega <omega> - Sets omega
241: Level: intermediate
243: .seealso: PCSORSetSymmetric(), PCSORSetIterations(), PCEisenstatSetOmega(), PCSORSetOmega()
244: @*/
245: PetscErrorCode PCSORGetOmega(PC pc,PetscReal *omega)
246: {
251: PetscUseMethod(pc,"PCSORGetOmega_C",(PC,PetscReal*),(pc,omega));
252: return(0);
253: }
255: /*@
256: PCSORGetIterations - Gets the number of inner iterations to
257: be used by the SOR preconditioner. The default is 1.
259: Logically Collective on PC
261: Input Parameter:
262: . pc - the preconditioner context
264: Output Parameters:
265: + lits - number of local iterations, smoothings over just variables on processor
266: - its - number of parallel iterations to use; each parallel iteration has lits local iterations
268: Options Database Key:
269: + -pc_sor_its <its> - Sets number of iterations
270: - -pc_sor_lits <lits> - Sets number of local iterations
272: Level: intermediate
274: Notes:
275: When run on one processor the number of smoothings is lits*its
277: .seealso: PCSORSetOmega(), PCSORSetSymmetric(), PCSORSetIterations()
278: @*/
279: PetscErrorCode PCSORGetIterations(PC pc,PetscInt *its,PetscInt *lits)
280: {
285: PetscUseMethod(pc,"PCSORGetIterations_C",(PC,PetscInt*,PetscInt*),(pc,its,lits));
286: return(0);
287: }
289: /*@
290: PCSORSetSymmetric - Sets the SOR preconditioner to use symmetric (SSOR),
291: backward, or forward relaxation. The local variants perform SOR on
292: each processor. By default forward relaxation is used.
294: Logically Collective on PC
296: Input Parameters:
297: + pc - the preconditioner context
298: - flag - one of the following
299: .vb
300: SOR_FORWARD_SWEEP
301: SOR_BACKWARD_SWEEP
302: SOR_SYMMETRIC_SWEEP
303: SOR_LOCAL_FORWARD_SWEEP
304: SOR_LOCAL_BACKWARD_SWEEP
305: SOR_LOCAL_SYMMETRIC_SWEEP
306: .ve
308: Options Database Keys:
309: + -pc_sor_symmetric - Activates symmetric version
310: . -pc_sor_backward - Activates backward version
311: . -pc_sor_local_forward - Activates local forward version
312: . -pc_sor_local_symmetric - Activates local symmetric version
313: - -pc_sor_local_backward - Activates local backward version
315: Notes:
316: To use the Eisenstat trick with SSOR, employ the PCEISENSTAT preconditioner,
317: which can be chosen with the option
318: . -pc_type eisenstat - Activates Eisenstat trick
320: Level: intermediate
322: .seealso: PCEisenstatSetOmega(), PCSORSetIterations(), PCSORSetOmega()
323: @*/
324: PetscErrorCode PCSORSetSymmetric(PC pc,MatSORType flag)
325: {
331: PetscTryMethod(pc,"PCSORSetSymmetric_C",(PC,MatSORType),(pc,flag));
332: return(0);
333: }
335: /*@
336: PCSORSetOmega - Sets the SOR relaxation coefficient, omega
337: (where omega = 1.0 by default).
339: Logically Collective on PC
341: Input Parameters:
342: + pc - the preconditioner context
343: - omega - relaxation coefficient (0 < omega < 2).
345: Options Database Key:
346: . -pc_sor_omega <omega> - Sets omega
348: Level: intermediate
350: Note:
351: If omega != 1, you will need to set the MAT_USE_INODES option to PETSC_FALSE on the matrix.
353: .seealso: PCSORSetSymmetric(), PCSORSetIterations(), PCEisenstatSetOmega(), MatSetOption()
354: @*/
355: PetscErrorCode PCSORSetOmega(PC pc,PetscReal omega)
356: {
362: PetscTryMethod(pc,"PCSORSetOmega_C",(PC,PetscReal),(pc,omega));
363: return(0);
364: }
366: /*@
367: PCSORSetIterations - Sets the number of inner iterations to
368: be used by the SOR preconditioner. The default is 1.
370: Logically Collective on PC
372: Input Parameters:
373: + pc - the preconditioner context
374: . lits - number of local iterations, smoothings over just variables on processor
375: - its - number of parallel iterations to use; each parallel iteration has lits local iterations
377: Options Database Key:
378: + -pc_sor_its <its> - Sets number of iterations
379: - -pc_sor_lits <lits> - Sets number of local iterations
381: Level: intermediate
383: Notes:
384: When run on one processor the number of smoothings is lits*its
386: .seealso: PCSORSetOmega(), PCSORSetSymmetric()
387: @*/
388: PetscErrorCode PCSORSetIterations(PC pc,PetscInt its,PetscInt lits)
389: {
395: PetscTryMethod(pc,"PCSORSetIterations_C",(PC,PetscInt,PetscInt),(pc,its,lits));
396: return(0);
397: }
399: /*MC
400: PCSOR - (S)SOR (successive over relaxation, Gauss-Seidel) preconditioning
402: Options Database Keys:
403: + -pc_sor_symmetric - Activates symmetric version
404: . -pc_sor_backward - Activates backward version
405: . -pc_sor_forward - Activates forward version
406: . -pc_sor_local_forward - Activates local forward version
407: . -pc_sor_local_symmetric - Activates local symmetric version (default version)
408: . -pc_sor_local_backward - Activates local backward version
409: . -pc_sor_omega <omega> - Sets omega
410: . -pc_sor_diagonal_shift <shift> - shift the diagonal entries; useful if the matrix has zeros on the diagonal
411: . -pc_sor_its <its> - Sets number of iterations (default 1)
412: - -pc_sor_lits <lits> - Sets number of local iterations (default 1)
414: Level: beginner
416: Notes:
417: Only implemented for the AIJ and SeqBAIJ matrix formats.
418: Not a true parallel SOR, in parallel this implementation corresponds to block
419: Jacobi with SOR on each block.
421: For AIJ matrix if a diagonal entry is zero (and the diagonal shift is zero) then by default the inverse of that
422: zero will be used and hence the KSPSolve() will terminate with KSP_DIVERGED_NANORIF. If the option
423: KSPSetErrorIfNotConverged() or -ksp_error_if_not_converged the code will terminate as soon as it detects the
424: zero pivot.
426: For SeqBAIJ matrices this implements point-block SOR, but the omega, its, lits options are not supported.
428: For SeqBAIJ the diagonal blocks are inverted using dense LU with partial pivoting. If a zero pivot is detected
429: the computation is stopped with an error
431: If used with KSPRICHARDSON and no monitors the convergence test is skipped to improve speed, thus it always iterates
432: the maximum number of iterations you've selected for KSP. It is usually used in this mode as a smoother for multigrid.
434: If omega != 1, you will need to set the MAT_USE_INODES option to PETSC_FALSE on the matrix.
436: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
437: PCSORSetIterations(), PCSORSetSymmetric(), PCSORSetOmega(), PCEISENSTAT, MatSetOption()
438: M*/
440: PETSC_EXTERN PetscErrorCode PCCreate_SOR(PC pc)
441: {
443: PC_SOR *jac;
446: PetscNewLog(pc,&jac);
448: pc->ops->apply = PCApply_SOR;
449: pc->ops->applytranspose = PCApplyTranspose_SOR;
450: pc->ops->applyrichardson = PCApplyRichardson_SOR;
451: pc->ops->setfromoptions = PCSetFromOptions_SOR;
452: pc->ops->setup = NULL;
453: pc->ops->view = PCView_SOR;
454: pc->ops->destroy = PCDestroy_SOR;
455: pc->data = (void*)jac;
456: jac->sym = SOR_LOCAL_SYMMETRIC_SWEEP;
457: jac->omega = 1.0;
458: jac->fshift = 0.0;
459: jac->its = 1;
460: jac->lits = 1;
462: PetscObjectComposeFunction((PetscObject)pc,"PCSORSetSymmetric_C",PCSORSetSymmetric_SOR);
463: PetscObjectComposeFunction((PetscObject)pc,"PCSORSetOmega_C",PCSORSetOmega_SOR);
464: PetscObjectComposeFunction((PetscObject)pc,"PCSORSetIterations_C",PCSORSetIterations_SOR);
465: PetscObjectComposeFunction((PetscObject)pc,"PCSORGetSymmetric_C",PCSORGetSymmetric_SOR);
466: PetscObjectComposeFunction((PetscObject)pc,"PCSORGetOmega_C",PCSORGetOmega_SOR);
467: PetscObjectComposeFunction((PetscObject)pc,"PCSORGetIterations_C",PCSORGetIterations_SOR);
468: return(0);
469: }