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: static int PCDestroy_SOR(PC pc)
16: {
17: PC_SOR *jac = (PC_SOR*)pc->data;
18: int ierr;
21: PetscFree(jac);
22: return(0);
23: }
25: static int PCApply_SOR(PC pc,Vec x,Vec y)
26: {
27: PC_SOR *jac = (PC_SOR*)pc->data;
28: int ierr,flag = jac->sym | SOR_ZERO_INITIAL_GUESS;
31: MatRelax(pc->pmat,x,jac->omega,(MatSORType)flag,0.0,jac->its,jac->lits,y);
32: return(0);
33: }
35: static int PCApplyRichardson_SOR(PC pc,Vec b,Vec y,Vec w,PetscReal rtol,PetscReal atol, PetscReal dtol,int its)
36: {
37: PC_SOR *jac = (PC_SOR*)pc->data;
38: int ierr;
41: PetscLogInfo(pc,"PCApplyRichardson_SOR: Warning, convergence critera ignored, using %d iterationsn",its);
42: MatRelax(pc->mat,b,jac->omega,(MatSORType)jac->sym,0.0,its,jac->lits,y);
43: return(0);
44: }
46: static int PCSetFromOptions_SOR(PC pc)
47: {
48: PC_SOR *jac = (PC_SOR*)pc->data;
49: int ierr;
50: PetscTruth flg;
53: PetscOptionsHead("(S)SOR options");
54: PetscOptionsReal("-pc_sor_omega","relaxation factor (0 < omega < 2)","PCSORSetOmega",jac->omega,&jac->omega,0);
55: PetscOptionsInt("-pc_sor_its","number of inner SOR iterations","PCSORSetIterations",jac->its,&jac->its,0);
56: PetscOptionsInt("-pc_sor_lits","number of local inner SOR iterations","PCSORSetIterations",jac->lits,&jac->lits,0);
57: PetscOptionsLogicalGroupBegin("-pc_sor_symmetric","SSOR, not SOR","PCSORSetSymmetric",&flg);
58: if (flg) {PCSORSetSymmetric(pc,SOR_SYMMETRIC_SWEEP);}
59: PetscOptionsLogicalGroup("-pc_sor_backward","use backward sweep instead of forward","PCSORSetSymmetric",&flg);
60: if (flg) {PCSORSetSymmetric(pc,SOR_BACKWARD_SWEEP);}
61: PetscOptionsLogicalGroup("-pc_sor_local_symmetric","use SSOR seperately on each processor","PCSORSetSymmetric",&flg);
62: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_SYMMETRIC_SWEEP);}
63: PetscOptionsLogicalGroup("-pc_sor_local_backward","use backward sweep locally","PCSORSetSymmetric",&flg);
64: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_BACKWARD_SWEEP);}
65: PetscOptionsLogicalGroupEnd("-pc_sor_local_forward","use forward sweep locally","PCSORSetSymmetric",&flg);
66: if (flg) {PCSORSetSymmetric(pc,SOR_LOCAL_FORWARD_SWEEP);}
67: PetscOptionsTail();
68: return(0);
69: }
71: static int PCView_SOR(PC pc,PetscViewer viewer)
72: {
73: PC_SOR *jac = (PC_SOR*)pc->data;
74: MatSORType sym = jac->sym;
75: char *sortype;
76: int ierr;
77: PetscTruth isascii;
80: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
81: if (isascii) {
82: if (sym & SOR_ZERO_INITIAL_GUESS) {PetscViewerASCIIPrintf(viewer," SOR: zero initial guessn");}
83: if (sym == SOR_APPLY_UPPER) sortype = "apply_upper";
84: else if (sym == SOR_APPLY_LOWER) sortype = "apply_lower";
85: else if (sym & SOR_EISENSTAT) sortype = "Eisenstat";
86: else if ((sym & SOR_SYMMETRIC_SWEEP) == SOR_SYMMETRIC_SWEEP)
87: sortype = "symmetric";
88: else if (sym & SOR_BACKWARD_SWEEP) sortype = "backward";
89: else if (sym & SOR_FORWARD_SWEEP) sortype = "forward";
90: else if ((sym & SOR_LOCAL_SYMMETRIC_SWEEP) == SOR_LOCAL_SYMMETRIC_SWEEP)
91: sortype = "local_symmetric";
92: else if (sym & SOR_LOCAL_FORWARD_SWEEP) sortype = "local_forward";
93: else if (sym & SOR_LOCAL_BACKWARD_SWEEP) sortype = "local_backward";
94: else sortype = "unknown";
95: PetscViewerASCIIPrintf(viewer," SOR: type = %s, iterations = %d, omega = %gn",sortype,jac->its,jac->omega);
96: } else {
97: SETERRQ1(1,"Viewer type %s not supported for PCSOR",((PetscObject)viewer)->type_name);
98: }
99: return(0);
100: }
103: /* ------------------------------------------------------------------------------*/
104: EXTERN_C_BEGIN
105: int PCSORSetSymmetric_SOR(PC pc,MatSORType flag)
106: {
107: PC_SOR *jac;
110: jac = (PC_SOR*)pc->data;
111: jac->sym = flag;
112: return(0);
113: }
114: EXTERN_C_END
116: EXTERN_C_BEGIN
117: int PCSORSetOmega_SOR(PC pc,PetscReal omega)
118: {
119: PC_SOR *jac;
122: if (omega >= 2.0 || omega <= 0.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Relaxation out of range");
123: jac = (PC_SOR*)pc->data;
124: jac->omega = omega;
125: return(0);
126: }
127: EXTERN_C_END
129: EXTERN_C_BEGIN
130: int PCSORSetIterations_SOR(PC pc,int its,int lits)
131: {
132: PC_SOR *jac;
135: jac = (PC_SOR*)pc->data;
136: jac->its = its;
137: jac->lits = lits;
138: return(0);
139: }
140: EXTERN_C_END
142: /* ------------------------------------------------------------------------------*/
143: /*@
144: PCSORSetSymmetric - Sets the SOR preconditioner to use symmetric (SSOR),
145: backward, or forward relaxation. The local variants perform SOR on
146: each processor. By default forward relaxation is used.
148: Collective on PC
150: Input Parameters:
151: + pc - the preconditioner context
152: - flag - one of the following
153: .vb
154: SOR_FORWARD_SWEEP
155: SOR_BACKWARD_SWEEP
156: SOR_SYMMETRIC_SWEEP
157: SOR_LOCAL_FORWARD_SWEEP
158: SOR_LOCAL_BACKWARD_SWEEP
159: SOR_LOCAL_SYMMETRIC_SWEEP
160: .ve
162: Options Database Keys:
163: + -pc_sor_symmetric - Activates symmetric version
164: . -pc_sor_backward - Activates backward version
165: . -pc_sor_local_forward - Activates local forward version
166: . -pc_sor_local_symmetric - Activates local symmetric version
167: - -pc_sor_local_backward - Activates local backward version
169: Notes:
170: To use the Eisenstat trick with SSOR, employ the PCEISENSTAT preconditioner,
171: which can be chosen with the option
172: . -pc_type eisenstat - Activates Eisenstat trick
174: Level: intermediate
176: .keywords: PC, SOR, SSOR, set, relaxation, sweep, forward, backward, symmetric
178: .seealso: PCEisenstatSetOmega(), PCSORSetIterations(), PCSORSetOmega()
179: @*/
180: int PCSORSetSymmetric(PC pc,MatSORType flag)
181: {
182: int ierr,(*f)(PC,MatSORType);
186: PetscObjectQueryFunction((PetscObject)pc,"PCSORSetSymmetric_C",(void (**)(void))&f);
187: if (f) {
188: (*f)(pc,flag);
189: }
190: return(0);
191: }
193: /*@
194: PCSORSetOmega - Sets the SOR relaxation coefficient, omega
195: (where omega = 1.0 by default).
197: Collective on PC
199: Input Parameters:
200: + pc - the preconditioner context
201: - omega - relaxation coefficient (0 < omega < 2).
203: Options Database Key:
204: . -pc_sor_omega <omega> - Sets omega
206: Level: intermediate
208: .keywords: PC, SOR, SSOR, set, relaxation, omega
210: .seealso: PCSORSetSymmetric(), PCSORSetIterations(), PCEisenstatSetOmega()
211: @*/
212: int PCSORSetOmega(PC pc,PetscReal omega)
213: {
214: int ierr,(*f)(PC,PetscReal);
217: PetscObjectQueryFunction((PetscObject)pc,"PCSORSetOmega_C",(void (**)(void))&f);
218: if (f) {
219: (*f)(pc,omega);
220: }
221: return(0);
222: }
224: /*@
225: PCSORSetIterations - Sets the number of inner iterations to
226: be used by the SOR preconditioner. The default is 1.
228: Collective on PC
230: Input Parameters:
231: + pc - the preconditioner context
232: . lits - number of local iterations, smoothings over just variables on processor
233: - its - number of parallel iterations to use; each parallel iteration has lits local iterations
235: Options Database Key:
236: + -pc_sor_its <its> - Sets number of iterations
237: - -pc_sor_lits <lits> - Sets number of local iterations
239: Level: intermediate
241: Notes: When run on one processor the number of smoothings is lits*its
243: .keywords: PC, SOR, SSOR, set, iterations
245: .seealso: PCSORSetOmega(), PCSORSetSymmetric()
246: @*/
247: int PCSORSetIterations(PC pc,int its,int lits)
248: {
249: int ierr,(*f)(PC,int,int);
253: PetscObjectQueryFunction((PetscObject)pc,"PCSORSetIterations_C",(void (**)(void))&f);
254: if (f) {
255: (*f)(pc,its,lits);
256: }
257: return(0);
258: }
260: EXTERN_C_BEGIN
261: int PCCreate_SOR(PC pc)
262: {
263: int ierr;
264: PC_SOR *jac;
267: PetscNew(PC_SOR,&jac);
268: PetscLogObjectMemory(pc,sizeof(PC_SOR));
270: pc->ops->apply = PCApply_SOR;
271: pc->ops->applyrichardson = PCApplyRichardson_SOR;
272: pc->ops->setfromoptions = PCSetFromOptions_SOR;
273: pc->ops->setup = 0;
274: pc->ops->view = PCView_SOR;
275: pc->ops->destroy = PCDestroy_SOR;
276: pc->data = (void*)jac;
277: jac->sym = SOR_FORWARD_SWEEP;
278: jac->omega = 1.0;
279: jac->its = 1;
280: jac->lits = 1;
282: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCSORSetSymmetric_C","PCSORSetSymmetric_SOR",
283: PCSORSetSymmetric_SOR);
284: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCSORSetOmega_C","PCSORSetOmega_SOR",
285: PCSORSetOmega_SOR);
286: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCSORSetIterations_C","PCSORSetIterations_SOR",
287: PCSORSetIterations_SOR);
289: return(0);
290: }
291: EXTERN_C_END