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: {
16: PetscObjectComposeFunction((PetscObject)pc, "PCSORSetSymmetric_C", NULL);
17: PetscObjectComposeFunction((PetscObject)pc, "PCSORSetOmega_C", NULL);
18: PetscObjectComposeFunction((PetscObject)pc, "PCSORSetIterations_C", NULL);
19: PetscObjectComposeFunction((PetscObject)pc, "PCSORGetSymmetric_C", NULL);
20: PetscObjectComposeFunction((PetscObject)pc, "PCSORGetOmega_C", NULL);
21: PetscObjectComposeFunction((PetscObject)pc, "PCSORGetIterations_C", NULL);
22: PetscFree(pc->data);
23: return 0;
24: }
26: static PetscErrorCode PCApply_SOR(PC pc, Vec x, Vec y)
27: {
28: PC_SOR *jac = (PC_SOR *)pc->data;
29: 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;
39: PetscInt flag = jac->sym | SOR_ZERO_INITIAL_GUESS;
40: PetscBool set, sym;
42: MatIsSymmetricKnown(pc->pmat, &set, &sym);
44: MatSOR(pc->pmat, x, jac->omega, (MatSORType)flag, jac->fshift, jac->its, jac->lits, y);
45: MatFactorGetError(pc->pmat, (MatFactorError *)&pc->failedreason);
46: return 0;
47: }
49: 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)
50: {
51: PC_SOR *jac = (PC_SOR *)pc->data;
52: MatSORType stype = jac->sym;
54: PetscInfo(pc, "Warning, convergence criteria ignored, using %" PetscInt_FMT " iterations\n", its);
55: if (guesszero) stype = (MatSORType)(stype | SOR_ZERO_INITIAL_GUESS);
56: MatSOR(pc->pmat, b, jac->omega, stype, jac->fshift, its * jac->its, jac->lits, y);
57: MatFactorGetError(pc->pmat, (MatFactorError *)&pc->failedreason);
58: *outits = its;
59: *reason = PCRICHARDSON_CONVERGED_ITS;
60: return 0;
61: }
63: PetscErrorCode PCSetFromOptions_SOR(PC pc, PetscOptionItems *PetscOptionsObject)
64: {
65: PC_SOR *jac = (PC_SOR *)pc->data;
66: PetscBool flg;
68: PetscOptionsHeadBegin(PetscOptionsObject, "(S)SOR options");
69: PetscOptionsReal("-pc_sor_omega", "relaxation factor (0 < omega < 2)", "PCSORSetOmega", jac->omega, &jac->omega, NULL);
70: PetscOptionsReal("-pc_sor_diagonal_shift", "Add to the diagonal entries", "", jac->fshift, &jac->fshift, NULL);
71: PetscOptionsInt("-pc_sor_its", "number of inner SOR iterations", "PCSORSetIterations", jac->its, &jac->its, NULL);
72: PetscOptionsInt("-pc_sor_lits", "number of local inner SOR iterations", "PCSORSetIterations", jac->lits, &jac->lits, NULL);
73: PetscOptionsBoolGroupBegin("-pc_sor_symmetric", "SSOR, not SOR", "PCSORSetSymmetric", &flg);
74: if (flg) PCSORSetSymmetric(pc, SOR_SYMMETRIC_SWEEP);
75: PetscOptionsBoolGroup("-pc_sor_backward", "use backward sweep instead of forward", "PCSORSetSymmetric", &flg);
76: if (flg) PCSORSetSymmetric(pc, SOR_BACKWARD_SWEEP);
77: PetscOptionsBoolGroup("-pc_sor_forward", "use forward sweep", "PCSORSetSymmetric", &flg);
78: if (flg) PCSORSetSymmetric(pc, SOR_FORWARD_SWEEP);
79: PetscOptionsBoolGroup("-pc_sor_local_symmetric", "use SSOR separately on each processor", "PCSORSetSymmetric", &flg);
80: if (flg) PCSORSetSymmetric(pc, SOR_LOCAL_SYMMETRIC_SWEEP);
81: PetscOptionsBoolGroup("-pc_sor_local_backward", "use backward sweep locally", "PCSORSetSymmetric", &flg);
82: if (flg) PCSORSetSymmetric(pc, SOR_LOCAL_BACKWARD_SWEEP);
83: PetscOptionsBoolGroupEnd("-pc_sor_local_forward", "use forward sweep locally", "PCSORSetSymmetric", &flg);
84: if (flg) PCSORSetSymmetric(pc, SOR_LOCAL_FORWARD_SWEEP);
85: PetscOptionsHeadEnd();
86: return 0;
87: }
89: PetscErrorCode PCView_SOR(PC pc, PetscViewer viewer)
90: {
91: PC_SOR *jac = (PC_SOR *)pc->data;
92: MatSORType sym = jac->sym;
93: const char *sortype;
94: PetscBool iascii;
96: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
97: if (iascii) {
98: if (sym & SOR_ZERO_INITIAL_GUESS) PetscViewerASCIIPrintf(viewer, " zero initial guess\n");
99: if (sym == SOR_APPLY_UPPER) sortype = "apply_upper";
100: else if (sym == SOR_APPLY_LOWER) sortype = "apply_lower";
101: else if (sym & SOR_EISENSTAT) sortype = "Eisenstat";
102: else if ((sym & SOR_SYMMETRIC_SWEEP) == SOR_SYMMETRIC_SWEEP) sortype = "symmetric";
103: else if (sym & SOR_BACKWARD_SWEEP) sortype = "backward";
104: else if (sym & SOR_FORWARD_SWEEP) sortype = "forward";
105: else if ((sym & SOR_LOCAL_SYMMETRIC_SWEEP) == SOR_LOCAL_SYMMETRIC_SWEEP) sortype = "local_symmetric";
106: else if (sym & SOR_LOCAL_FORWARD_SWEEP) sortype = "local_forward";
107: else if (sym & SOR_LOCAL_BACKWARD_SWEEP) sortype = "local_backward";
108: else sortype = "unknown";
109: PetscViewerASCIIPrintf(viewer, " type = %s, iterations = %" PetscInt_FMT ", local iterations = %" PetscInt_FMT ", omega = %g\n", sortype, jac->its, jac->lits, (double)jac->omega);
110: }
111: return 0;
112: }
114: static PetscErrorCode PCSORSetSymmetric_SOR(PC pc, MatSORType flag)
115: {
116: PC_SOR *jac = (PC_SOR *)pc->data;
118: jac->sym = flag;
119: return 0;
120: }
122: static PetscErrorCode PCSORSetOmega_SOR(PC pc, PetscReal omega)
123: {
124: PC_SOR *jac = (PC_SOR *)pc->data;
127: jac->omega = omega;
128: return 0;
129: }
131: static PetscErrorCode PCSORSetIterations_SOR(PC pc, PetscInt its, PetscInt lits)
132: {
133: PC_SOR *jac = (PC_SOR *)pc->data;
135: jac->its = its;
136: jac->lits = lits;
137: return 0;
138: }
140: static PetscErrorCode PCSORGetSymmetric_SOR(PC pc, MatSORType *flag)
141: {
142: PC_SOR *jac = (PC_SOR *)pc->data;
144: *flag = jac->sym;
145: return 0;
146: }
148: static PetscErrorCode PCSORGetOmega_SOR(PC pc, PetscReal *omega)
149: {
150: PC_SOR *jac = (PC_SOR *)pc->data;
152: *omega = jac->omega;
153: return 0;
154: }
156: static PetscErrorCode PCSORGetIterations_SOR(PC pc, PetscInt *its, PetscInt *lits)
157: {
158: PC_SOR *jac = (PC_SOR *)pc->data;
160: if (its) *its = jac->its;
161: if (lits) *lits = jac->lits;
162: return 0;
163: }
165: /*@
166: PCSORGetSymmetric - Gets the form the SOR preconditioner is using; backward, or forward relaxation. The local variants perform SOR on
167: each processor. By default forward relaxation is used.
169: Logically Collective
171: Input Parameter:
172: . pc - the preconditioner context
174: Output Parameter:
175: . flag - one of the following
176: .vb
177: SOR_FORWARD_SWEEP
178: SOR_BACKWARD_SWEEP
179: SOR_SYMMETRIC_SWEEP
180: SOR_LOCAL_FORWARD_SWEEP
181: SOR_LOCAL_BACKWARD_SWEEP
182: SOR_LOCAL_SYMMETRIC_SWEEP
183: .ve
185: Options Database Keys:
186: + -pc_sor_symmetric - Activates symmetric version
187: . -pc_sor_backward - Activates backward version
188: . -pc_sor_local_forward - Activates local forward version
189: . -pc_sor_local_symmetric - Activates local symmetric version
190: - -pc_sor_local_backward - Activates local backward version
192: Note:
193: To use the Eisenstat trick with SSOR, employ the `PCEISENSTAT` preconditioner,
194: which can be chosen with the option
195: . -pc_type eisenstat - Activates Eisenstat trick
197: Level: intermediate
199: .seealso: `PCSOR`, `PCEisenstatSetOmega()`, `PCSORSetIterations()`, `PCSORSetOmega()`, `PCSORSetSymmetric()`
200: @*/
201: PetscErrorCode PCSORGetSymmetric(PC pc, MatSORType *flag)
202: {
204: PetscUseMethod(pc, "PCSORGetSymmetric_C", (PC, MatSORType *), (pc, flag));
205: return 0;
206: }
208: /*@
209: PCSORGetOmega - Gets the SOR relaxation coefficient, omega
210: (where omega = 1.0 by default).
212: Logically Collective
214: Input Parameter:
215: . pc - the preconditioner context
217: Output Parameter:
218: . omega - relaxation coefficient (0 < omega < 2).
220: Options Database Key:
221: . -pc_sor_omega <omega> - Sets omega
223: Level: intermediate
225: .seealso: `PCSOR`, `PCSORSetSymmetric()`, `PCSORSetIterations()`, `PCEisenstatSetOmega()`, `PCSORSetOmega()`
226: @*/
227: PetscErrorCode PCSORGetOmega(PC pc, PetscReal *omega)
228: {
230: PetscUseMethod(pc, "PCSORGetOmega_C", (PC, PetscReal *), (pc, omega));
231: return 0;
232: }
234: /*@
235: PCSORGetIterations - Gets the number of inner iterations to
236: be used by the SOR preconditioner. The default is 1.
238: Logically Collective
240: Input Parameter:
241: . pc - the preconditioner context
243: Output Parameters:
244: + lits - number of local iterations, smoothings over just variables on processor
245: - its - number of parallel iterations to use; each parallel iteration has lits local iterations
247: Options Database Keys:
248: + -pc_sor_its <its> - Sets number of iterations
249: - -pc_sor_lits <lits> - Sets number of local iterations
251: Level: intermediate
253: Note:
254: When run on one processor the number of smoothings is lits*its
256: .seealso: `PCSOR`, `PCSORSetOmega()`, `PCSORSetSymmetric()`, `PCSORSetIterations()`
257: @*/
258: PetscErrorCode PCSORGetIterations(PC pc, PetscInt *its, PetscInt *lits)
259: {
261: PetscUseMethod(pc, "PCSORGetIterations_C", (PC, PetscInt *, PetscInt *), (pc, its, lits));
262: return 0;
263: }
265: /*@
266: PCSORSetSymmetric - Sets the SOR preconditioner to use symmetric (SSOR),
267: backward, or forward relaxation. The local variants perform SOR on
268: each processor. By default forward relaxation is used.
270: Logically Collective
272: Input Parameters:
273: + pc - the preconditioner context
274: - flag - one of the following
275: .vb
276: SOR_FORWARD_SWEEP
277: SOR_BACKWARD_SWEEP
278: SOR_SYMMETRIC_SWEEP
279: SOR_LOCAL_FORWARD_SWEEP
280: SOR_LOCAL_BACKWARD_SWEEP
281: SOR_LOCAL_SYMMETRIC_SWEEP
282: .ve
284: Options Database Keys:
285: + -pc_sor_symmetric - Activates symmetric version
286: . -pc_sor_backward - Activates backward version
287: . -pc_sor_local_forward - Activates local forward version
288: . -pc_sor_local_symmetric - Activates local symmetric version
289: - -pc_sor_local_backward - Activates local backward version
291: Note:
292: To use the Eisenstat trick with SSOR, employ the PCEISENSTAT preconditioner,
293: which can be chosen with the option
294: . -pc_type eisenstat - Activates Eisenstat trick
296: Level: intermediate
298: .seealso: `PCSOR`, `PCEisenstatSetOmega()`, `PCSORSetIterations()`, `PCSORSetOmega()`
299: @*/
300: PetscErrorCode PCSORSetSymmetric(PC pc, MatSORType flag)
301: {
304: PetscTryMethod(pc, "PCSORSetSymmetric_C", (PC, MatSORType), (pc, flag));
305: return 0;
306: }
308: /*@
309: PCSORSetOmega - Sets the SOR relaxation coefficient, omega
310: (where omega = 1.0 by default).
312: Logically Collective
314: Input Parameters:
315: + pc - the preconditioner context
316: - omega - relaxation coefficient (0 < omega < 2).
318: Options Database Key:
319: . -pc_sor_omega <omega> - Sets omega
321: Level: intermediate
323: Note:
324: If omega != 1, you will need to set the `MAT_USE_INODE`S option to `PETSC_FALSE` on the matrix.
326: .seealso: `PCSOR`, `PCSORSetSymmetric()`, `PCSORSetIterations()`, `PCEisenstatSetOmega()`, `MatSetOption()`
327: @*/
328: PetscErrorCode PCSORSetOmega(PC pc, PetscReal omega)
329: {
332: PetscTryMethod(pc, "PCSORSetOmega_C", (PC, PetscReal), (pc, omega));
333: return 0;
334: }
336: /*@
337: PCSORSetIterations - Sets the number of inner iterations to
338: be used by the SOR preconditioner. The default is 1.
340: Logically Collective
342: Input Parameters:
343: + pc - the preconditioner context
344: . lits - number of local iterations, smoothings over just variables on processor
345: - its - number of parallel iterations to use; each parallel iteration has lits local iterations
347: Options Database Keys:
348: + -pc_sor_its <its> - Sets number of iterations
349: - -pc_sor_lits <lits> - Sets number of local iterations
351: Level: intermediate
353: Note:
354: When run on one processor the number of smoothings is lits*its
356: .seealso: `PCSOR`, `PCSORSetOmega()`, `PCSORSetSymmetric()`
357: @*/
358: PetscErrorCode PCSORSetIterations(PC pc, PetscInt its, PetscInt lits)
359: {
362: PetscTryMethod(pc, "PCSORSetIterations_C", (PC, PetscInt, PetscInt), (pc, its, lits));
363: return 0;
364: }
366: /*MC
367: PCSOR - (S)SOR (successive over relaxation, Gauss-Seidel) preconditioning
369: Options Database Keys:
370: + -pc_sor_symmetric - Activates symmetric version
371: . -pc_sor_backward - Activates backward version
372: . -pc_sor_forward - Activates forward version
373: . -pc_sor_local_forward - Activates local forward version
374: . -pc_sor_local_symmetric - Activates local symmetric version (default version)
375: . -pc_sor_local_backward - Activates local backward version
376: . -pc_sor_omega <omega> - Sets omega
377: . -pc_sor_diagonal_shift <shift> - shift the diagonal entries; useful if the matrix has zeros on the diagonal
378: . -pc_sor_its <its> - Sets number of iterations (default 1)
379: - -pc_sor_lits <lits> - Sets number of local iterations (default 1)
381: Level: beginner
383: Notes:
384: Only implemented for the `MATAIJ` and `MATSEQBAIJ` matrix formats.
386: Not a true parallel SOR, in parallel this implementation corresponds to block
387: Jacobi with SOR on each block.
389: For `MATAIJ` matrices if a diagonal entry is zero (and the diagonal shift is zero) then by default the inverse of that
390: zero will be used and hence the `KSPSolve()` will terminate with `KSP_DIVERGED_NANORINF`. If the option
391: `KSPSetErrorIfNotConverged()` or -ksp_error_if_not_converged the code will terminate as soon as it detects the
392: zero pivot.
394: For `MATSEQBAIJ` matrices this implements point-block SOR, but the omega, its, lits options are not supported.
396: For `MATSEQBAIJ` the diagonal blocks are inverted using dense LU with partial pivoting. If a zero pivot is detected
397: the computation is stopped with an error
399: If used with `KSPRICHARDSON` and no monitors the convergence test is skipped to improve speed, thus it always iterates
400: the maximum number of iterations you've selected for `KSP`. It is usually used in this mode as a smoother for multigrid.
402: If omega != 1, you will need to set the `MAT_USE_INODES` option to `PETSC_FALSE` on the matrix.
404: .seealso: `PCCreate()`, `PCSetType()`, `PCType`, `PC`, `PCJACOBI`,
405: `PCSORSetIterations()`, `PCSORSetSymmetric()`, `PCSORSetOmega()`, `PCEISENSTAT`, `MatSetOption()`
406: M*/
408: PETSC_EXTERN PetscErrorCode PCCreate_SOR(PC pc)
409: {
410: PC_SOR *jac;
412: PetscNew(&jac);
414: pc->ops->apply = PCApply_SOR;
415: pc->ops->applytranspose = PCApplyTranspose_SOR;
416: pc->ops->applyrichardson = PCApplyRichardson_SOR;
417: pc->ops->setfromoptions = PCSetFromOptions_SOR;
418: pc->ops->setup = NULL;
419: pc->ops->view = PCView_SOR;
420: pc->ops->destroy = PCDestroy_SOR;
421: pc->data = (void *)jac;
422: jac->sym = SOR_LOCAL_SYMMETRIC_SWEEP;
423: jac->omega = 1.0;
424: jac->fshift = 0.0;
425: jac->its = 1;
426: jac->lits = 1;
428: PetscObjectComposeFunction((PetscObject)pc, "PCSORSetSymmetric_C", PCSORSetSymmetric_SOR);
429: PetscObjectComposeFunction((PetscObject)pc, "PCSORSetOmega_C", PCSORSetOmega_SOR);
430: PetscObjectComposeFunction((PetscObject)pc, "PCSORSetIterations_C", PCSORSetIterations_SOR);
431: PetscObjectComposeFunction((PetscObject)pc, "PCSORGetSymmetric_C", PCSORGetSymmetric_SOR);
432: PetscObjectComposeFunction((PetscObject)pc, "PCSORGetOmega_C", PCSORGetOmega_SOR);
433: PetscObjectComposeFunction((PetscObject)pc, "PCSORGetIterations_C", PCSORGetIterations_SOR);
434: return 0;
435: }