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