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: PetscFree(pc->data);
17: return 0;
18: }
20: static PetscErrorCode PCApply_SOR(PC pc,Vec x,Vec y)
21: {
22: PC_SOR *jac = (PC_SOR*)pc->data;
23: PetscInt flag = jac->sym | SOR_ZERO_INITIAL_GUESS;
25: MatSOR(pc->pmat,x,jac->omega,(MatSORType)flag,jac->fshift,jac->its,jac->lits,y);
26: MatFactorGetError(pc->pmat,(MatFactorError*)&pc->failedreason);
27: return 0;
28: }
30: static PetscErrorCode PCApplyTranspose_SOR(PC pc,Vec x,Vec y)
31: {
32: PC_SOR *jac = (PC_SOR*)pc->data;
33: PetscInt flag = jac->sym | SOR_ZERO_INITIAL_GUESS;
34: PetscBool set,sym;
36: MatIsSymmetricKnown(pc->pmat,&set,&sym);
38: MatSOR(pc->pmat,x,jac->omega,(MatSORType)flag,jac->fshift,jac->its,jac->lits,y);
39: MatFactorGetError(pc->pmat,(MatFactorError*)&pc->failedreason);
40: return 0;
41: }
43: 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)
44: {
45: PC_SOR *jac = (PC_SOR*)pc->data;
46: MatSORType stype = jac->sym;
48: PetscInfo(pc,"Warning, convergence critera ignored, using %D iterations\n",its);
49: if (guesszero) stype = (MatSORType) (stype | SOR_ZERO_INITIAL_GUESS);
50: MatSOR(pc->pmat,b,jac->omega,stype,jac->fshift,its*jac->its,jac->lits,y);
51: MatFactorGetError(pc->pmat,(MatFactorError*)&pc->failedreason);
52: *outits = its;
53: *reason = PCRICHARDSON_CONVERGED_ITS;
54: return 0;
55: }
57: PetscErrorCode PCSetFromOptions_SOR(PetscOptionItems *PetscOptionsObject,PC pc)
58: {
59: PC_SOR *jac = (PC_SOR*)pc->data;
60: PetscBool flg;
62: PetscOptionsHead(PetscOptionsObject,"(S)SOR options");
63: PetscOptionsReal("-pc_sor_omega","relaxation factor (0 < omega < 2)","PCSORSetOmega",jac->omega,&jac->omega,NULL);
64: PetscOptionsReal("-pc_sor_diagonal_shift","Add to the diagonal entries","",jac->fshift,&jac->fshift,NULL);
65: PetscOptionsInt("-pc_sor_its","number of inner SOR iterations","PCSORSetIterations",jac->its,&jac->its,NULL);
66: PetscOptionsInt("-pc_sor_lits","number of local inner SOR iterations","PCSORSetIterations",jac->lits,&jac->lits,NULL);
67: PetscOptionsBoolGroupBegin("-pc_sor_symmetric","SSOR, not SOR","PCSORSetSymmetric",&flg);
68: if (flg) PCSORSetSymmetric(pc,SOR_SYMMETRIC_SWEEP);
69: PetscOptionsBoolGroup("-pc_sor_backward","use backward sweep instead of forward","PCSORSetSymmetric",&flg);
70: if (flg) PCSORSetSymmetric(pc,SOR_BACKWARD_SWEEP);
71: PetscOptionsBoolGroup("-pc_sor_forward","use forward sweep","PCSORSetSymmetric",&flg);
72: if (flg) PCSORSetSymmetric(pc,SOR_FORWARD_SWEEP);
73: PetscOptionsBoolGroup("-pc_sor_local_symmetric","use SSOR separately on each processor","PCSORSetSymmetric",&flg);
74: if (flg) PCSORSetSymmetric(pc,SOR_LOCAL_SYMMETRIC_SWEEP);
75: PetscOptionsBoolGroup("-pc_sor_local_backward","use backward sweep locally","PCSORSetSymmetric",&flg);
76: if (flg) PCSORSetSymmetric(pc,SOR_LOCAL_BACKWARD_SWEEP);
77: PetscOptionsBoolGroupEnd("-pc_sor_local_forward","use forward sweep locally","PCSORSetSymmetric",&flg);
78: if (flg) PCSORSetSymmetric(pc,SOR_LOCAL_FORWARD_SWEEP);
79: PetscOptionsTail();
80: return 0;
81: }
83: PetscErrorCode PCView_SOR(PC pc,PetscViewer viewer)
84: {
85: PC_SOR *jac = (PC_SOR*)pc->data;
86: MatSORType sym = jac->sym;
87: const char *sortype;
88: PetscBool iascii;
90: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
91: if (iascii) {
92: if (sym & SOR_ZERO_INITIAL_GUESS) PetscViewerASCIIPrintf(viewer," zero initial guess\n");
93: if (sym == SOR_APPLY_UPPER) sortype = "apply_upper";
94: else if (sym == SOR_APPLY_LOWER) sortype = "apply_lower";
95: else if (sym & SOR_EISENSTAT) sortype = "Eisenstat";
96: else if ((sym & SOR_SYMMETRIC_SWEEP) == SOR_SYMMETRIC_SWEEP) sortype = "symmetric";
97: else if (sym & SOR_BACKWARD_SWEEP) sortype = "backward";
98: else if (sym & SOR_FORWARD_SWEEP) sortype = "forward";
99: else if ((sym & SOR_LOCAL_SYMMETRIC_SWEEP) == SOR_LOCAL_SYMMETRIC_SWEEP) sortype = "local_symmetric";
100: else if (sym & SOR_LOCAL_FORWARD_SWEEP) sortype = "local_forward";
101: else if (sym & SOR_LOCAL_BACKWARD_SWEEP) sortype = "local_backward";
102: else sortype = "unknown";
103: PetscViewerASCIIPrintf(viewer," type = %s, iterations = %D, local iterations = %D, omega = %g\n",sortype,jac->its,jac->lits,(double)jac->omega);
104: }
105: return 0;
106: }
108: /* ------------------------------------------------------------------------------*/
109: static PetscErrorCode PCSORSetSymmetric_SOR(PC pc,MatSORType flag)
110: {
111: PC_SOR *jac = (PC_SOR*)pc->data;
113: jac->sym = flag;
114: return 0;
115: }
117: static PetscErrorCode PCSORSetOmega_SOR(PC pc,PetscReal omega)
118: {
119: PC_SOR *jac = (PC_SOR*)pc->data;
122: jac->omega = omega;
123: return 0;
124: }
126: static PetscErrorCode PCSORSetIterations_SOR(PC pc,PetscInt its,PetscInt lits)
127: {
128: PC_SOR *jac = (PC_SOR*)pc->data;
130: jac->its = its;
131: jac->lits = lits;
132: return 0;
133: }
135: static PetscErrorCode PCSORGetSymmetric_SOR(PC pc,MatSORType *flag)
136: {
137: PC_SOR *jac = (PC_SOR*)pc->data;
139: *flag = jac->sym;
140: return 0;
141: }
143: static PetscErrorCode PCSORGetOmega_SOR(PC pc,PetscReal *omega)
144: {
145: PC_SOR *jac = (PC_SOR*)pc->data;
147: *omega = jac->omega;
148: return 0;
149: }
151: static PetscErrorCode PCSORGetIterations_SOR(PC pc,PetscInt *its,PetscInt *lits)
152: {
153: PC_SOR *jac = (PC_SOR*)pc->data;
155: if (its) *its = jac->its;
156: if (lits) *lits = jac->lits;
157: return 0;
158: }
160: /* ------------------------------------------------------------------------------*/
161: /*@
162: PCSORGetSymmetric - Gets the form the SOR preconditioner is using; backward, or forward relaxation. The local variants perform SOR on
163: each processor. By default forward relaxation is used.
165: Logically Collective on PC
167: Input Parameter:
168: . pc - the preconditioner context
170: Output Parameter:
171: . flag - one of the following
172: .vb
173: SOR_FORWARD_SWEEP
174: SOR_BACKWARD_SWEEP
175: SOR_SYMMETRIC_SWEEP
176: SOR_LOCAL_FORWARD_SWEEP
177: SOR_LOCAL_BACKWARD_SWEEP
178: SOR_LOCAL_SYMMETRIC_SWEEP
179: .ve
181: Options Database Keys:
182: + -pc_sor_symmetric - Activates symmetric version
183: . -pc_sor_backward - Activates backward version
184: . -pc_sor_local_forward - Activates local forward version
185: . -pc_sor_local_symmetric - Activates local symmetric version
186: - -pc_sor_local_backward - Activates local backward version
188: Notes:
189: To use the Eisenstat trick with SSOR, employ the PCEISENSTAT preconditioner,
190: which can be chosen with the option
191: . -pc_type eisenstat - Activates Eisenstat trick
193: Level: intermediate
195: .seealso: PCEisenstatSetOmega(), PCSORSetIterations(), PCSORSetOmega(), PCSORSetSymmetric()
196: @*/
197: PetscErrorCode PCSORGetSymmetric(PC pc,MatSORType *flag)
198: {
200: PetscUseMethod(pc,"PCSORGetSymmetric_C",(PC,MatSORType*),(pc,flag));
201: return 0;
202: }
204: /*@
205: PCSORGetOmega - Gets the SOR relaxation coefficient, omega
206: (where omega = 1.0 by default).
208: Logically Collective on PC
210: Input Parameter:
211: . pc - the preconditioner context
213: Output Parameter:
214: . omega - relaxation coefficient (0 < omega < 2).
216: Options Database Key:
217: . -pc_sor_omega <omega> - Sets omega
219: Level: intermediate
221: .seealso: PCSORSetSymmetric(), PCSORSetIterations(), PCEisenstatSetOmega(), PCSORSetOmega()
222: @*/
223: PetscErrorCode PCSORGetOmega(PC pc,PetscReal *omega)
224: {
226: PetscUseMethod(pc,"PCSORGetOmega_C",(PC,PetscReal*),(pc,omega));
227: return 0;
228: }
230: /*@
231: PCSORGetIterations - Gets the number of inner iterations to
232: be used by the SOR preconditioner. The default is 1.
234: Logically Collective on PC
236: Input Parameter:
237: . pc - the preconditioner context
239: Output Parameters:
240: + lits - number of local iterations, smoothings over just variables on processor
241: - its - number of parallel iterations to use; each parallel iteration has lits local iterations
243: Options Database Key:
244: + -pc_sor_its <its> - Sets number of iterations
245: - -pc_sor_lits <lits> - Sets number of local iterations
247: Level: intermediate
249: Notes:
250: When run on one processor the number of smoothings is lits*its
252: .seealso: PCSORSetOmega(), PCSORSetSymmetric(), PCSORSetIterations()
253: @*/
254: PetscErrorCode PCSORGetIterations(PC pc,PetscInt *its,PetscInt *lits)
255: {
257: PetscUseMethod(pc,"PCSORGetIterations_C",(PC,PetscInt*,PetscInt*),(pc,its,lits));
258: return 0;
259: }
261: /*@
262: PCSORSetSymmetric - Sets the SOR preconditioner to use symmetric (SSOR),
263: backward, or forward relaxation. The local variants perform SOR on
264: each processor. By default forward relaxation is used.
266: Logically Collective on PC
268: Input Parameters:
269: + pc - the preconditioner context
270: - flag - one of the following
271: .vb
272: SOR_FORWARD_SWEEP
273: SOR_BACKWARD_SWEEP
274: SOR_SYMMETRIC_SWEEP
275: SOR_LOCAL_FORWARD_SWEEP
276: SOR_LOCAL_BACKWARD_SWEEP
277: SOR_LOCAL_SYMMETRIC_SWEEP
278: .ve
280: Options Database Keys:
281: + -pc_sor_symmetric - Activates symmetric version
282: . -pc_sor_backward - Activates backward version
283: . -pc_sor_local_forward - Activates local forward version
284: . -pc_sor_local_symmetric - Activates local symmetric version
285: - -pc_sor_local_backward - Activates local backward version
287: Notes:
288: To use the Eisenstat trick with SSOR, employ the PCEISENSTAT preconditioner,
289: which can be chosen with the option
290: . -pc_type eisenstat - Activates Eisenstat trick
292: Level: intermediate
294: .seealso: PCEisenstatSetOmega(), PCSORSetIterations(), PCSORSetOmega()
295: @*/
296: PetscErrorCode PCSORSetSymmetric(PC pc,MatSORType flag)
297: {
300: PetscTryMethod(pc,"PCSORSetSymmetric_C",(PC,MatSORType),(pc,flag));
301: return 0;
302: }
304: /*@
305: PCSORSetOmega - Sets the SOR relaxation coefficient, omega
306: (where omega = 1.0 by default).
308: Logically Collective on PC
310: Input Parameters:
311: + pc - the preconditioner context
312: - omega - relaxation coefficient (0 < omega < 2).
314: Options Database Key:
315: . -pc_sor_omega <omega> - Sets omega
317: Level: intermediate
319: Note:
320: If omega != 1, you will need to set the MAT_USE_INODES option to PETSC_FALSE on the matrix.
322: .seealso: PCSORSetSymmetric(), PCSORSetIterations(), PCEisenstatSetOmega(), MatSetOption()
323: @*/
324: PetscErrorCode PCSORSetOmega(PC pc,PetscReal omega)
325: {
328: PetscTryMethod(pc,"PCSORSetOmega_C",(PC,PetscReal),(pc,omega));
329: return 0;
330: }
332: /*@
333: PCSORSetIterations - Sets the number of inner iterations to
334: be used by the SOR preconditioner. The default is 1.
336: Logically Collective on PC
338: Input Parameters:
339: + pc - the preconditioner context
340: . lits - number of local iterations, smoothings over just variables on processor
341: - its - number of parallel iterations to use; each parallel iteration has lits local iterations
343: Options Database Key:
344: + -pc_sor_its <its> - Sets number of iterations
345: - -pc_sor_lits <lits> - Sets number of local iterations
347: Level: intermediate
349: Notes:
350: When run on one processor the number of smoothings is lits*its
352: .seealso: PCSORSetOmega(), PCSORSetSymmetric()
353: @*/
354: PetscErrorCode PCSORSetIterations(PC pc,PetscInt its,PetscInt lits)
355: {
358: PetscTryMethod(pc,"PCSORSetIterations_C",(PC,PetscInt,PetscInt),(pc,its,lits));
359: return 0;
360: }
362: /*MC
363: PCSOR - (S)SOR (successive over relaxation, Gauss-Seidel) preconditioning
365: Options Database Keys:
366: + -pc_sor_symmetric - Activates symmetric version
367: . -pc_sor_backward - Activates backward version
368: . -pc_sor_forward - Activates forward version
369: . -pc_sor_local_forward - Activates local forward version
370: . -pc_sor_local_symmetric - Activates local symmetric version (default version)
371: . -pc_sor_local_backward - Activates local backward version
372: . -pc_sor_omega <omega> - Sets omega
373: . -pc_sor_diagonal_shift <shift> - shift the diagonal entries; useful if the matrix has zeros on the diagonal
374: . -pc_sor_its <its> - Sets number of iterations (default 1)
375: - -pc_sor_lits <lits> - Sets number of local iterations (default 1)
377: Level: beginner
379: Notes:
380: Only implemented for the AIJ and SeqBAIJ matrix formats.
381: Not a true parallel SOR, in parallel this implementation corresponds to block
382: Jacobi with SOR on each block.
384: For AIJ matrix if a diagonal entry is zero (and the diagonal shift is zero) then by default the inverse of that
385: zero will be used and hence the KSPSolve() will terminate with KSP_DIVERGED_NANORIF. If the option
386: KSPSetErrorIfNotConverged() or -ksp_error_if_not_converged the code will terminate as soon as it detects the
387: zero pivot.
389: For SeqBAIJ matrices this implements point-block SOR, but the omega, its, lits options are not supported.
391: For SeqBAIJ the diagonal blocks are inverted using dense LU with partial pivoting. If a zero pivot is detected
392: the computation is stopped with an error
394: If used with KSPRICHARDSON and no monitors the convergence test is skipped to improve speed, thus it always iterates
395: the maximum number of iterations you've selected for KSP. It is usually used in this mode as a smoother for multigrid.
397: If omega != 1, you will need to set the MAT_USE_INODES option to PETSC_FALSE on the matrix.
399: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
400: PCSORSetIterations(), PCSORSetSymmetric(), PCSORSetOmega(), PCEISENSTAT, MatSetOption()
401: M*/
403: PETSC_EXTERN PetscErrorCode PCCreate_SOR(PC pc)
404: {
405: PC_SOR *jac;
407: PetscNewLog(pc,&jac);
409: pc->ops->apply = PCApply_SOR;
410: pc->ops->applytranspose = PCApplyTranspose_SOR;
411: pc->ops->applyrichardson = PCApplyRichardson_SOR;
412: pc->ops->setfromoptions = PCSetFromOptions_SOR;
413: pc->ops->setup = NULL;
414: pc->ops->view = PCView_SOR;
415: pc->ops->destroy = PCDestroy_SOR;
416: pc->data = (void*)jac;
417: jac->sym = SOR_LOCAL_SYMMETRIC_SWEEP;
418: jac->omega = 1.0;
419: jac->fshift = 0.0;
420: jac->its = 1;
421: jac->lits = 1;
423: PetscObjectComposeFunction((PetscObject)pc,"PCSORSetSymmetric_C",PCSORSetSymmetric_SOR);
424: PetscObjectComposeFunction((PetscObject)pc,"PCSORSetOmega_C",PCSORSetOmega_SOR);
425: PetscObjectComposeFunction((PetscObject)pc,"PCSORSetIterations_C",PCSORSetIterations_SOR);
426: PetscObjectComposeFunction((PetscObject)pc,"PCSORGetSymmetric_C",PCSORGetSymmetric_SOR);
427: PetscObjectComposeFunction((PetscObject)pc,"PCSORGetOmega_C",PCSORGetOmega_SOR);
428: PetscObjectComposeFunction((PetscObject)pc,"PCSORGetIterations_C",PCSORGetIterations_SOR);
429: return 0;
430: }