Actual source code: sor.c
petsc-3.11.4 2019-09-28
2: /*
3: Defines a (S)SOR preconditioner for any Mat implementation
4: */
5: #include <petsc/private/pcimpl.h>
7: typedef struct {
8: PetscInt its; /* inner iterations, number of sweeps */
9: PetscInt lits; /* local inner iterations, number of sweeps applied by the local matrix mat->A */
10: MatSORType sym; /* forward, reverse, symmetric etc. */
11: PetscReal omega;
12: PetscReal fshift;
13: } PC_SOR;
15: static PetscErrorCode PCDestroy_SOR(PC pc)
16: {
20: PetscFree(pc->data);
21: return(0);
22: }
24: static PetscErrorCode PCApply_SOR(PC pc,Vec x,Vec y)
25: {
26: PC_SOR *jac = (PC_SOR*)pc->data;
28: PetscInt flag = jac->sym | SOR_ZERO_INITIAL_GUESS;
31: MatSOR(pc->pmat,x,jac->omega,(MatSORType)flag,jac->fshift,jac->its,jac->lits,y);
32: MatFactorGetError(pc->pmat,(MatFactorError*)&pc->failedreason);
33: return(0);
34: }
36: static PetscErrorCode PCApplyTranspose_SOR(PC pc,Vec x,Vec y)
37: {
38: PC_SOR *jac = (PC_SOR*)pc->data;
40: PetscInt flag = jac->sym | SOR_ZERO_INITIAL_GUESS;
41: PetscBool set,sym;
44: MatIsSymmetricKnown(pc->pmat,&set,&sym);
45: 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");
46: MatSOR(pc->pmat,x,jac->omega,(MatSORType)flag,jac->fshift,jac->its,jac->lits,y);
47: MatFactorGetError(pc->pmat,(MatFactorError*)&pc->failedreason);
48: return(0);
49: }
51: 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)
52: {
53: PC_SOR *jac = (PC_SOR*)pc->data;
55: MatSORType stype = jac->sym;
58: PetscInfo1(pc,"Warning, convergence critera ignored, using %D iterations\n",its);
59: if (guesszero) stype = (MatSORType) (stype | SOR_ZERO_INITIAL_GUESS);
60: MatSOR(pc->pmat,b,jac->omega,stype,jac->fshift,its*jac->its,jac->lits,y);
61: MatFactorGetError(pc->pmat,(MatFactorError*)&pc->failedreason);
62: *outits = its;
63: *reason = PCRICHARDSON_CONVERGED_ITS;
64: return(0);
65: }
67: PetscErrorCode PCSetFromOptions_SOR(PetscOptionItems *PetscOptionsObject,PC pc)
68: {
69: PC_SOR *jac = (PC_SOR*)pc->data;
71: PetscBool flg;
74: PetscOptionsHead(PetscOptionsObject,"(S)SOR options");
75: PetscOptionsReal("-pc_sor_omega","relaxation factor (0 < omega < 2)","PCSORSetOmega",jac->omega,&jac->omega,NULL);
76: PetscOptionsReal("-pc_sor_diagonal_shift","Add to the diagonal entries","",jac->fshift,&jac->fshift,NULL);
77: PetscOptionsInt("-pc_sor_its","number of inner SOR iterations","PCSORSetIterations",jac->its,&jac->its,NULL);
78: PetscOptionsInt("-pc_sor_lits","number of local inner SOR iterations","PCSORSetIterations",jac->lits,&jac->lits,NULL);
79: PetscOptionsBoolGroupBegin("-pc_sor_symmetric","SSOR, not SOR","PCSORSetSymmetric",&flg);
80: if (flg) {PCSORSetSymmetric(pc,SOR_SYMMETRIC_SWEEP);}
81: PetscOptionsBoolGroup("-pc_sor_backward","use backward sweep instead of forward","PCSORSetSymmetric",&flg);
82: if (flg) {PCSORSetSymmetric(pc,SOR_BACKWARD_SWEEP);}
83: PetscOptionsBoolGroup("-pc_sor_forward","use forward sweep","PCSORSetSymmetric",&flg);
84: if (flg) {PCSORSetSymmetric(pc,SOR_FORWARD_SWEEP);}
85: PetscOptionsBoolGroup("-pc_sor_local_symmetric","use SSOR separately on each processor","PCSORSetSymmetric",&flg);
86: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_SYMMETRIC_SWEEP);}
87: PetscOptionsBoolGroup("-pc_sor_local_backward","use backward sweep locally","PCSORSetSymmetric",&flg);
88: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_BACKWARD_SWEEP);}
89: PetscOptionsBoolGroupEnd("-pc_sor_local_forward","use forward sweep locally","PCSORSetSymmetric",&flg);
90: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_FORWARD_SWEEP);}
91: PetscOptionsTail();
92: return(0);
93: }
95: PetscErrorCode PCView_SOR(PC pc,PetscViewer viewer)
96: {
97: PC_SOR *jac = (PC_SOR*)pc->data;
98: MatSORType sym = jac->sym;
99: const char *sortype;
101: PetscBool iascii;
104: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
105: if (iascii) {
106: if (sym & SOR_ZERO_INITIAL_GUESS) {PetscViewerASCIIPrintf(viewer," zero initial guess\n");}
107: if (sym == SOR_APPLY_UPPER) sortype = "apply_upper";
108: else if (sym == SOR_APPLY_LOWER) sortype = "apply_lower";
109: else if (sym & SOR_EISENSTAT) sortype = "Eisenstat";
110: else if ((sym & SOR_SYMMETRIC_SWEEP) == SOR_SYMMETRIC_SWEEP) sortype = "symmetric";
111: else if (sym & SOR_BACKWARD_SWEEP) sortype = "backward";
112: else if (sym & SOR_FORWARD_SWEEP) sortype = "forward";
113: else if ((sym & SOR_LOCAL_SYMMETRIC_SWEEP) == SOR_LOCAL_SYMMETRIC_SWEEP) sortype = "local_symmetric";
114: else if (sym & SOR_LOCAL_FORWARD_SWEEP) sortype = "local_forward";
115: else if (sym & SOR_LOCAL_BACKWARD_SWEEP) sortype = "local_backward";
116: else sortype = "unknown";
117: PetscViewerASCIIPrintf(viewer," type = %s, iterations = %D, local iterations = %D, omega = %g\n",sortype,jac->its,jac->lits,(double)jac->omega);
118: }
119: return(0);
120: }
123: /* ------------------------------------------------------------------------------*/
124: static PetscErrorCode PCSORSetSymmetric_SOR(PC pc,MatSORType flag)
125: {
126: PC_SOR *jac = (PC_SOR*)pc->data;
129: jac->sym = flag;
130: return(0);
131: }
133: static PetscErrorCode PCSORSetOmega_SOR(PC pc,PetscReal omega)
134: {
135: PC_SOR *jac = (PC_SOR*)pc->data;
138: if (omega >= 2.0 || omega <= 0.0) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Relaxation out of range");
139: jac->omega = omega;
140: return(0);
141: }
143: static PetscErrorCode PCSORSetIterations_SOR(PC pc,PetscInt its,PetscInt lits)
144: {
145: PC_SOR *jac = (PC_SOR*)pc->data;
148: jac->its = its;
149: jac->lits = lits;
150: return(0);
151: }
153: static PetscErrorCode PCSORGetSymmetric_SOR(PC pc,MatSORType *flag)
154: {
155: PC_SOR *jac = (PC_SOR*)pc->data;
158: *flag = jac->sym;
159: return(0);
160: }
162: static PetscErrorCode PCSORGetOmega_SOR(PC pc,PetscReal *omega)
163: {
164: PC_SOR *jac = (PC_SOR*)pc->data;
167: *omega = jac->omega;
168: return(0);
169: }
171: static PetscErrorCode PCSORGetIterations_SOR(PC pc,PetscInt *its,PetscInt *lits)
172: {
173: PC_SOR *jac = (PC_SOR*)pc->data;
176: if (its) *its = jac->its;
177: if (lits) *lits = jac->lits;
178: return(0);
179: }
181: /* ------------------------------------------------------------------------------*/
182: /*@
183: PCSORGetSymmetric - Gets the form the SOR preconditioner is using; backward, or forward relaxation. The local variants perform SOR on
184: each processor. By default forward relaxation is used.
186: Logically Collective on PC
188: Input Parameter:
189: . pc - the preconditioner context
191: Output Parameter:
192: . flag - one of the following
193: .vb
194: SOR_FORWARD_SWEEP
195: SOR_BACKWARD_SWEEP
196: SOR_SYMMETRIC_SWEEP
197: SOR_LOCAL_FORWARD_SWEEP
198: SOR_LOCAL_BACKWARD_SWEEP
199: SOR_LOCAL_SYMMETRIC_SWEEP
200: .ve
202: Options Database Keys:
203: + -pc_sor_symmetric - Activates symmetric version
204: . -pc_sor_backward - Activates backward version
205: . -pc_sor_local_forward - Activates local forward version
206: . -pc_sor_local_symmetric - Activates local symmetric version
207: - -pc_sor_local_backward - Activates local backward version
209: Notes:
210: To use the Eisenstat trick with SSOR, employ the PCEISENSTAT preconditioner,
211: which can be chosen with the option
212: . -pc_type eisenstat - Activates Eisenstat trick
214: Level: intermediate
216: .keywords: PC, SOR, SSOR, set, relaxation, sweep, forward, backward, symmetric
218: .seealso: PCEisenstatSetOmega(), PCSORSetIterations(), PCSORSetOmega(), PCSORSetSymmetric()
219: @*/
220: PetscErrorCode PCSORGetSymmetric(PC pc,MatSORType *flag)
221: {
226: PetscUseMethod(pc,"PCSORGetSymmetric_C",(PC,MatSORType*),(pc,flag));
227: return(0);
228: }
230: /*@
231: PCSORGetOmega - Gets the SOR relaxation coefficient, omega
232: (where omega = 1.0 by default).
234: Logically Collective on PC
236: Input Parameter:
237: . pc - the preconditioner context
239: Output Parameter:
240: . omega - relaxation coefficient (0 < omega < 2).
242: Options Database Key:
243: . -pc_sor_omega <omega> - Sets omega
245: Level: intermediate
247: .keywords: PC, SOR, SSOR, set, relaxation, omega
249: .seealso: PCSORSetSymmetric(), PCSORSetIterations(), PCEisenstatSetOmega(), PCSORSetOmega()
250: @*/
251: PetscErrorCode PCSORGetOmega(PC pc,PetscReal *omega)
252: {
257: PetscUseMethod(pc,"PCSORGetOmega_C",(PC,PetscReal*),(pc,omega));
258: return(0);
259: }
261: /*@
262: PCSORGetIterations - Gets the number of inner iterations to
263: be used by the SOR preconditioner. The default is 1.
265: Logically Collective on PC
267: Input Parameter:
268: . pc - the preconditioner context
270: Output Parameter:
271: + lits - number of local iterations, smoothings over just variables on processor
272: - its - number of parallel iterations to use; each parallel iteration has lits local iterations
274: Options Database Key:
275: + -pc_sor_its <its> - Sets number of iterations
276: - -pc_sor_lits <lits> - Sets number of local iterations
278: Level: intermediate
280: Notes:
281: When run on one processor the number of smoothings is lits*its
283: .keywords: PC, SOR, SSOR, set, iterations
285: .seealso: PCSORSetOmega(), PCSORSetSymmetric(), PCSORSetIterations()
286: @*/
287: PetscErrorCode PCSORGetIterations(PC pc,PetscInt *its,PetscInt *lits)
288: {
293: PetscUseMethod(pc,"PCSORGetIterations_C",(PC,PetscInt*,PetscInt*),(pc,its,lits));
294: return(0);
295: }
297: /*@
298: PCSORSetSymmetric - Sets the SOR preconditioner to use symmetric (SSOR),
299: backward, or forward relaxation. The local variants perform SOR on
300: each processor. By default forward relaxation is used.
302: Logically Collective on PC
304: Input Parameters:
305: + pc - the preconditioner context
306: - flag - one of the following
307: .vb
308: SOR_FORWARD_SWEEP
309: SOR_BACKWARD_SWEEP
310: SOR_SYMMETRIC_SWEEP
311: SOR_LOCAL_FORWARD_SWEEP
312: SOR_LOCAL_BACKWARD_SWEEP
313: SOR_LOCAL_SYMMETRIC_SWEEP
314: .ve
316: Options Database Keys:
317: + -pc_sor_symmetric - Activates symmetric version
318: . -pc_sor_backward - Activates backward version
319: . -pc_sor_local_forward - Activates local forward version
320: . -pc_sor_local_symmetric - Activates local symmetric version
321: - -pc_sor_local_backward - Activates local backward version
323: Notes:
324: To use the Eisenstat trick with SSOR, employ the PCEISENSTAT preconditioner,
325: which can be chosen with the option
326: . -pc_type eisenstat - Activates Eisenstat trick
328: Level: intermediate
330: .keywords: PC, SOR, SSOR, set, relaxation, sweep, forward, backward, symmetric
332: .seealso: PCEisenstatSetOmega(), PCSORSetIterations(), PCSORSetOmega()
333: @*/
334: PetscErrorCode PCSORSetSymmetric(PC pc,MatSORType flag)
335: {
341: PetscTryMethod(pc,"PCSORSetSymmetric_C",(PC,MatSORType),(pc,flag));
342: return(0);
343: }
345: /*@
346: PCSORSetOmega - Sets the SOR relaxation coefficient, omega
347: (where omega = 1.0 by default).
349: Logically Collective on PC
351: Input Parameters:
352: + pc - the preconditioner context
353: - omega - relaxation coefficient (0 < omega < 2).
355: Options Database Key:
356: . -pc_sor_omega <omega> - Sets omega
358: Level: intermediate
360: .keywords: PC, SOR, SSOR, set, relaxation, omega
362: .seealso: PCSORSetSymmetric(), PCSORSetIterations(), PCEisenstatSetOmega()
363: @*/
364: PetscErrorCode PCSORSetOmega(PC pc,PetscReal omega)
365: {
371: PetscTryMethod(pc,"PCSORSetOmega_C",(PC,PetscReal),(pc,omega));
372: return(0);
373: }
375: /*@
376: PCSORSetIterations - Sets the number of inner iterations to
377: be used by the SOR preconditioner. The default is 1.
379: Logically Collective on PC
381: Input Parameters:
382: + pc - the preconditioner context
383: . lits - number of local iterations, smoothings over just variables on processor
384: - its - number of parallel iterations to use; each parallel iteration has lits local iterations
386: Options Database Key:
387: + -pc_sor_its <its> - Sets number of iterations
388: - -pc_sor_lits <lits> - Sets number of local iterations
390: Level: intermediate
392: Notes:
393: When run on one processor the number of smoothings is lits*its
395: .keywords: PC, SOR, SSOR, set, iterations
397: .seealso: PCSORSetOmega(), PCSORSetSymmetric()
398: @*/
399: PetscErrorCode PCSORSetIterations(PC pc,PetscInt its,PetscInt lits)
400: {
406: PetscTryMethod(pc,"PCSORSetIterations_C",(PC,PetscInt,PetscInt),(pc,its,lits));
407: return(0);
408: }
410: /*MC
411: PCSOR - (S)SOR (successive over relaxation, Gauss-Seidel) preconditioning
413: Options Database Keys:
414: + -pc_sor_symmetric - Activates symmetric version
415: . -pc_sor_backward - Activates backward version
416: . -pc_sor_forward - Activates forward version
417: . -pc_sor_local_forward - Activates local forward version
418: . -pc_sor_local_symmetric - Activates local symmetric version (default version)
419: . -pc_sor_local_backward - Activates local backward version
420: . -pc_sor_omega <omega> - Sets omega
421: . -pc_sor_diagonal_shift <shift> - shift the diagonal entries; useful if the matrix has zeros on the diagonal
422: . -pc_sor_its <its> - Sets number of iterations (default 1)
423: - -pc_sor_lits <lits> - Sets number of local iterations (default 1)
425: Level: beginner
427: Concepts: SOR, preconditioners, Gauss-Seidel
429: Notes:
430: Only implemented for the AIJ and SeqBAIJ matrix formats.
431: Not a true parallel SOR, in parallel this implementation corresponds to block
432: Jacobi with SOR on each block.
434: For AIJ matrix if a diagonal entry is zero (and the diagonal shift is zero) then by default the inverse of that
435: zero will be used and hence the KSPSolve() will terminate with KSP_DIVERGED_NANORIF. If the option
436: KSPSetErrorIfNotConverged() or -ksp_error_if_not_converged the code will terminate as soon as it detects the
437: zero pivot.
439: For SeqBAIJ matrices this implements point-block SOR, but the omega, its, lits options are not supported.
441: For SeqBAIJ the diagonal blocks are inverted using dense LU with partial pivoting. If a zero pivot is detected
442: the computation is stopped with an error
444: If used with KSPRICHARDSON and no monitors the convergence test is skipped to improve speed, thus it always iterates
445: the maximum number of iterations you've selected for KSP. It is usually used in this mode as a smoother for multigrid.
447: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
448: PCSORSetIterations(), PCSORSetSymmetric(), PCSORSetOmega(), PCEISENSTAT
449: M*/
451: PETSC_EXTERN PetscErrorCode PCCreate_SOR(PC pc)
452: {
454: PC_SOR *jac;
457: PetscNewLog(pc,&jac);
459: pc->ops->apply = PCApply_SOR;
460: pc->ops->applytranspose = PCApplyTranspose_SOR;
461: pc->ops->applyrichardson = PCApplyRichardson_SOR;
462: pc->ops->setfromoptions = PCSetFromOptions_SOR;
463: pc->ops->setup = 0;
464: pc->ops->view = PCView_SOR;
465: pc->ops->destroy = PCDestroy_SOR;
466: pc->data = (void*)jac;
467: jac->sym = SOR_LOCAL_SYMMETRIC_SWEEP;
468: jac->omega = 1.0;
469: jac->fshift = 0.0;
470: jac->its = 1;
471: jac->lits = 1;
473: PetscObjectComposeFunction((PetscObject)pc,"PCSORSetSymmetric_C",PCSORSetSymmetric_SOR);
474: PetscObjectComposeFunction((PetscObject)pc,"PCSORSetOmega_C",PCSORSetOmega_SOR);
475: PetscObjectComposeFunction((PetscObject)pc,"PCSORSetIterations_C",PCSORSetIterations_SOR);
476: PetscObjectComposeFunction((PetscObject)pc,"PCSORGetSymmetric_C",PCSORGetSymmetric_SOR);
477: PetscObjectComposeFunction((PetscObject)pc,"PCSORGetOmega_C",PCSORGetOmega_SOR);
478: PetscObjectComposeFunction((PetscObject)pc,"PCSORGetIterations_C",PCSORGetIterations_SOR);
479: return(0);
480: }