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 static PetscErrorCode PCDeflationGetType_Deflation(PC pc,PCDeflationType *type) 94 { 95 PC_Deflation *def = (PC_Deflation*)pc->data; 96 97 PetscFunctionBegin; 98 if (def->init) { 99 *type = PC_DEFLATION_INIT; 100 } else if (def->pre) { 101 *type = PC_DEFLATION_PRE; 102 } else { 103 *type = PC_DEFLATION_POST; 104 } 105 PetscFunctionReturn(0); 106 } 107 108 /* -------------------------------------------------------------------------- */ 109 /* 110 PCSetUp_Deflation - Prepares for the use of the Deflation preconditioner 111 by setting data structures and options. 112 113 Input Parameter: 114 . pc - the preconditioner context 115 116 Application Interface Routine: PCSetUp() 117 118 Notes: 119 The interface routine PCSetUp() is not usually called directly by 120 the user, but instead is called by PCApply() if necessary. 121 */ 122 static PetscErrorCode PCSetUp_Deflation(PC pc) 123 { 124 PC_Deflation *jac = (PC_Deflation*)pc->data; 125 Vec diag,diagsqrt; 126 PetscErrorCode ierr; 127 PetscInt n,i; 128 PetscScalar *x; 129 PetscBool zeroflag = PETSC_FALSE; 130 131 PetscFunctionBegin; 132 /* 133 For most preconditioners the code would begin here something like 134 135 if (pc->setupcalled == 0) { allocate space the first time this is ever called 136 ierr = MatCreateVecs(pc->mat,&jac->diag);CHKERRQ(ierr); 137 PetscLogObjectParent((PetscObject)pc,(PetscObject)jac->diag); 138 } 139 140 But for this preconditioner we want to support use of both the matrix' diagonal 141 elements (for left or right preconditioning) and square root of diagonal elements 142 (for symmetric preconditioning). Hence we do not allocate space here, since we 143 don't know at this point which will be needed (diag and/or diagsqrt) until the user 144 applies the preconditioner, and we don't want to allocate BOTH unless we need 145 them both. Thus, the diag and diagsqrt are allocated in PCSetUp_Deflation_NonSymmetric() 146 and PCSetUp_Deflation_Symmetric(), respectively. 147 */ 148 149 /* 150 Here we set up the preconditioner; that is, we copy the diagonal values from 151 the matrix and put them into a format to make them quick to apply as a preconditioner. 152 */ 153 if (zeroflag) { 154 ierr = PetscInfo(pc,"Zero detected in diagonal of matrix, using 1 at those locations\n");CHKERRQ(ierr); 155 } 156 PetscFunctionReturn(0); 157 } 158 /* -------------------------------------------------------------------------- */ 159 /* 160 PCApply_Deflation - Applies the Deflation preconditioner to a vector. 161 162 Input Parameters: 163 . pc - the preconditioner context 164 . x - input vector 165 166 Output Parameter: 167 . y - output vector 168 169 Application Interface Routine: PCApply() 170 */ 171 static PetscErrorCode PCApply_Deflation(PC pc,Vec x,Vec y) 172 { 173 PC_Deflation *jac = (PC_Deflation*)pc->data; 174 PetscErrorCode ierr; 175 176 PetscFunctionBegin; 177 PetscFunctionReturn(0); 178 } 179 /* -------------------------------------------------------------------------- */ 180 static PetscErrorCode PCReset_Deflation(PC pc) 181 { 182 PC_Deflation *jac = (PC_Deflation*)pc->data; 183 PetscErrorCode ierr; 184 185 PetscFunctionBegin; 186 PetscFunctionReturn(0); 187 } 188 189 /* 190 PCDestroy_Deflation - Destroys the private context for the Deflation preconditioner 191 that was created with PCCreate_Deflation(). 192 193 Input Parameter: 194 . pc - the preconditioner context 195 196 Application Interface Routine: PCDestroy() 197 */ 198 static PetscErrorCode PCDestroy_Deflation(PC pc) 199 { 200 PetscErrorCode ierr; 201 202 PetscFunctionBegin; 203 ierr = PCReset_Deflation(pc);CHKERRQ(ierr); 204 205 /* 206 Free the private data structure that was hanging off the PC 207 */ 208 ierr = PetscFree(pc->data);CHKERRQ(ierr); 209 PetscFunctionReturn(0); 210 } 211 212 static PetscErrorCode PCSetFromOptions_Deflation(PetscOptionItems *PetscOptionsObject,PC pc) 213 { 214 PC_Deflation *jac = (PC_Deflation*)pc->data; 215 PetscErrorCode ierr; 216 PetscBool flg; 217 PCDeflationType deflt,type; 218 219 PetscFunctionBegin; 220 ierr = PCDeflationGetType(pc,&deflt);CHKERRQ(ierr); 221 ierr = PetscOptionsHead(PetscOptionsObject,"Deflation options");CHKERRQ(ierr); 222 ierr = PetscOptionsEnum("-pc_jacobi_type","How to construct diagonal matrix","PCDeflationSetType",PCDeflationTypes,(PetscEnum)deflt,(PetscEnum*)&type,&flg);CHKERRQ(ierr); 223 if (flg) { 224 ierr = PCDeflationSetType(pc,type);CHKERRQ(ierr); 225 } 226 ierr = PetscOptionsTail();CHKERRQ(ierr); 227 PetscFunctionReturn(0); 228 } 229 230 /* -------------------------------------------------------------------------- */ 231 /* 232 PCCreate_Deflation - Creates a Deflation preconditioner context, PC_DeflationPC_Deflation, 233 and sets this as the private data within the generic preconditioning 234 context, PC, that was created within PCCreate(). 235 236 Input Parameter: 237 . pc - the preconditioner context 238 239 Application Interface Routine: PCCreate() 240 */ 241 242 /*MC 243 PCDEFLATION - Deflation (i.e. diagonal scaling preconditioning) 244 245 Options Database Key: 246 + -pc_jacobi_type <diagonal,rowmax,rowsum> - approach for forming the preconditioner 247 - -pc_jacobi_abs - use the absolute value of the diagonal entry 248 249 Level: beginner 250 251 Concepts: Deflation, diagonal scaling, preconditioners 252 253 Notes: 254 By using KSPSetPCSide(ksp,PC_SYMMETRIC) or -ksp_pc_side symmetric 255 can scale each side of the matrix by the square root of the diagonal entries. 256 257 Zero entries along the diagonal are replaced with the value 1.0 258 259 See PCPBDEFLATION for a point-block Deflation preconditioner 260 261 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 262 PCDeflationSetType(), PCDeflationSetUseAbs(), PCDeflationGetUseAbs(), PCPBDEFLATION 263 M*/ 264 265 PETSC_EXTERN PetscErrorCode PCCreate_Deflation(PC pc) 266 { 267 PC_Deflation *def; 268 PetscErrorCode ierr; 269 270 PetscFunctionBegin; 271 /* 272 Creates the private data structure for this preconditioner and 273 attach it to the PC object. 274 */ 275 ierr = PetscNewLog(pc,&def);CHKERRQ(ierr); 276 pc->data = (void*)def; 277 278 /* 279 Initialize the pointers to vectors to ZERO; these will be used to store 280 diagonal entries of the matrix for fast preconditioner application. 281 */ 282 //jac->diag = 0; 283 //jac->diagsqrt = 0; 284 //jac->userowmax = PETSC_FALSE; 285 //jac->userowsum = PETSC_FALSE; 286 //jac->useabs = PETSC_FALSE; 287 288 /* 289 Set the pointers for the functions that are provided above. 290 Now when the user-level routines (such as PCApply(), PCDestroy(), etc.) 291 are called, they will automatically call these functions. Note we 292 choose not to provide a couple of these functions since they are 293 not needed. 294 */ 295 pc->ops->apply = PCApply_Deflation; 296 pc->ops->applytranspose = PCApply_Deflation; 297 pc->ops->setup = PCSetUp_Deflation; 298 pc->ops->reset = PCReset_Deflation; 299 pc->ops->destroy = PCDestroy_Deflation; 300 pc->ops->setfromoptions = PCSetFromOptions_Deflation; 301 pc->ops->view = 0; 302 pc->ops->applyrichardson = 0; 303 pc->ops->applysymmetricleft = 0; 304 pc->ops->applysymmetricright = 0; 305 306 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCDeflationSetType_C",PCDeflationSetType_Deflation);CHKERRQ(ierr); 307 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCDeflationGetType_C",PCDeflationGetType_Deflation);CHKERRQ(ierr); 308 PetscFunctionReturn(0); 309 } 310 311 /*@ 312 PCDeflationSetType - Causes the deflation preconditioner to use only a special 313 initial gues or pre/post solve solution update 314 315 Logically Collective on PC 316 317 Input Parameters: 318 + pc - the preconditioner context 319 - type - PC_DEFLATION_PRE, PC_DEFLATION_INIT, PC_DEFLATION_POST 320 321 Options Database Key: 322 . -pc_deflation_type <pre,init,post> 323 324 Level: intermediate 325 326 Concepts: Deflation preconditioner 327 328 .seealso: PCDeflationGetType() 329 @*/ 330 PetscErrorCode PCDeflationSetType(PC pc,PCDeflationType type) 331 { 332 PetscErrorCode ierr; 333 334 PetscFunctionBegin; 335 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 336 ierr = PetscTryMethod(pc,"PCDeflationSetType_C",(PC,PCDeflationType),(pc,type));CHKERRQ(ierr); 337 PetscFunctionReturn(0); 338 } 339 340 /*@ 341 PCDeflationGetType - Gets how the diagonal matrix is produced for the preconditioner 342 343 Not Collective on PC 344 345 Input Parameter: 346 . pc - the preconditioner context 347 348 Output Parameter: 349 - type - PC_DEFLATION_PRE, PC_DEFLATION_INIT, PC_DEFLATION_POST 350 351 Level: intermediate 352 353 Concepts: Deflation preconditioner 354 355 .seealso: PCDeflationSetType() 356 @*/ 357 PetscErrorCode PCDeflationGetType(PC pc,PCDeflationType *type) 358 { 359 PetscErrorCode ierr; 360 361 PetscFunctionBegin; 362 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 363 ierr = PetscUseMethod(pc,"PCDeflationGetType_C",(PC,PCDeflationType*),(pc,type));CHKERRQ(ierr); 364 PetscFunctionReturn(0); 365 } 366