Actual source code: sor.c
1: #define PETSCKSP_DLL
3: /*
4: Defines a (S)SOR preconditioner for any Mat implementation
5: */
6: #include private/pcimpl.h
8: typedef struct {
9: PetscInt its; /* inner iterations, number of sweeps */
10: PetscInt lits; /* local inner iterations, number of sweeps applied by the local matrix mat->A */
11: MatSORType sym; /* forward, reverse, symmetric etc. */
12: PetscReal omega;
13: PetscReal fshift;
14: } PC_SOR;
18: static PetscErrorCode PCDestroy_SOR(PC pc)
19: {
20: PC_SOR *jac = (PC_SOR*)pc->data;
24: PetscFree(jac);
25: return(0);
26: }
30: static PetscErrorCode PCApply_SOR(PC pc,Vec x,Vec y)
31: {
32: PC_SOR *jac = (PC_SOR*)pc->data;
34: PetscInt flag = jac->sym | SOR_ZERO_INITIAL_GUESS;
37: MatSOR(pc->pmat,x,jac->omega,(MatSORType)flag,jac->fshift,jac->its,jac->lits,y);
38: return(0);
39: }
43: static PetscErrorCode PCApplyRichardson_SOR(PC pc,Vec b,Vec y,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt its,PetscTruth guesszero,PetscInt *outits,PCRichardsonConvergedReason *reason)
44: {
45: PC_SOR *jac = (PC_SOR*)pc->data;
47: MatSORType stype = jac->sym;
50: PetscInfo1(pc,"Warning, convergence critera ignored, using %D iterations\n",its);
51: if (guesszero) {
52: stype = (MatSORType) (stype | SOR_ZERO_INITIAL_GUESS);
53: }
54: MatSOR(pc->pmat,b,jac->omega,stype,jac->fshift,its*jac->its,jac->lits,y);
55: *outits = its;
56: *reason = PCRICHARDSON_CONVERGED_ITS;
57: return(0);
58: }
62: PetscErrorCode PCSetFromOptions_SOR(PC pc)
63: {
64: PC_SOR *jac = (PC_SOR*)pc->data;
66: PetscTruth flg;
69: PetscOptionsHead("(S)SOR options");
70: PetscOptionsReal("-pc_sor_omega","relaxation factor (0 < omega < 2)","PCSORSetOmega",jac->omega,&jac->omega,0);
71: PetscOptionsReal("-pc_sor_diagonal_shift","Add to the diagonal entries","",jac->fshift,&jac->fshift,0);
72: PetscOptionsInt("-pc_sor_its","number of inner SOR iterations","PCSORSetIterations",jac->its,&jac->its,0);
73: PetscOptionsInt("-pc_sor_lits","number of local inner SOR iterations","PCSORSetIterations",jac->lits,&jac->lits,0);
74: PetscOptionsTruthGroupBegin("-pc_sor_symmetric","SSOR, not SOR","PCSORSetSymmetric",&flg);
75: if (flg) {PCSORSetSymmetric(pc,SOR_SYMMETRIC_SWEEP);}
76: PetscOptionsTruthGroup("-pc_sor_backward","use backward sweep instead of forward","PCSORSetSymmetric",&flg);
77: if (flg) {PCSORSetSymmetric(pc,SOR_BACKWARD_SWEEP);}
78: PetscOptionsTruthGroup("-pc_sor_forward","use forward sweep","PCSORSetSymmetric",&flg);
79: if (flg) {PCSORSetSymmetric(pc,SOR_FORWARD_SWEEP);}
80: PetscOptionsTruthGroup("-pc_sor_local_symmetric","use SSOR separately on each processor","PCSORSetSymmetric",&flg);
81: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_SYMMETRIC_SWEEP);}
82: PetscOptionsTruthGroup("-pc_sor_local_backward","use backward sweep locally","PCSORSetSymmetric",&flg);
83: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_BACKWARD_SWEEP);}
84: PetscOptionsTruthGroupEnd("-pc_sor_local_forward","use forward sweep locally","PCSORSetSymmetric",&flg);
85: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_FORWARD_SWEEP);}
86: PetscOptionsTail();
87: return(0);
88: }
92: PetscErrorCode PCView_SOR(PC pc,PetscViewer viewer)
93: {
94: PC_SOR *jac = (PC_SOR*)pc->data;
95: MatSORType sym = jac->sym;
96: const char *sortype;
98: PetscTruth iascii;
101: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
102: if (iascii) {
103: if (sym & SOR_ZERO_INITIAL_GUESS) {PetscViewerASCIIPrintf(viewer," SOR: zero initial guess\n");}
104: if (sym == SOR_APPLY_UPPER) sortype = "apply_upper";
105: else if (sym == SOR_APPLY_LOWER) sortype = "apply_lower";
106: else if (sym & SOR_EISENSTAT) sortype = "Eisenstat";
107: else if ((sym & SOR_SYMMETRIC_SWEEP) == SOR_SYMMETRIC_SWEEP)
108: sortype = "symmetric";
109: else if (sym & SOR_BACKWARD_SWEEP) sortype = "backward";
110: else if (sym & SOR_FORWARD_SWEEP) sortype = "forward";
111: else if ((sym & SOR_LOCAL_SYMMETRIC_SWEEP) == SOR_LOCAL_SYMMETRIC_SWEEP)
112: 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," SOR: type = %s, iterations = %D, local iterations = %D, omega = %G\n",sortype,jac->its,jac->lits,jac->omega);
117: } else {
118: SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PCSOR",((PetscObject)viewer)->type_name);
119: }
120: return(0);
121: }
124: /* ------------------------------------------------------------------------------*/
128: PetscErrorCode PCSORSetSymmetric_SOR(PC pc,MatSORType flag)
129: {
130: PC_SOR *jac;
133: jac = (PC_SOR*)pc->data;
134: jac->sym = flag;
135: return(0);
136: }
142: PetscErrorCode PCSORSetOmega_SOR(PC pc,PetscReal omega)
143: {
144: PC_SOR *jac;
147: if (omega >= 2.0 || omega <= 0.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Relaxation out of range");
148: jac = (PC_SOR*)pc->data;
149: jac->omega = omega;
150: return(0);
151: }
157: PetscErrorCode PCSORSetIterations_SOR(PC pc,PetscInt its,PetscInt lits)
158: {
159: PC_SOR *jac;
162: jac = (PC_SOR*)pc->data;
163: jac->its = its;
164: jac->lits = lits;
165: return(0);
166: }
169: /* ------------------------------------------------------------------------------*/
172: /*@
173: PCSORSetSymmetric - Sets the SOR preconditioner to use symmetric (SSOR),
174: backward, or forward relaxation. The local variants perform SOR on
175: each processor. By default forward relaxation is used.
177: Collective on PC
179: Input Parameters:
180: + pc - the preconditioner context
181: - flag - one of the following
182: .vb
183: SOR_FORWARD_SWEEP
184: SOR_BACKWARD_SWEEP
185: SOR_SYMMETRIC_SWEEP
186: SOR_LOCAL_FORWARD_SWEEP
187: SOR_LOCAL_BACKWARD_SWEEP
188: SOR_LOCAL_SYMMETRIC_SWEEP
189: .ve
191: Options Database Keys:
192: + -pc_sor_symmetric - Activates symmetric version
193: . -pc_sor_backward - Activates backward version
194: . -pc_sor_local_forward - Activates local forward version
195: . -pc_sor_local_symmetric - Activates local symmetric version
196: - -pc_sor_local_backward - Activates local backward version
198: Notes:
199: To use the Eisenstat trick with SSOR, employ the PCEISENSTAT preconditioner,
200: which can be chosen with the option
201: . -pc_type eisenstat - Activates Eisenstat trick
203: Level: intermediate
205: .keywords: PC, SOR, SSOR, set, relaxation, sweep, forward, backward, symmetric
207: .seealso: PCEisenstatSetOmega(), PCSORSetIterations(), PCSORSetOmega()
208: @*/
209: PetscErrorCode PCSORSetSymmetric(PC pc,MatSORType flag)
210: {
211: PetscErrorCode ierr,(*f)(PC,MatSORType);
215: PetscObjectQueryFunction((PetscObject)pc,"PCSORSetSymmetric_C",(void (**)(void))&f);
216: if (f) {
217: (*f)(pc,flag);
218: }
219: return(0);
220: }
224: /*@
225: PCSORSetOmega - Sets the SOR relaxation coefficient, omega
226: (where omega = 1.0 by default).
228: Collective on PC
230: Input Parameters:
231: + pc - the preconditioner context
232: - omega - relaxation coefficient (0 < omega < 2).
234: Options Database Key:
235: . -pc_sor_omega <omega> - Sets omega
237: Level: intermediate
239: .keywords: PC, SOR, SSOR, set, relaxation, omega
241: .seealso: PCSORSetSymmetric(), PCSORSetIterations(), PCEisenstatSetOmega()
242: @*/
243: PetscErrorCode PCSORSetOmega(PC pc,PetscReal omega)
244: {
245: PetscErrorCode ierr,(*f)(PC,PetscReal);
248: PetscObjectQueryFunction((PetscObject)pc,"PCSORSetOmega_C",(void (**)(void))&f);
249: if (f) {
250: (*f)(pc,omega);
251: }
252: return(0);
253: }
257: /*@
258: PCSORSetIterations - Sets the number of inner iterations to
259: be used by the SOR preconditioner. The default is 1.
261: Collective on PC
263: Input Parameters:
264: + pc - the preconditioner context
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: When run on one processor the number of smoothings is lits*its
276: .keywords: PC, SOR, SSOR, set, iterations
278: .seealso: PCSORSetOmega(), PCSORSetSymmetric()
279: @*/
280: PetscErrorCode PCSORSetIterations(PC pc,PetscInt its,PetscInt lits)
281: {
282: PetscErrorCode ierr,(*f)(PC,PetscInt,PetscInt);
286: PetscObjectQueryFunction((PetscObject)pc,"PCSORSetIterations_C",(void (**)(void))&f);
287: if (f) {
288: (*f)(pc,its,lits);
289: }
290: return(0);
291: }
293: /*MC
294: PCSOR - (S)SOR (successive over relaxation, Gauss-Seidel) preconditioning
296: Options Database Keys:
297: + -pc_sor_symmetric - Activates symmetric version
298: . -pc_sor_backward - Activates backward version
299: . -pc_sor_forward - Activates forward version
300: . -pc_sor_local_forward - Activates local forward version
301: . -pc_sor_local_symmetric - Activates local symmetric version (default version)
302: . -pc_sor_local_backward - Activates local backward version
303: . -pc_sor_omega <omega> - Sets omega
304: . -pc_sor_its <its> - Sets number of iterations (default 1)
305: - -pc_sor_lits <lits> - Sets number of local iterations (default 1)
307: Level: beginner
309: Concepts: SOR, preconditioners, Gauss-Seidel
311: Notes: Only implemented for the AIJ and SeqBAIJ matrix formats.
312: Not a true parallel SOR, in parallel this implementation corresponds to block
313: Jacobi with SOR on each block.
315: For SeqBAIJ matrices this implements point-block SOR, but the omega, its, lits options are not supported.
317: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
318: PCSORSetIterations(), PCSORSetSymmetric(), PCSORSetOmega(), PCEISENSTAT
319: M*/
324: PetscErrorCode PCCreate_SOR(PC pc)
325: {
327: PC_SOR *jac;
330: PetscNewLog(pc,PC_SOR,&jac);
332: pc->ops->apply = PCApply_SOR;
333: pc->ops->applyrichardson = PCApplyRichardson_SOR;
334: pc->ops->setfromoptions = PCSetFromOptions_SOR;
335: pc->ops->setup = 0;
336: pc->ops->view = PCView_SOR;
337: pc->ops->destroy = PCDestroy_SOR;
338: pc->data = (void*)jac;
339: jac->sym = SOR_LOCAL_SYMMETRIC_SWEEP;
340: jac->omega = 1.0;
341: jac->fshift = 0.0;
342: jac->its = 1;
343: jac->lits = 1;
345: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCSORSetSymmetric_C","PCSORSetSymmetric_SOR",
346: PCSORSetSymmetric_SOR);
347: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCSORSetOmega_C","PCSORSetOmega_SOR",
348: PCSORSetOmega_SOR);
349: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCSORSetIterations_C","PCSORSetIterations_SOR",
350: PCSORSetIterations_SOR);
352: return(0);
353: }