Actual source code: sor.c
1: /*$Id: sor.c,v 1.103 2001/08/21 21:03:17 bsmith Exp $*/
3: /*
4: Defines a (S)SOR preconditioner for any Mat implementation
5: */
6: #include src/sles/pc/pcimpl.h
8: typedef struct {
9: int its; /* inner iterations, number of sweeps */
10: int 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: } PC_SOR;
15: #undef __FUNCT__
17: static int PCDestroy_SOR(PC pc)
18: {
19: PC_SOR *jac = (PC_SOR*)pc->data;
20: int ierr;
23: PetscFree(jac);
24: return(0);
25: }
27: #undef __FUNCT__
29: static int PCApply_SOR(PC pc,Vec x,Vec y)
30: {
31: PC_SOR *jac = (PC_SOR*)pc->data;
32: int ierr,flag = jac->sym | SOR_ZERO_INITIAL_GUESS;
35: MatRelax(pc->pmat,x,jac->omega,(MatSORType)flag,0.0,jac->its,jac->lits,y);
36: return(0);
37: }
39: #undef __FUNCT__
41: static int PCApplyRichardson_SOR(PC pc,Vec b,Vec y,Vec w,PetscReal rtol,PetscReal atol, PetscReal dtol,int its)
42: {
43: PC_SOR *jac = (PC_SOR*)pc->data;
44: int ierr;
47: PetscLogInfo(pc,"PCApplyRichardson_SOR: Warning, convergence critera ignored, using %d iterationsn",its);
48: MatRelax(pc->mat,b,jac->omega,(MatSORType)jac->sym,0.0,its,jac->lits,y);
49: return(0);
50: }
52: #undef __FUNCT__
54: static int PCSetFromOptions_SOR(PC pc)
55: {
56: PC_SOR *jac = (PC_SOR*)pc->data;
57: int ierr;
58: PetscTruth flg;
61: PetscOptionsHead("(S)SOR options");
62: PetscOptionsReal("-pc_sor_omega","relaxation factor (0 < omega < 2)","PCSORSetOmega",jac->omega,&jac->omega,0);
63: PetscOptionsInt("-pc_sor_its","number of inner SOR iterations","PCSORSetIterations",jac->its,&jac->its,0);
64: PetscOptionsInt("-pc_sor_lits","number of local inner SOR iterations","PCSORSetIterations",jac->lits,&jac->lits,0);
65: PetscOptionsLogicalGroupBegin("-pc_sor_symmetric","SSOR, not SOR","PCSORSetSymmetric",&flg);
66: if (flg) {PCSORSetSymmetric(pc,SOR_SYMMETRIC_SWEEP);}
67: PetscOptionsLogicalGroup("-pc_sor_backward","use backward sweep instead of forward","PCSORSetSymmetric",&flg);
68: if (flg) {PCSORSetSymmetric(pc,SOR_BACKWARD_SWEEP);}
69: PetscOptionsLogicalGroup("-pc_sor_local_symmetric","use SSOR seperately on each processor","PCSORSetSymmetric",&flg);
70: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_SYMMETRIC_SWEEP);}
71: PetscOptionsLogicalGroup("-pc_sor_local_backward","use backward sweep locally","PCSORSetSymmetric",&flg);
72: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_BACKWARD_SWEEP);}
73: PetscOptionsLogicalGroupEnd("-pc_sor_local_forward","use forward sweep locally","PCSORSetSymmetric",&flg);
74: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_FORWARD_SWEEP);}
75: PetscOptionsTail();
76: return(0);
77: }
79: #undef __FUNCT__
81: static int PCView_SOR(PC pc,PetscViewer viewer)
82: {
83: PC_SOR *jac = (PC_SOR*)pc->data;
84: MatSORType sym = jac->sym;
85: char *sortype;
86: int ierr;
87: PetscTruth isascii;
90: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
91: if (isascii) {
92: if (sym & SOR_ZERO_INITIAL_GUESS) {PetscViewerASCIIPrintf(viewer," SOR: zero initial guessn");}
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)
97: sortype = "symmetric";
98: else if (sym & SOR_BACKWARD_SWEEP) sortype = "backward";
99: else if (sym & SOR_FORWARD_SWEEP) sortype = "forward";
100: else if ((sym & SOR_LOCAL_SYMMETRIC_SWEEP) == SOR_LOCAL_SYMMETRIC_SWEEP)
101: sortype = "local_symmetric";
102: else if (sym & SOR_LOCAL_FORWARD_SWEEP) sortype = "local_forward";
103: else if (sym & SOR_LOCAL_BACKWARD_SWEEP) sortype = "local_backward";
104: else sortype = "unknown";
105: PetscViewerASCIIPrintf(viewer," SOR: type = %s, iterations = %d, omega = %gn",sortype,jac->its,jac->omega);
106: } else {
107: SETERRQ1(1,"Viewer type %s not supported for PCSOR",((PetscObject)viewer)->type_name);
108: }
109: return(0);
110: }
113: /* ------------------------------------------------------------------------------*/
114: EXTERN_C_BEGIN
115: #undef __FUNCT__
117: int PCSORSetSymmetric_SOR(PC pc,MatSORType flag)
118: {
119: PC_SOR *jac;
122: jac = (PC_SOR*)pc->data;
123: jac->sym = flag;
124: return(0);
125: }
126: EXTERN_C_END
128: EXTERN_C_BEGIN
129: #undef __FUNCT__
131: int PCSORSetOmega_SOR(PC pc,PetscReal omega)
132: {
133: PC_SOR *jac;
136: if (omega >= 2.0 || omega <= 0.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Relaxation out of range");
137: jac = (PC_SOR*)pc->data;
138: jac->omega = omega;
139: return(0);
140: }
141: EXTERN_C_END
143: EXTERN_C_BEGIN
144: #undef __FUNCT__
146: int PCSORSetIterations_SOR(PC pc,int its,int lits)
147: {
148: PC_SOR *jac;
151: jac = (PC_SOR*)pc->data;
152: jac->its = its;
153: jac->lits = lits;
154: return(0);
155: }
156: EXTERN_C_END
158: /* ------------------------------------------------------------------------------*/
159: #undef __FUNCT__
161: /*@
162: PCSORSetSymmetric - Sets the SOR preconditioner to use symmetric (SSOR),
163: backward, or forward relaxation. The local variants perform SOR on
164: each processor. By default forward relaxation is used.
166: Collective on PC
168: Input Parameters:
169: + pc - the preconditioner context
170: - flag - one of the following
171: .vb
172: SOR_FORWARD_SWEEP
173: SOR_BACKWARD_SWEEP
174: SOR_SYMMETRIC_SWEEP
175: SOR_LOCAL_FORWARD_SWEEP
176: SOR_LOCAL_BACKWARD_SWEEP
177: SOR_LOCAL_SYMMETRIC_SWEEP
178: .ve
180: Options Database Keys:
181: + -pc_sor_symmetric - Activates symmetric version
182: . -pc_sor_backward - Activates backward version
183: . -pc_sor_local_forward - Activates local forward version
184: . -pc_sor_local_symmetric - Activates local symmetric version
185: - -pc_sor_local_backward - Activates local backward version
187: Notes:
188: To use the Eisenstat trick with SSOR, employ the PCEISENSTAT preconditioner,
189: which can be chosen with the option
190: . -pc_type eisenstat - Activates Eisenstat trick
192: Level: intermediate
194: .keywords: PC, SOR, SSOR, set, relaxation, sweep, forward, backward, symmetric
196: .seealso: PCEisenstatSetOmega(), PCSORSetIterations(), PCSORSetOmega()
197: @*/
198: int PCSORSetSymmetric(PC pc,MatSORType flag)
199: {
200: int ierr,(*f)(PC,MatSORType);
204: PetscObjectQueryFunction((PetscObject)pc,"PCSORSetSymmetric_C",(void (**)(void))&f);
205: if (f) {
206: (*f)(pc,flag);
207: }
208: return(0);
209: }
211: #undef __FUNCT__
213: /*@
214: PCSORSetOmega - Sets the SOR relaxation coefficient, omega
215: (where omega = 1.0 by default).
217: Collective on PC
219: Input Parameters:
220: + pc - the preconditioner context
221: - omega - relaxation coefficient (0 < omega < 2).
223: Options Database Key:
224: . -pc_sor_omega <omega> - Sets omega
226: Level: intermediate
228: .keywords: PC, SOR, SSOR, set, relaxation, omega
230: .seealso: PCSORSetSymmetric(), PCSORSetIterations(), PCEisenstatSetOmega()
231: @*/
232: int PCSORSetOmega(PC pc,PetscReal omega)
233: {
234: int ierr,(*f)(PC,PetscReal);
237: PetscObjectQueryFunction((PetscObject)pc,"PCSORSetOmega_C",(void (**)(void))&f);
238: if (f) {
239: (*f)(pc,omega);
240: }
241: return(0);
242: }
244: #undef __FUNCT__
246: /*@
247: PCSORSetIterations - Sets the number of inner iterations to
248: be used by the SOR preconditioner. The default is 1.
250: Collective on PC
252: Input Parameters:
253: + pc - the preconditioner context
254: . lits - number of local iterations, smoothings over just variables on processor
255: - its - number of parallel iterations to use; each parallel iteration has lits local iterations
257: Options Database Key:
258: + -pc_sor_its <its> - Sets number of iterations
259: - -pc_sor_lits <lits> - Sets number of local iterations
261: Level: intermediate
263: Notes: When run on one processor the number of smoothings is lits*its
265: .keywords: PC, SOR, SSOR, set, iterations
267: .seealso: PCSORSetOmega(), PCSORSetSymmetric()
268: @*/
269: int PCSORSetIterations(PC pc,int its,int lits)
270: {
271: int ierr,(*f)(PC,int,int);
275: PetscObjectQueryFunction((PetscObject)pc,"PCSORSetIterations_C",(void (**)(void))&f);
276: if (f) {
277: (*f)(pc,its,lits);
278: }
279: return(0);
280: }
282: /*MC
283: PCSOR - (S)SOR (successive over relaxation, Gauss-Seidel) preconditioning
285: Options Database Keys:
286: + -pc_sor_symmetric - Activates symmetric version
287: . -pc_sor_backward - Activates backward version
288: . -pc_sor_local_forward - Activates local forward version
289: . -pc_sor_local_symmetric - Activates local symmetric version
290: . -pc_sor_local_backward - Activates local backward version
291: . -pc_sor_omega <omega> - Sets omega
292: . -pc_sor_its <its> - Sets number of iterations
293: - -pc_sor_lits <lits> - Sets number of local iterations
295: Level: beginner
297: Concepts: SOR, preconditioners, Gauss-Seidel
299: Notes: Only implemented for the AIJ matrix format.
300: Not a true parallel SOR, in parallel this implementation corresponds to block
301: Jacobi with SOR on each block.
303: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
304: PCSORSetIterations(), PCSORSetSymmetric(), PCSORSetOmega(), PCEISENSTAT
305: M*/
307: EXTERN_C_BEGIN
308: #undef __FUNCT__
310: int PCCreate_SOR(PC pc)
311: {
312: int ierr;
313: PC_SOR *jac;
316: PetscNew(PC_SOR,&jac);
317: PetscLogObjectMemory(pc,sizeof(PC_SOR));
319: pc->ops->apply = PCApply_SOR;
320: pc->ops->applyrichardson = PCApplyRichardson_SOR;
321: pc->ops->setfromoptions = PCSetFromOptions_SOR;
322: pc->ops->setup = 0;
323: pc->ops->view = PCView_SOR;
324: pc->ops->destroy = PCDestroy_SOR;
325: pc->data = (void*)jac;
326: jac->sym = SOR_FORWARD_SWEEP;
327: jac->omega = 1.0;
328: jac->its = 1;
329: jac->lits = 1;
331: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCSORSetSymmetric_C","PCSORSetSymmetric_SOR",
332: PCSORSetSymmetric_SOR);
333: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCSORSetOmega_C","PCSORSetOmega_SOR",
334: PCSORSetOmega_SOR);
335: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCSORSetIterations_C","PCSORSetIterations_SOR",
336: PCSORSetIterations_SOR);
338: return(0);
339: }
340: EXTERN_C_END