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