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