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