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