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