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