1 2 /* -------------------------------------------------------------------- 3 4 This file implements a Deflation preconditioner in PETSc as part of PC. 5 You can use this as a starting point for implementing your own 6 preconditioner that is not provided with PETSc. (You might also consider 7 just using PCSHELL) 8 9 The following basic routines are required for each preconditioner. 10 PCCreate_XXX() - Creates a preconditioner context 11 PCSetFromOptions_XXX() - Sets runtime options 12 PCApply_XXX() - Applies the preconditioner 13 PCDestroy_XXX() - Destroys the preconditioner context 14 where the suffix "_XXX" denotes a particular implementation, in 15 this case we use _Deflation (e.g., PCCreate_Deflation, PCApply_Deflation). 16 These routines are actually called via the common user interface 17 routines PCCreate(), PCSetFromOptions(), PCApply(), and PCDestroy(), 18 so the application code interface remains identical for all 19 preconditioners. 20 21 Another key routine is: 22 PCSetUp_XXX() - Prepares for the use of a preconditioner 23 by setting data structures and options. The interface routine PCSetUp() 24 is not usually called directly by the user, but instead is called by 25 PCApply() if necessary. 26 27 Additional basic routines are: 28 PCView_XXX() - Prints details of runtime options that 29 have actually been used. 30 These are called by application codes via the interface routines 31 PCView(). 32 33 The various types of solvers (preconditioners, Krylov subspace methods, 34 nonlinear solvers, timesteppers) are all organized similarly, so the 35 above description applies to these categories also. One exception is 36 that the analogues of PCApply() for these components are KSPSolve(), 37 SNESSolve(), and TSSolve(). 38 39 Additional optional functionality unique to preconditioners is left and 40 right symmetric preconditioner application via PCApplySymmetricLeft() 41 and PCApplySymmetricRight(). The Deflation implementation is 42 PCApplySymmetricLeftOrRight_Deflation(). 43 44 -------------------------------------------------------------------- */ 45 46 /* 47 Include files needed for the Deflation preconditioner: 48 pcimpl.h - private include file intended for use by all preconditioners 49 */ 50 51 #include <petsc/private/pcimpl.h> /*I "petscpc.h" I*/ 52 53 const char *const PCDeflationTypes[] = {"INIT","PRE","POST","PCDeflationType","PC_DEFLATION_",0}; 54 55 /* 56 Private context (data structure) for the deflation preconditioner. 57 */ 58 typedef struct { 59 PetscBool init; /* do only init step - error correction of direction is omitted */ 60 PetscBool pre; /* start with x0 being the solution in the deflation space */ 61 PetscBool correct; /* add CP (Qr) correction to descent direction */ 62 PetscBool truenorm; 63 PetscBool adaptiveconv; 64 PetscReal adaptiveconst; 65 PetscInt reductionfact; 66 Mat W,Wt,AW,WtAW; /* deflation space, coarse problem mats */ 67 KSP WtAWinv; /* deflation coarse problem */ 68 Vec *work; 69 70 PCDEFLATIONSpaceType spacetype; 71 PetscInt spacesize; 72 PetscInt nestedlvl; 73 PetscInt maxnestedlvl; 74 PetscBool extendsp; 75 } PC_Deflation; 76 77 static PetscErrorCode PCDeflationSetType_Deflation(PC pc,PCDeflationType type) 78 { 79 PC_Deflation *def = (PC_Deflation*)pc->data; 80 81 PetscFunctionBegin; 82 def->init = PETSC_FALSE; 83 def->pre = PETSC_FALSE; 84 if (type == PC_DEFLATION_INIT) { 85 def->init = PETSC_TRUE; 86 def->pre = PETSC_TRUE; 87 } else if (type == PC_DEFLATION_PRE) { 88 def->pre = PETSC_TRUE; 89 } 90 PetscFunctionReturn(0); 91 } 92 93 /*@ 94 PCDeflationSetType - Causes the deflation preconditioner to use only a special 95 initial gues or pre/post solve solution update 96 97 Logically Collective on PC 98 99 Input Parameters: 100 + pc - the preconditioner context 101 - type - PC_DEFLATION_PRE, PC_DEFLATION_INIT, PC_DEFLATION_POST 102 103 Options Database Key: 104 . -pc_deflation_type <pre,init,post> 105 106 Level: intermediate 107 108 Concepts: Deflation preconditioner 109 110 .seealso: PCDeflationGetType() 111 @*/ 112 PetscErrorCode PCDeflationSetType(PC pc,PCDeflationType type) 113 { 114 PetscErrorCode ierr; 115 116 PetscFunctionBegin; 117 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 118 ierr = PetscTryMethod(pc,"PCDeflationSetType_C",(PC,PCDeflationType),(pc,type));CHKERRQ(ierr); 119 PetscFunctionReturn(0); 120 } 121 122 static PetscErrorCode PCDeflationGetType_Deflation(PC pc,PCDeflationType *type) 123 { 124 PC_Deflation *def = (PC_Deflation*)pc->data; 125 126 PetscFunctionBegin; 127 if (def->init) { 128 *type = PC_DEFLATION_INIT; 129 } else if (def->pre) { 130 *type = PC_DEFLATION_PRE; 131 } else { 132 *type = PC_DEFLATION_POST; 133 } 134 PetscFunctionReturn(0); 135 } 136 137 /*@ 138 PCDeflationGetType - Gets how the diagonal matrix is produced for the preconditioner 139 140 Not Collective on PC 141 142 Input Parameter: 143 . pc - the preconditioner context 144 145 Output Parameter: 146 - type - PC_DEFLATION_PRE, PC_DEFLATION_INIT, PC_DEFLATION_POST 147 148 Level: intermediate 149 150 Concepts: Deflation preconditioner 151 152 .seealso: PCDeflationSetType() 153 @*/ 154 PetscErrorCode PCDeflationGetType(PC pc,PCDeflationType *type) 155 { 156 PetscErrorCode ierr; 157 158 PetscFunctionBegin; 159 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 160 ierr = PetscUseMethod(pc,"PCDeflationGetType_C",(PC,PCDeflationType*),(pc,type));CHKERRQ(ierr); 161 PetscFunctionReturn(0); 162 } 163 164 static PetscErrorCode PCDeflationSetSpace_Deflation(PC pc,Mat W,PetscBool transpose) 165 { 166 PC_Deflation *def = (PC_Deflation*)pc->data; 167 PetscErrorCode ierr; 168 169 PetscFunctionBegin; 170 if (transpose) { 171 def->Wt = W; 172 def->W = NULL; 173 } else { 174 def->W = W; 175 } 176 ierr = PetscObjectReference((PetscObject)W);CHKERRQ(ierr); 177 PetscFunctionReturn(0); 178 } 179 180 /* TODO create PCDeflationSetSpaceTranspose? */ 181 /*@ 182 PCDeflationSetSpace - Set deflation space matrix (or its transpose). 183 184 Logically Collective on PC 185 186 Input Parameters: 187 + pc - the preconditioner context 188 . W - deflation matrix 189 - tranpose - indicates that W is an explicit transpose of the deflation matrix 190 191 Level: intermediate 192 193 .seealso: PCDEFLATION 194 @*/ 195 PetscErrorCode PCDeflationSetSpace(PC pc,Mat W,PetscBool transpose) 196 { 197 PetscErrorCode ierr; 198 199 PetscFunctionBegin; 200 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 201 PetscValidHeaderSpecific(W,MAT_CLASSID,2); 202 PetscValidLogicalCollectiveBool(pc,transpose,3); 203 ierr = PetscTryMethod(pc,"PCDeflationSetSpace_C",(PC,Mat,PetscBool),(pc,W,transpose));CHKERRQ(ierr); 204 PetscFunctionReturn(0); 205 } 206 207 208 /* -------------------------------------------------------------------------- */ 209 /* 210 PCSetUp_Deflation - Prepares for the use of the Deflation preconditioner 211 by setting data structures and options. 212 213 Input Parameter: 214 . pc - the preconditioner context 215 216 Application Interface Routine: PCSetUp() 217 218 Notes: 219 The interface routine PCSetUp() is not usually called directly by 220 the user, but instead is called by PCApply() if necessary. 221 */ 222 static PetscErrorCode PCSetUp_Deflation(PC pc) 223 { 224 PC_Deflation *jac = (PC_Deflation*)pc->data; 225 Vec diag,diagsqrt; 226 PetscErrorCode ierr; 227 PetscInt n,i; 228 PetscScalar *x; 229 PetscBool zeroflag = PETSC_FALSE; 230 231 PetscFunctionBegin; 232 /* 233 For most preconditioners the code would begin here something like 234 235 if (pc->setupcalled == 0) { allocate space the first time this is ever called 236 ierr = MatCreateVecs(pc->mat,&jac->diag);CHKERRQ(ierr); 237 PetscLogObjectParent((PetscObject)pc,(PetscObject)jac->diag); 238 } 239 240 But for this preconditioner we want to support use of both the matrix' diagonal 241 elements (for left or right preconditioning) and square root of diagonal elements 242 (for symmetric preconditioning). Hence we do not allocate space here, since we 243 don't know at this point which will be needed (diag and/or diagsqrt) until the user 244 applies the preconditioner, and we don't want to allocate BOTH unless we need 245 them both. Thus, the diag and diagsqrt are allocated in PCSetUp_Deflation_NonSymmetric() 246 and PCSetUp_Deflation_Symmetric(), respectively. 247 */ 248 249 /* 250 Here we set up the preconditioner; that is, we copy the diagonal values from 251 the matrix and put them into a format to make them quick to apply as a preconditioner. 252 */ 253 if (zeroflag) { 254 ierr = PetscInfo(pc,"Zero detected in diagonal of matrix, using 1 at those locations\n");CHKERRQ(ierr); 255 } 256 PetscFunctionReturn(0); 257 } 258 /* -------------------------------------------------------------------------- */ 259 /* 260 PCApply_Deflation - Applies the Deflation preconditioner to a vector. 261 262 Input Parameters: 263 . pc - the preconditioner context 264 . x - input vector 265 266 Output Parameter: 267 . y - output vector 268 269 Application Interface Routine: PCApply() 270 */ 271 static PetscErrorCode PCApply_Deflation(PC pc,Vec x,Vec y) 272 { 273 PC_Deflation *jac = (PC_Deflation*)pc->data; 274 PetscErrorCode ierr; 275 276 PetscFunctionBegin; 277 PetscFunctionReturn(0); 278 } 279 /* -------------------------------------------------------------------------- */ 280 static PetscErrorCode PCReset_Deflation(PC pc) 281 { 282 PC_Deflation *jac = (PC_Deflation*)pc->data; 283 PetscErrorCode ierr; 284 285 PetscFunctionBegin; 286 PetscFunctionReturn(0); 287 } 288 289 /* 290 PCDestroy_Deflation - Destroys the private context for the Deflation preconditioner 291 that was created with PCCreate_Deflation(). 292 293 Input Parameter: 294 . pc - the preconditioner context 295 296 Application Interface Routine: PCDestroy() 297 */ 298 static PetscErrorCode PCDestroy_Deflation(PC pc) 299 { 300 PetscErrorCode ierr; 301 302 PetscFunctionBegin; 303 ierr = PCReset_Deflation(pc);CHKERRQ(ierr); 304 305 /* 306 Free the private data structure that was hanging off the PC 307 */ 308 ierr = PetscFree(pc->data);CHKERRQ(ierr); 309 PetscFunctionReturn(0); 310 } 311 312 static PetscErrorCode PCSetFromOptions_Deflation(PetscOptionItems *PetscOptionsObject,PC pc) 313 { 314 PC_Deflation *jac = (PC_Deflation*)pc->data; 315 PetscErrorCode ierr; 316 PetscBool flg; 317 PCDeflationType deflt,type; 318 319 PetscFunctionBegin; 320 ierr = PCDeflationGetType(pc,&deflt);CHKERRQ(ierr); 321 ierr = PetscOptionsHead(PetscOptionsObject,"Deflation options");CHKERRQ(ierr); 322 ierr = PetscOptionsEnum("-pc_jacobi_type","How to construct diagonal matrix","PCDeflationSetType",PCDeflationTypes,(PetscEnum)deflt,(PetscEnum*)&type,&flg);CHKERRQ(ierr); 323 if (flg) { 324 ierr = PCDeflationSetType(pc,type);CHKERRQ(ierr); 325 } 326 ierr = PetscOptionsTail();CHKERRQ(ierr); 327 PetscFunctionReturn(0); 328 } 329 330 /*MC 331 PCDEFLATION - Deflation preconditioner shifts part of the spectrum to zero (deflates) 332 or to a predefined value 333 334 Options Database Key: 335 + -pc_deflation_type <init,pre,post> - selects approach to deflation (default: pre) 336 - -pc_jacobi_abs - use the absolute value of the diagonal entry 337 338 Level: beginner 339 340 Concepts: Deflation, preconditioners 341 342 Notes: 343 todo 344 345 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 346 PCDeflationSetType(), PCDeflationSetSpace() 347 M*/ 348 349 PETSC_EXTERN PetscErrorCode PCCreate_Deflation(PC pc) 350 { 351 PC_Deflation *def; 352 PetscErrorCode ierr; 353 354 PetscFunctionBegin; 355 ierr = PetscNewLog(pc,&def);CHKERRQ(ierr); 356 pc->data = (void*)def; 357 358 def->init = PETSC_FALSE; 359 def->pre = PETSC_TRUE; 360 def->correct = PETSC_FALSE; 361 def->truenorm = PETSC_TRUE; 362 def->reductionfact = -1; 363 def->spacetype = PC_DEFLATION_SPACE_HAAR; 364 def->spacesize = 1; 365 def->extendsp = PETSC_FALSE; 366 def->nestedlvl = 0; 367 def->maxnestedlvl = 0; 368 def->adaptiveconv = PETSC_FALSE; 369 def->adaptiveconst = 1.0; 370 371 /* 372 Set the pointers for the functions that are provided above. 373 Now when the user-level routines (such as PCApply(), PCDestroy(), etc.) 374 are called, they will automatically call these functions. Note we 375 choose not to provide a couple of these functions since they are 376 not needed. 377 */ 378 pc->ops->apply = PCApply_Deflation; 379 pc->ops->applytranspose = PCApply_Deflation; 380 pc->ops->setup = PCSetUp_Deflation; 381 pc->ops->reset = PCReset_Deflation; 382 pc->ops->destroy = PCDestroy_Deflation; 383 pc->ops->setfromoptions = PCSetFromOptions_Deflation; 384 pc->ops->view = 0; 385 pc->ops->applyrichardson = 0; 386 pc->ops->applysymmetricleft = 0; 387 pc->ops->applysymmetricright = 0; 388 389 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCDeflationSetType_C",PCDeflationSetType_Deflation);CHKERRQ(ierr); 390 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCDeflationGetType_C",PCDeflationGetType_Deflation);CHKERRQ(ierr); 391 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCDeflationSetSpace_C",PCDeflationSetSpace_Deflation);CHKERRQ(ierr); 392 PetscFunctionReturn(0); 393 } 394 395