1 /*
2 Defines a (S)SOR preconditioner for any Mat implementation
3 */
4 #include <petsc/private/pcimpl.h> /*I "petscpc.h" I*/
5
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;
13
PCDestroy_SOR(PC pc)14 static PetscErrorCode PCDestroy_SOR(PC pc)
15 {
16 PetscFunctionBegin;
17 PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSORSetSymmetric_C", NULL));
18 PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSORSetOmega_C", NULL));
19 PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSORSetIterations_C", NULL));
20 PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSORGetSymmetric_C", NULL));
21 PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSORGetOmega_C", NULL));
22 PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSORGetIterations_C", NULL));
23 PetscCall(PetscFree(pc->data));
24 PetscFunctionReturn(PETSC_SUCCESS);
25 }
26
PCApply_SOR(PC pc,Vec x,Vec y)27 static PetscErrorCode PCApply_SOR(PC pc, Vec x, Vec y)
28 {
29 PC_SOR *jac = (PC_SOR *)pc->data;
30 PetscInt flag = jac->sym | SOR_ZERO_INITIAL_GUESS;
31
32 PetscFunctionBegin;
33 PetscCall(MatSOR(pc->pmat, x, jac->omega, (MatSORType)flag, jac->fshift, jac->its, jac->lits, y));
34 PetscCall(MatFactorGetError(pc->pmat, (MatFactorError *)&pc->failedreason));
35 PetscFunctionReturn(PETSC_SUCCESS);
36 }
37
PCApplyTranspose_SOR(PC pc,Vec x,Vec y)38 static PetscErrorCode PCApplyTranspose_SOR(PC pc, Vec x, Vec y)
39 {
40 PC_SOR *jac = (PC_SOR *)pc->data;
41 PetscInt flag = jac->sym | SOR_ZERO_INITIAL_GUESS;
42 PetscBool set, sym;
43
44 PetscFunctionBegin;
45 PetscCall(MatIsSymmetricKnown(pc->pmat, &set, &sym));
46 PetscCheck(set && sym && (jac->sym == SOR_SYMMETRIC_SWEEP || jac->sym == SOR_LOCAL_SYMMETRIC_SWEEP), PetscObjectComm((PetscObject)pc), PETSC_ERR_SUP, "Can only apply transpose of SOR if matrix is symmetric and sweep is symmetric");
47 PetscCall(MatSOR(pc->pmat, x, jac->omega, (MatSORType)flag, jac->fshift, jac->its, jac->lits, y));
48 PetscCall(MatFactorGetError(pc->pmat, (MatFactorError *)&pc->failedreason));
49 PetscFunctionReturn(PETSC_SUCCESS);
50 }
51
PCApplyRichardson_SOR(PC pc,Vec b,Vec y,Vec w,PetscReal rtol,PetscReal abstol,PetscReal dtol,PetscInt its,PetscBool guesszero,PetscInt * outits,PCRichardsonConvergedReason * reason)52 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)
53 {
54 PC_SOR *jac = (PC_SOR *)pc->data;
55 MatSORType stype = jac->sym;
56
57 PetscFunctionBegin;
58 if (guesszero) stype = (MatSORType)(stype | SOR_ZERO_INITIAL_GUESS);
59 PetscCall(MatSOR(pc->pmat, b, jac->omega, stype, jac->fshift, its * jac->its, jac->lits, y));
60 PetscCall(MatFactorGetError(pc->pmat, (MatFactorError *)&pc->failedreason));
61 *outits = its;
62 *reason = PCRICHARDSON_CONVERGED_ITS;
63 PetscFunctionReturn(PETSC_SUCCESS);
64 }
65
PCSetFromOptions_SOR(PC pc,PetscOptionItems PetscOptionsObject)66 static PetscErrorCode PCSetFromOptions_SOR(PC pc, PetscOptionItems PetscOptionsObject)
67 {
68 PC_SOR *jac = (PC_SOR *)pc->data;
69 PetscBool flg;
70 PetscReal omega;
71
72 PetscFunctionBegin;
73 PetscOptionsHeadBegin(PetscOptionsObject, "(S)SOR options");
74 PetscCall(PetscOptionsReal("-pc_sor_omega", "relaxation factor (0 < omega < 2)", "PCSORSetOmega", jac->omega, &omega, &flg));
75 if (flg) PetscCall(PCSORSetOmega(pc, omega));
76 PetscCall(PetscOptionsReal("-pc_sor_diagonal_shift", "Add to the diagonal entries", "", jac->fshift, &jac->fshift, NULL));
77 PetscCall(PetscOptionsInt("-pc_sor_its", "number of inner SOR iterations", "PCSORSetIterations", jac->its, &jac->its, NULL));
78 PetscCall(PetscOptionsInt("-pc_sor_lits", "number of local inner SOR iterations", "PCSORSetIterations", jac->lits, &jac->lits, NULL));
79 PetscCall(PetscOptionsBoolGroupBegin("-pc_sor_symmetric", "SSOR, not SOR", "PCSORSetSymmetric", &flg));
80 if (flg) PetscCall(PCSORSetSymmetric(pc, SOR_SYMMETRIC_SWEEP));
81 PetscCall(PetscOptionsBoolGroup("-pc_sor_backward", "use backward sweep instead of forward", "PCSORSetSymmetric", &flg));
82 if (flg) PetscCall(PCSORSetSymmetric(pc, SOR_BACKWARD_SWEEP));
83 PetscCall(PetscOptionsBoolGroup("-pc_sor_forward", "use forward sweep", "PCSORSetSymmetric", &flg));
84 if (flg) PetscCall(PCSORSetSymmetric(pc, SOR_FORWARD_SWEEP));
85 PetscCall(PetscOptionsBoolGroup("-pc_sor_local_symmetric", "use SSOR separately on each processor", "PCSORSetSymmetric", &flg));
86 if (flg) PetscCall(PCSORSetSymmetric(pc, SOR_LOCAL_SYMMETRIC_SWEEP));
87 PetscCall(PetscOptionsBoolGroup("-pc_sor_local_backward", "use backward sweep locally", "PCSORSetSymmetric", &flg));
88 if (flg) PetscCall(PCSORSetSymmetric(pc, SOR_LOCAL_BACKWARD_SWEEP));
89 PetscCall(PetscOptionsBoolGroupEnd("-pc_sor_local_forward", "use forward sweep locally", "PCSORSetSymmetric", &flg));
90 if (flg) PetscCall(PCSORSetSymmetric(pc, SOR_LOCAL_FORWARD_SWEEP));
91 PetscOptionsHeadEnd();
92 PetscFunctionReturn(PETSC_SUCCESS);
93 }
94
PCView_SOR(PC pc,PetscViewer viewer)95 static PetscErrorCode PCView_SOR(PC pc, PetscViewer viewer)
96 {
97 PC_SOR *jac = (PC_SOR *)pc->data;
98 MatSORType sym = jac->sym;
99 const char *sortype;
100 PetscBool isascii;
101
102 PetscFunctionBegin;
103 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
104 if (isascii) {
105 if (sym & SOR_ZERO_INITIAL_GUESS) PetscCall(PetscViewerASCIIPrintf(viewer, " zero initial guess\n"));
106 if (sym == SOR_APPLY_UPPER) sortype = "apply_upper";
107 else if (sym == SOR_APPLY_LOWER) sortype = "apply_lower";
108 else if (sym & SOR_EISENSTAT) sortype = "Eisenstat";
109 else if ((sym & SOR_SYMMETRIC_SWEEP) == SOR_SYMMETRIC_SWEEP) sortype = "symmetric";
110 else if (sym & SOR_BACKWARD_SWEEP) sortype = "backward";
111 else if (sym & SOR_FORWARD_SWEEP) sortype = "forward";
112 else if ((sym & SOR_LOCAL_SYMMETRIC_SWEEP) == SOR_LOCAL_SYMMETRIC_SWEEP) 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 PetscCall(PetscViewerASCIIPrintf(viewer, " type = %s, iterations = %" PetscInt_FMT ", local iterations = %" PetscInt_FMT ", omega = %g\n", sortype, jac->its, jac->lits, (double)jac->omega));
117 }
118 PetscFunctionReturn(PETSC_SUCCESS);
119 }
120
PCSORSetSymmetric_SOR(PC pc,MatSORType flag)121 static PetscErrorCode PCSORSetSymmetric_SOR(PC pc, MatSORType flag)
122 {
123 PC_SOR *jac = (PC_SOR *)pc->data;
124
125 PetscFunctionBegin;
126 jac->sym = flag;
127 PetscFunctionReturn(PETSC_SUCCESS);
128 }
129
PCSORSetOmega_SOR(PC pc,PetscReal omega)130 static PetscErrorCode PCSORSetOmega_SOR(PC pc, PetscReal omega)
131 {
132 PC_SOR *jac = (PC_SOR *)pc->data;
133
134 PetscFunctionBegin;
135 PetscCheck(omega > 0.0 && omega < 2.0, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "Relaxation out of range");
136 jac->omega = omega;
137 PetscFunctionReturn(PETSC_SUCCESS);
138 }
139
PCSORSetIterations_SOR(PC pc,PetscInt its,PetscInt lits)140 static PetscErrorCode PCSORSetIterations_SOR(PC pc, PetscInt its, PetscInt lits)
141 {
142 PC_SOR *jac = (PC_SOR *)pc->data;
143
144 PetscFunctionBegin;
145 jac->its = its;
146 jac->lits = lits;
147 PetscFunctionReturn(PETSC_SUCCESS);
148 }
149
PCSORGetSymmetric_SOR(PC pc,MatSORType * flag)150 static PetscErrorCode PCSORGetSymmetric_SOR(PC pc, MatSORType *flag)
151 {
152 PC_SOR *jac = (PC_SOR *)pc->data;
153
154 PetscFunctionBegin;
155 *flag = jac->sym;
156 PetscFunctionReturn(PETSC_SUCCESS);
157 }
158
PCSORGetOmega_SOR(PC pc,PetscReal * omega)159 static PetscErrorCode PCSORGetOmega_SOR(PC pc, PetscReal *omega)
160 {
161 PC_SOR *jac = (PC_SOR *)pc->data;
162
163 PetscFunctionBegin;
164 *omega = jac->omega;
165 PetscFunctionReturn(PETSC_SUCCESS);
166 }
167
PCSORGetIterations_SOR(PC pc,PetscInt * its,PetscInt * lits)168 static PetscErrorCode PCSORGetIterations_SOR(PC pc, PetscInt *its, PetscInt *lits)
169 {
170 PC_SOR *jac = (PC_SOR *)pc->data;
171
172 PetscFunctionBegin;
173 if (its) *its = jac->its;
174 if (lits) *lits = jac->lits;
175 PetscFunctionReturn(PETSC_SUCCESS);
176 }
177
178 /*@
179 PCSORGetSymmetric - Gets the form the SOR preconditioner is using; backward, or forward relaxation. The local variants perform SOR on
180 each processor. By default forward relaxation is used.
181
182 Logically Collective
183
184 Input Parameter:
185 . pc - the preconditioner context
186
187 Output Parameter:
188 . flag - one of the following
189 .vb
190 SOR_FORWARD_SWEEP
191 SOR_BACKWARD_SWEEP
192 SOR_SYMMETRIC_SWEEP
193 SOR_LOCAL_FORWARD_SWEEP
194 SOR_LOCAL_BACKWARD_SWEEP
195 SOR_LOCAL_SYMMETRIC_SWEEP
196 .ve
197
198 Options Database Keys:
199 + -pc_sor_symmetric - Activates symmetric version
200 . -pc_sor_backward - Activates backward version
201 . -pc_sor_local_forward - Activates local forward version
202 . -pc_sor_local_symmetric - Activates local symmetric version
203 - -pc_sor_local_backward - Activates local backward version
204
205 Note:
206 To use the Eisenstat trick with SSOR, employ the `PCEISENSTAT` preconditioner,
207 which can be chosen with the option
208 . -pc_type eisenstat - Activates Eisenstat trick
209
210 Level: intermediate
211
212 .seealso: [](ch_ksp), `PCSOR`, `PCEisenstatSetOmega()`, `PCSORSetIterations()`, `PCSORSetOmega()`, `PCSORSetSymmetric()`
213 @*/
PCSORGetSymmetric(PC pc,MatSORType * flag)214 PetscErrorCode PCSORGetSymmetric(PC pc, MatSORType *flag)
215 {
216 PetscFunctionBegin;
217 PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
218 PetscUseMethod(pc, "PCSORGetSymmetric_C", (PC, MatSORType *), (pc, flag));
219 PetscFunctionReturn(PETSC_SUCCESS);
220 }
221
222 /*@
223 PCSORGetOmega - Gets the SOR relaxation coefficient, omega
224 (where omega = 1.0 by default).
225
226 Logically Collective
227
228 Input Parameter:
229 . pc - the preconditioner context
230
231 Output Parameter:
232 . omega - relaxation coefficient (0 < omega < 2).
233
234 Options Database Key:
235 . -pc_sor_omega <omega> - Sets omega
236
237 Level: intermediate
238
239 .seealso: [](ch_ksp), `PCSOR`, `PCSORSetSymmetric()`, `PCSORSetIterations()`, `PCEisenstatSetOmega()`, `PCSORSetOmega()`
240 @*/
PCSORGetOmega(PC pc,PetscReal * omega)241 PetscErrorCode PCSORGetOmega(PC pc, PetscReal *omega)
242 {
243 PetscFunctionBegin;
244 PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
245 PetscUseMethod(pc, "PCSORGetOmega_C", (PC, PetscReal *), (pc, omega));
246 PetscFunctionReturn(PETSC_SUCCESS);
247 }
248
249 /*@
250 PCSORGetIterations - Gets the number of inner iterations to
251 be used by the SOR preconditioner. The default is 1.
252
253 Logically Collective
254
255 Input Parameter:
256 . pc - the preconditioner context
257
258 Output Parameters:
259 + lits - number of local iterations, smoothings over just variables on processor
260 - its - number of parallel iterations to use; each parallel iteration has lits local iterations
261
262 Options Database Keys:
263 + -pc_sor_its <its> - Sets number of iterations
264 - -pc_sor_lits <lits> - Sets number of local iterations
265
266 Level: intermediate
267
268 Note:
269 When run on one processor the number of smoothings is lits*its
270
271 .seealso: [](ch_ksp), `PCSOR`, `PCSORSetOmega()`, `PCSORSetSymmetric()`, `PCSORSetIterations()`
272 @*/
PCSORGetIterations(PC pc,PetscInt * its,PetscInt * lits)273 PetscErrorCode PCSORGetIterations(PC pc, PetscInt *its, PetscInt *lits)
274 {
275 PetscFunctionBegin;
276 PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
277 PetscUseMethod(pc, "PCSORGetIterations_C", (PC, PetscInt *, PetscInt *), (pc, its, lits));
278 PetscFunctionReturn(PETSC_SUCCESS);
279 }
280
281 /*@
282 PCSORSetSymmetric - Sets the SOR preconditioner to use symmetric (SSOR),
283 backward, or forward relaxation. The local variants perform SOR on
284 each processor. By default forward relaxation is used.
285
286 Logically Collective
287
288 Input Parameters:
289 + pc - the preconditioner context
290 - flag - one of the following
291 .vb
292 SOR_FORWARD_SWEEP
293 SOR_BACKWARD_SWEEP
294 SOR_SYMMETRIC_SWEEP
295 SOR_LOCAL_FORWARD_SWEEP
296 SOR_LOCAL_BACKWARD_SWEEP
297 SOR_LOCAL_SYMMETRIC_SWEEP
298 .ve
299
300 Options Database Keys:
301 + -pc_sor_symmetric - Activates symmetric version
302 . -pc_sor_backward - Activates backward version
303 . -pc_sor_local_forward - Activates local forward version
304 . -pc_sor_local_symmetric - Activates local symmetric version
305 - -pc_sor_local_backward - Activates local backward version
306
307 Note:
308 To use the Eisenstat trick with SSOR, employ the PCEISENSTAT preconditioner,
309 which can be chosen with the option
310 . -pc_type eisenstat - Activates Eisenstat trick
311
312 Level: intermediate
313
314 .seealso: [](ch_ksp), `PCSOR`, `PCEisenstatSetOmega()`, `PCSORSetIterations()`, `PCSORSetOmega()`
315 @*/
PCSORSetSymmetric(PC pc,MatSORType flag)316 PetscErrorCode PCSORSetSymmetric(PC pc, MatSORType flag)
317 {
318 PetscFunctionBegin;
319 PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
320 PetscValidLogicalCollectiveEnum(pc, flag, 2);
321 PetscTryMethod(pc, "PCSORSetSymmetric_C", (PC, MatSORType), (pc, flag));
322 PetscFunctionReturn(PETSC_SUCCESS);
323 }
324
325 /*@
326 PCSORSetOmega - Sets the SOR relaxation coefficient, omega
327 (where omega = 1.0 by default).
328
329 Logically Collective
330
331 Input Parameters:
332 + pc - the preconditioner context
333 - omega - relaxation coefficient (0 < omega < 2).
334
335 Options Database Key:
336 . -pc_sor_omega <omega> - Sets omega
337
338 Level: intermediate
339
340 Note:
341 If omega != 1, you will need to set the `MAT_USE_INODE`S option to `PETSC_FALSE` on the matrix.
342
343 .seealso: [](ch_ksp), `PCSOR`, `PCSORSetSymmetric()`, `PCSORSetIterations()`, `PCEisenstatSetOmega()`, `MatSetOption()`
344 @*/
PCSORSetOmega(PC pc,PetscReal omega)345 PetscErrorCode PCSORSetOmega(PC pc, PetscReal omega)
346 {
347 PetscFunctionBegin;
348 PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
349 PetscValidLogicalCollectiveReal(pc, omega, 2);
350 PetscTryMethod(pc, "PCSORSetOmega_C", (PC, PetscReal), (pc, omega));
351 PetscFunctionReturn(PETSC_SUCCESS);
352 }
353
354 /*@
355 PCSORSetIterations - Sets the number of inner iterations to
356 be used by the SOR preconditioner. The default is 1.
357
358 Logically Collective
359
360 Input Parameters:
361 + pc - the preconditioner context
362 . lits - number of local iterations, smoothings over just variables on processor
363 - its - number of parallel iterations to use; each parallel iteration has lits local iterations
364
365 Options Database Keys:
366 + -pc_sor_its <its> - Sets number of iterations
367 - -pc_sor_lits <lits> - Sets number of local iterations
368
369 Level: intermediate
370
371 Note:
372 When run on one processor the number of smoothings is lits*its
373
374 .seealso: [](ch_ksp), `PCSOR`, `PCSORSetOmega()`, `PCSORSetSymmetric()`
375 @*/
PCSORSetIterations(PC pc,PetscInt its,PetscInt lits)376 PetscErrorCode PCSORSetIterations(PC pc, PetscInt its, PetscInt lits)
377 {
378 PetscFunctionBegin;
379 PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
380 PetscValidLogicalCollectiveInt(pc, its, 2);
381 PetscTryMethod(pc, "PCSORSetIterations_C", (PC, PetscInt, PetscInt), (pc, its, lits));
382 PetscFunctionReturn(PETSC_SUCCESS);
383 }
384
385 /*MC
386 PCSOR - (S)SOR (successive over relaxation, Gauss-Seidel) preconditioning
387
388 Options Database Keys:
389 + -pc_sor_symmetric - Activates symmetric version
390 . -pc_sor_backward - Activates backward version
391 . -pc_sor_forward - Activates forward version
392 . -pc_sor_local_forward - Activates local forward version
393 . -pc_sor_local_symmetric - Activates local symmetric version (default version)
394 . -pc_sor_local_backward - Activates local backward version
395 . -pc_sor_omega <omega> - Sets omega
396 . -pc_sor_diagonal_shift <shift> - shift the diagonal entries; useful if the matrix has zeros on the diagonal
397 . -pc_sor_its <its> - Sets number of iterations (default 1)
398 - -pc_sor_lits <lits> - Sets number of local iterations (default 1)
399
400 Level: beginner
401
402 Notes:
403 Only implemented for the `MATAIJ` and `MATSEQBAIJ` matrix formats.
404
405 Not a true parallel SOR, in parallel this implementation corresponds to block
406 Jacobi with SOR on each block.
407
408 For `MATAIJ` matrices if a diagonal entry is zero (and the diagonal shift is zero) then by default the inverse of that
409 zero will be used and hence the `KSPSolve()` will terminate with `KSP_DIVERGED_NANORINF`. If the option
410 `KSPSetErrorIfNotConverged()` or -ksp_error_if_not_converged the code will terminate as soon as it detects the
411 zero pivot.
412
413 For `MATSEQBAIJ` matrices this implements point-block SOR, but the omega, its, lits options are not supported.
414
415 For `MATSEQBAIJ` the diagonal blocks are inverted using dense LU with partial pivoting. If a zero pivot is detected
416 the computation is stopped with an error
417
418 If used with `KSPRICHARDSON` and no monitors the convergence test is skipped to improve speed, thus it always iterates
419 the maximum number of iterations you've selected for `KSP`. It is usually used in this mode as a smoother for multigrid.
420
421 If omega != 1, you will need to set the `MAT_USE_INODES` option to `PETSC_FALSE` on the matrix.
422
423 .seealso: [](ch_ksp), `PCCreate()`, `PCSetType()`, `PCType`, `PC`, `PCJACOBI`,
424 `PCSORSetIterations()`, `PCSORSetSymmetric()`, `PCSORSetOmega()`, `PCEISENSTAT`, `MatSetOption()`
425 M*/
426
PCCreate_SOR(PC pc)427 PETSC_EXTERN PetscErrorCode PCCreate_SOR(PC pc)
428 {
429 PC_SOR *jac;
430
431 PetscFunctionBegin;
432 PetscCall(PetscNew(&jac));
433
434 pc->ops->apply = PCApply_SOR;
435 pc->ops->applytranspose = PCApplyTranspose_SOR;
436 pc->ops->applyrichardson = PCApplyRichardson_SOR;
437 pc->ops->setfromoptions = PCSetFromOptions_SOR;
438 pc->ops->setup = NULL;
439 pc->ops->view = PCView_SOR;
440 pc->ops->destroy = PCDestroy_SOR;
441 pc->data = (void *)jac;
442 jac->sym = SOR_LOCAL_SYMMETRIC_SWEEP;
443 jac->omega = 1.0;
444 jac->fshift = 0.0;
445 jac->its = 1;
446 jac->lits = 1;
447
448 PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSORSetSymmetric_C", PCSORSetSymmetric_SOR));
449 PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSORSetOmega_C", PCSORSetOmega_SOR));
450 PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSORSetIterations_C", PCSORSetIterations_SOR));
451 PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSORGetSymmetric_C", PCSORGetSymmetric_SOR));
452 PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSORGetOmega_C", PCSORGetOmega_SOR));
453 PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCSORGetIterations_C", PCSORGetIterations_SOR));
454 PetscFunctionReturn(PETSC_SUCCESS);
455 }
456