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