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 34 PetscFunctionBegin; 35 ierr = MatSOR(pc->pmat,x,jac->omega,(MatSORType)flag,jac->fshift,jac->its,jac->lits,y);CHKERRQ(ierr); 36 PetscFunctionReturn(0); 37 } 38 39 #undef __FUNCT__ 40 #define __FUNCT__ "PCApplyRichardson_SOR" 41 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) 42 { 43 PC_SOR *jac = (PC_SOR*)pc->data; 44 PetscErrorCode ierr; 45 MatSORType stype = jac->sym; 46 47 PetscFunctionBegin; 48 ierr = PetscInfo1(pc,"Warning, convergence critera ignored, using %D iterations\n",its);CHKERRQ(ierr); 49 if (guesszero) stype = (MatSORType) (stype | SOR_ZERO_INITIAL_GUESS); 50 ierr = MatSOR(pc->pmat,b,jac->omega,stype,jac->fshift,its*jac->its,jac->lits,y);CHKERRQ(ierr); 51 *outits = its; 52 *reason = PCRICHARDSON_CONVERGED_ITS; 53 PetscFunctionReturn(0); 54 } 55 56 #undef __FUNCT__ 57 #define __FUNCT__ "PCSetFromOptions_SOR" 58 PetscErrorCode PCSetFromOptions_SOR(PC pc) 59 { 60 PC_SOR *jac = (PC_SOR*)pc->data; 61 PetscErrorCode ierr; 62 PetscBool flg; 63 64 PetscFunctionBegin; 65 ierr = PetscOptionsHead("(S)SOR options");CHKERRQ(ierr); 66 ierr = PetscOptionsReal("-pc_sor_omega","relaxation factor (0 < omega < 2)","PCSORSetOmega",jac->omega,&jac->omega,0);CHKERRQ(ierr); 67 ierr = PetscOptionsReal("-pc_sor_diagonal_shift","Add to the diagonal entries","",jac->fshift,&jac->fshift,0);CHKERRQ(ierr); 68 ierr = PetscOptionsInt("-pc_sor_its","number of inner SOR iterations","PCSORSetIterations",jac->its,&jac->its,0);CHKERRQ(ierr); 69 ierr = PetscOptionsInt("-pc_sor_lits","number of local inner SOR iterations","PCSORSetIterations",jac->lits,&jac->lits,0);CHKERRQ(ierr); 70 ierr = PetscOptionsBoolGroupBegin("-pc_sor_symmetric","SSOR, not SOR","PCSORSetSymmetric",&flg);CHKERRQ(ierr); 71 if (flg) {ierr = PCSORSetSymmetric(pc,SOR_SYMMETRIC_SWEEP);CHKERRQ(ierr);} 72 ierr = PetscOptionsBoolGroup("-pc_sor_backward","use backward sweep instead of forward","PCSORSetSymmetric",&flg);CHKERRQ(ierr); 73 if (flg) {ierr = PCSORSetSymmetric(pc,SOR_BACKWARD_SWEEP);CHKERRQ(ierr);} 74 ierr = PetscOptionsBoolGroup("-pc_sor_forward","use forward sweep","PCSORSetSymmetric",&flg);CHKERRQ(ierr); 75 if (flg) {ierr = PCSORSetSymmetric(pc,SOR_FORWARD_SWEEP);CHKERRQ(ierr);} 76 ierr = PetscOptionsBoolGroup("-pc_sor_local_symmetric","use SSOR separately on each processor","PCSORSetSymmetric",&flg);CHKERRQ(ierr); 77 if (flg) {ierr = PCSORSetSymmetric(pc,SOR_LOCAL_SYMMETRIC_SWEEP);CHKERRQ(ierr);} 78 ierr = PetscOptionsBoolGroup("-pc_sor_local_backward","use backward sweep locally","PCSORSetSymmetric",&flg);CHKERRQ(ierr); 79 if (flg) {ierr = PCSORSetSymmetric(pc,SOR_LOCAL_BACKWARD_SWEEP);CHKERRQ(ierr);} 80 ierr = PetscOptionsBoolGroupEnd("-pc_sor_local_forward","use forward sweep locally","PCSORSetSymmetric",&flg);CHKERRQ(ierr); 81 if (flg) {ierr = PCSORSetSymmetric(pc,SOR_LOCAL_FORWARD_SWEEP);CHKERRQ(ierr);} 82 ierr = PetscOptionsTail();CHKERRQ(ierr); 83 PetscFunctionReturn(0); 84 } 85 86 #undef __FUNCT__ 87 #define __FUNCT__ "PCView_SOR" 88 PetscErrorCode PCView_SOR(PC pc,PetscViewer viewer) 89 { 90 PC_SOR *jac = (PC_SOR*)pc->data; 91 MatSORType sym = jac->sym; 92 const char *sortype; 93 PetscErrorCode ierr; 94 PetscBool iascii; 95 96 PetscFunctionBegin; 97 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 98 if (iascii) { 99 if (sym & SOR_ZERO_INITIAL_GUESS) {ierr = PetscViewerASCIIPrintf(viewer," SOR: zero initial guess\n");CHKERRQ(ierr);} 100 if (sym == SOR_APPLY_UPPER) sortype = "apply_upper"; 101 else if (sym == SOR_APPLY_LOWER) sortype = "apply_lower"; 102 else if (sym & SOR_EISENSTAT) sortype = "Eisenstat"; 103 else if ((sym & SOR_SYMMETRIC_SWEEP) == SOR_SYMMETRIC_SWEEP) sortype = "symmetric"; 104 else if (sym & SOR_BACKWARD_SWEEP) sortype = "backward"; 105 else if (sym & SOR_FORWARD_SWEEP) sortype = "forward"; 106 else if ((sym & SOR_LOCAL_SYMMETRIC_SWEEP) == SOR_LOCAL_SYMMETRIC_SWEEP) sortype = "local_symmetric"; 107 else if (sym & SOR_LOCAL_FORWARD_SWEEP) sortype = "local_forward"; 108 else if (sym & SOR_LOCAL_BACKWARD_SWEEP) sortype = "local_backward"; 109 else sortype = "unknown"; 110 ierr = PetscViewerASCIIPrintf(viewer," SOR: type = %s, iterations = %D, local iterations = %D, omega = %g\n",sortype,jac->its,jac->lits,(double)jac->omega);CHKERRQ(ierr); 111 } 112 PetscFunctionReturn(0); 113 } 114 115 116 /* ------------------------------------------------------------------------------*/ 117 #undef __FUNCT__ 118 #define __FUNCT__ "PCSORSetSymmetric_SOR" 119 static PetscErrorCode PCSORSetSymmetric_SOR(PC pc,MatSORType flag) 120 { 121 PC_SOR *jac; 122 123 PetscFunctionBegin; 124 jac = (PC_SOR*)pc->data; 125 jac->sym = flag; 126 PetscFunctionReturn(0); 127 } 128 129 #undef __FUNCT__ 130 #define __FUNCT__ "PCSORSetOmega_SOR" 131 static PetscErrorCode PCSORSetOmega_SOR(PC pc,PetscReal omega) 132 { 133 PC_SOR *jac; 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 = (PC_SOR*)pc->data; 138 jac->omega = omega; 139 PetscFunctionReturn(0); 140 } 141 142 #undef __FUNCT__ 143 #define __FUNCT__ "PCSORSetIterations_SOR" 144 static PetscErrorCode PCSORSetIterations_SOR(PC pc,PetscInt its,PetscInt lits) 145 { 146 PC_SOR *jac; 147 148 PetscFunctionBegin; 149 jac = (PC_SOR*)pc->data; 150 jac->its = its; 151 jac->lits = lits; 152 PetscFunctionReturn(0); 153 } 154 155 /* ------------------------------------------------------------------------------*/ 156 #undef __FUNCT__ 157 #define __FUNCT__ "PCSORSetSymmetric" 158 /*@ 159 PCSORSetSymmetric - Sets the SOR preconditioner to use symmetric (SSOR), 160 backward, or forward relaxation. The local variants perform SOR on 161 each processor. By default forward relaxation is used. 162 163 Logically Collective on PC 164 165 Input Parameters: 166 + pc - the preconditioner context 167 - flag - one of the following 168 .vb 169 SOR_FORWARD_SWEEP 170 SOR_BACKWARD_SWEEP 171 SOR_SYMMETRIC_SWEEP 172 SOR_LOCAL_FORWARD_SWEEP 173 SOR_LOCAL_BACKWARD_SWEEP 174 SOR_LOCAL_SYMMETRIC_SWEEP 175 .ve 176 177 Options Database Keys: 178 + -pc_sor_symmetric - Activates symmetric version 179 . -pc_sor_backward - Activates backward version 180 . -pc_sor_local_forward - Activates local forward version 181 . -pc_sor_local_symmetric - Activates local symmetric version 182 - -pc_sor_local_backward - Activates local backward version 183 184 Notes: 185 To use the Eisenstat trick with SSOR, employ the PCEISENSTAT preconditioner, 186 which can be chosen with the option 187 . -pc_type eisenstat - Activates Eisenstat trick 188 189 Level: intermediate 190 191 .keywords: PC, SOR, SSOR, set, relaxation, sweep, forward, backward, symmetric 192 193 .seealso: PCEisenstatSetOmega(), PCSORSetIterations(), PCSORSetOmega() 194 @*/ 195 PetscErrorCode PCSORSetSymmetric(PC pc,MatSORType flag) 196 { 197 PetscErrorCode ierr; 198 199 PetscFunctionBegin; 200 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 201 PetscValidLogicalCollectiveEnum(pc,flag,2); 202 ierr = PetscTryMethod(pc,"PCSORSetSymmetric_C",(PC,MatSORType),(pc,flag));CHKERRQ(ierr); 203 PetscFunctionReturn(0); 204 } 205 206 #undef __FUNCT__ 207 #define __FUNCT__ "PCSORSetOmega" 208 /*@ 209 PCSORSetOmega - Sets the SOR relaxation coefficient, omega 210 (where omega = 1.0 by default). 211 212 Logically Collective on PC 213 214 Input Parameters: 215 + pc - the preconditioner context 216 - omega - relaxation coefficient (0 < omega < 2). 217 218 Options Database Key: 219 . -pc_sor_omega <omega> - Sets omega 220 221 Level: intermediate 222 223 .keywords: PC, SOR, SSOR, set, relaxation, omega 224 225 .seealso: PCSORSetSymmetric(), PCSORSetIterations(), PCEisenstatSetOmega() 226 @*/ 227 PetscErrorCode PCSORSetOmega(PC pc,PetscReal omega) 228 { 229 PetscErrorCode ierr; 230 231 PetscFunctionBegin; 232 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 233 PetscValidLogicalCollectiveReal(pc,omega,2); 234 ierr = PetscTryMethod(pc,"PCSORSetOmega_C",(PC,PetscReal),(pc,omega));CHKERRQ(ierr); 235 PetscFunctionReturn(0); 236 } 237 238 #undef __FUNCT__ 239 #define __FUNCT__ "PCSORSetIterations" 240 /*@ 241 PCSORSetIterations - Sets the number of inner iterations to 242 be used by the SOR preconditioner. The default is 1. 243 244 Logically Collective on PC 245 246 Input Parameters: 247 + pc - the preconditioner context 248 . lits - number of local iterations, smoothings over just variables on processor 249 - its - number of parallel iterations to use; each parallel iteration has lits local iterations 250 251 Options Database Key: 252 + -pc_sor_its <its> - Sets number of iterations 253 - -pc_sor_lits <lits> - Sets number of local iterations 254 255 Level: intermediate 256 257 Notes: When run on one processor the number of smoothings is lits*its 258 259 .keywords: PC, SOR, SSOR, set, iterations 260 261 .seealso: PCSORSetOmega(), PCSORSetSymmetric() 262 @*/ 263 PetscErrorCode PCSORSetIterations(PC pc,PetscInt its,PetscInt lits) 264 { 265 PetscErrorCode ierr; 266 267 PetscFunctionBegin; 268 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 269 PetscValidLogicalCollectiveInt(pc,its,2); 270 ierr = PetscTryMethod(pc,"PCSORSetIterations_C",(PC,PetscInt,PetscInt),(pc,its,lits));CHKERRQ(ierr); 271 PetscFunctionReturn(0); 272 } 273 274 /*MC 275 PCSOR - (S)SOR (successive over relaxation, Gauss-Seidel) preconditioning 276 277 Options Database Keys: 278 + -pc_sor_symmetric - Activates symmetric version 279 . -pc_sor_backward - Activates backward version 280 . -pc_sor_forward - Activates forward version 281 . -pc_sor_local_forward - Activates local forward version 282 . -pc_sor_local_symmetric - Activates local symmetric version (default version) 283 . -pc_sor_local_backward - Activates local backward version 284 . -pc_sor_omega <omega> - Sets omega 285 . -pc_sor_its <its> - Sets number of iterations (default 1) 286 - -pc_sor_lits <lits> - Sets number of local iterations (default 1) 287 288 Level: beginner 289 290 Concepts: SOR, preconditioners, Gauss-Seidel 291 292 Notes: Only implemented for the AIJ and SeqBAIJ matrix formats. 293 Not a true parallel SOR, in parallel this implementation corresponds to block 294 Jacobi with SOR on each block. 295 296 For SeqBAIJ matrices this implements point-block SOR, but the omega, its, lits options are not supported. 297 298 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 299 PCSORSetIterations(), PCSORSetSymmetric(), PCSORSetOmega(), PCEISENSTAT 300 M*/ 301 302 #undef __FUNCT__ 303 #define __FUNCT__ "PCCreate_SOR" 304 PETSC_EXTERN PetscErrorCode PCCreate_SOR(PC pc) 305 { 306 PetscErrorCode ierr; 307 PC_SOR *jac; 308 309 PetscFunctionBegin; 310 ierr = PetscNewLog(pc,&jac);CHKERRQ(ierr); 311 312 pc->ops->apply = PCApply_SOR; 313 pc->ops->applyrichardson = PCApplyRichardson_SOR; 314 pc->ops->setfromoptions = PCSetFromOptions_SOR; 315 pc->ops->setup = 0; 316 pc->ops->view = PCView_SOR; 317 pc->ops->destroy = PCDestroy_SOR; 318 pc->data = (void*)jac; 319 jac->sym = SOR_LOCAL_SYMMETRIC_SWEEP; 320 jac->omega = 1.0; 321 jac->fshift = 0.0; 322 jac->its = 1; 323 jac->lits = 1; 324 325 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCSORSetSymmetric_C",PCSORSetSymmetric_SOR);CHKERRQ(ierr); 326 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCSORSetOmega_C",PCSORSetOmega_SOR);CHKERRQ(ierr); 327 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCSORSetIterations_C",PCSORSetIterations_SOR);CHKERRQ(ierr); 328 PetscFunctionReturn(0); 329 } 330 331 332 333 334 335