1 #define PETSCKSP_DLL 2 3 #include "../src/ksp/pc/impls/mg/mgimpl.h" /*I "petscksp.h" I*/ 4 /*I "petscmg.h" I*/ 5 6 #undef __FUNCT__ 7 #define __FUNCT__ "PCMGDefaultResidual" 8 /*@C 9 PCMGDefaultResidual - Default routine to calculate the residual. 10 11 Collective on Mat and Vec 12 13 Input Parameters: 14 + mat - the matrix 15 . b - the right-hand-side 16 - x - the approximate solution 17 18 Output Parameter: 19 . r - location to store the residual 20 21 Level: advanced 22 23 .keywords: MG, default, multigrid, residual 24 25 .seealso: PCMGSetResidual() 26 @*/ 27 PetscErrorCode PETSCKSP_DLLEXPORT PCMGDefaultResidual(Mat mat,Vec b,Vec x,Vec r) 28 { 29 PetscErrorCode ierr; 30 31 PetscFunctionBegin; 32 ierr = MatMult(mat,x,r);CHKERRQ(ierr); 33 ierr = VecAYPX(r,-1.0,b);CHKERRQ(ierr); 34 PetscFunctionReturn(0); 35 } 36 37 /* ---------------------------------------------------------------------------*/ 38 39 #undef __FUNCT__ 40 #define __FUNCT__ "PCMGGetCoarseSolve" 41 /*@ 42 PCMGGetCoarseSolve - Gets the solver context to be used on the coarse grid. 43 44 Not Collective 45 46 Input Parameter: 47 . pc - the multigrid context 48 49 Output Parameter: 50 . ksp - the coarse grid solver context 51 52 Level: advanced 53 54 .keywords: MG, multigrid, get, coarse grid 55 @*/ 56 PetscErrorCode PETSCKSP_DLLEXPORT PCMGGetCoarseSolve(PC pc,KSP *ksp) 57 { 58 PC_MG *mg = (PC_MG*)pc->data; 59 PC_MG_Levels **mglevels = mg->levels; 60 61 PetscFunctionBegin; 62 *ksp = mglevels[0]->smoothd; 63 PetscFunctionReturn(0); 64 } 65 66 #undef __FUNCT__ 67 #define __FUNCT__ "PCMGSetResidual" 68 /*@C 69 PCMGSetResidual - Sets the function to be used to calculate the residual 70 on the lth level. 71 72 Collective on PC and Mat 73 74 Input Parameters: 75 + pc - the multigrid context 76 . l - the level (0 is coarsest) to supply 77 . residual - function used to form residual (usually PCMGDefaultResidual) 78 - mat - matrix associated with residual 79 80 Level: advanced 81 82 .keywords: MG, set, multigrid, residual, level 83 84 .seealso: PCMGDefaultResidual() 85 @*/ 86 PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetResidual(PC pc,PetscInt l,PetscErrorCode (*residual)(Mat,Vec,Vec,Vec),Mat mat) 87 { 88 PC_MG *mg = (PC_MG*)pc->data; 89 PC_MG_Levels **mglevels = mg->levels; 90 91 PetscFunctionBegin; 92 if (!mglevels) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling"); 93 mglevels[l]->residual = residual; 94 mglevels[l]->A = mat; 95 PetscFunctionReturn(0); 96 } 97 98 #undef __FUNCT__ 99 #define __FUNCT__ "PCMGSetInterpolation" 100 /*@ 101 PCMGSetInterpolation - Sets the function to be used to calculate the 102 interpolation from l-1 to the lth level 103 104 Collective on PC and Mat 105 106 Input Parameters: 107 + pc - the multigrid context 108 . mat - the interpolation operator 109 - l - the level (0 is coarsest) to supply [do not supply 0] 110 111 Level: advanced 112 113 Notes: 114 Usually this is the same matrix used also to set the restriction 115 for the same level. 116 117 One can pass in the interpolation matrix or its transpose; PETSc figures 118 out from the matrix size which one it is. 119 120 .keywords: multigrid, set, interpolate, level 121 122 .seealso: PCMGSetRestriction() 123 @*/ 124 PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetInterpolation(PC pc,PetscInt l,Mat mat) 125 { 126 PC_MG *mg = (PC_MG*)pc->data; 127 PC_MG_Levels **mglevels = mg->levels; 128 PetscErrorCode ierr; 129 130 PetscFunctionBegin; 131 if (!mglevels) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling"); 132 if (!l) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Do not set interpolation routine for coarsest level"); 133 ierr = PetscObjectReference((PetscObject)mat);CHKERRQ(ierr); 134 if (mglevels[l]->interpolate) {ierr = MatDestroy(mglevels[l]->interpolate);CHKERRQ(ierr);} 135 mglevels[l]->interpolate = mat; 136 PetscFunctionReturn(0); 137 } 138 139 #undef __FUNCT__ 140 #define __FUNCT__ "PCMGSetRestriction" 141 /*@ 142 PCMGSetRestriction - Sets the function to be used to restrict vector 143 from level l to l-1. 144 145 Collective on PC and Mat 146 147 Input Parameters: 148 + pc - the multigrid context 149 . mat - the restriction matrix 150 - l - the level (0 is coarsest) to supply [Do not supply 0] 151 152 Level: advanced 153 154 Notes: 155 Usually this is the same matrix used also to set the interpolation 156 for the same level. 157 158 One can pass in the interpolation matrix or its transpose; PETSc figures 159 out from the matrix size which one it is. 160 161 If you do not set this, the transpose of the Mat set with PCMGSetInterpolation() 162 is used. 163 164 .keywords: MG, set, multigrid, restriction, level 165 166 .seealso: PCMGSetInterpolation() 167 @*/ 168 PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetRestriction(PC pc,PetscInt l,Mat mat) 169 { 170 PetscErrorCode ierr; 171 PC_MG *mg = (PC_MG*)pc->data; 172 PC_MG_Levels **mglevels = mg->levels; 173 174 PetscFunctionBegin; 175 if (!mglevels) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling"); 176 if (!l) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Do not set restriction routine for coarsest level"); 177 ierr = PetscObjectReference((PetscObject)mat);CHKERRQ(ierr); 178 if (mglevels[l]->restrct) {ierr = MatDestroy(mglevels[l]->restrct);CHKERRQ(ierr);} 179 mglevels[l]->restrct = mat; 180 PetscFunctionReturn(0); 181 } 182 183 #undef __FUNCT__ 184 #define __FUNCT__ "PCMGGetSmoother" 185 /*@ 186 PCMGGetSmoother - Gets the KSP context to be used as smoother for 187 both pre- and post-smoothing. Call both PCMGGetSmootherUp() and 188 PCMGGetSmootherDown() to use different functions for pre- and 189 post-smoothing. 190 191 Not Collective, KSP returned is parallel if PC is 192 193 Input Parameters: 194 + pc - the multigrid context 195 - l - the level (0 is coarsest) to supply 196 197 Ouput Parameters: 198 . ksp - the smoother 199 200 Level: advanced 201 202 .keywords: MG, get, multigrid, level, smoother, pre-smoother, post-smoother 203 204 .seealso: PCMGGetSmootherUp(), PCMGGetSmootherDown() 205 @*/ 206 PetscErrorCode PETSCKSP_DLLEXPORT PCMGGetSmoother(PC pc,PetscInt l,KSP *ksp) 207 { 208 PC_MG *mg = (PC_MG*)pc->data; 209 PC_MG_Levels **mglevels = mg->levels; 210 211 PetscFunctionBegin; 212 *ksp = mglevels[l]->smoothd; 213 PetscFunctionReturn(0); 214 } 215 216 #undef __FUNCT__ 217 #define __FUNCT__ "PCMGGetSmootherUp" 218 /*@ 219 PCMGGetSmootherUp - Gets the KSP context to be used as smoother after 220 coarse grid correction (post-smoother). 221 222 Not Collective, KSP returned is parallel if PC is 223 224 Input Parameters: 225 + pc - the multigrid context 226 - l - the level (0 is coarsest) to supply 227 228 Ouput Parameters: 229 . ksp - the smoother 230 231 Level: advanced 232 233 .keywords: MG, multigrid, get, smoother, up, post-smoother, level 234 235 .seealso: PCMGGetSmootherUp(), PCMGGetSmootherDown() 236 @*/ 237 PetscErrorCode PETSCKSP_DLLEXPORT PCMGGetSmootherUp(PC pc,PetscInt l,KSP *ksp) 238 { 239 PC_MG *mg = (PC_MG*)pc->data; 240 PC_MG_Levels **mglevels = mg->levels; 241 PetscErrorCode ierr; 242 const char *prefix; 243 MPI_Comm comm; 244 245 PetscFunctionBegin; 246 /* 247 This is called only if user wants a different pre-smoother from post. 248 Thus we check if a different one has already been allocated, 249 if not we allocate it. 250 */ 251 if (!l) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_OUTOFRANGE,"There is no such thing as a up smoother on the coarse grid"); 252 if (mglevels[l]->smoothu == mglevels[l]->smoothd) { 253 ierr = PetscObjectGetComm((PetscObject)mglevels[l]->smoothd,&comm);CHKERRQ(ierr); 254 ierr = KSPGetOptionsPrefix(mglevels[l]->smoothd,&prefix);CHKERRQ(ierr); 255 ierr = KSPCreate(comm,&mglevels[l]->smoothu);CHKERRQ(ierr); 256 ierr = PetscObjectIncrementTabLevel((PetscObject)mglevels[l]->smoothu,(PetscObject)pc,mglevels[0]->levels-l);CHKERRQ(ierr); 257 ierr = KSPSetTolerances(mglevels[l]->smoothu,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,1);CHKERRQ(ierr); 258 ierr = KSPSetOptionsPrefix(mglevels[l]->smoothu,prefix);CHKERRQ(ierr); 259 ierr = PetscLogObjectParent(pc,mglevels[l]->smoothu);CHKERRQ(ierr); 260 } 261 if (ksp) *ksp = mglevels[l]->smoothu; 262 PetscFunctionReturn(0); 263 } 264 265 #undef __FUNCT__ 266 #define __FUNCT__ "PCMGGetSmootherDown" 267 /*@ 268 PCMGGetSmootherDown - Gets the KSP context to be used as smoother before 269 coarse grid correction (pre-smoother). 270 271 Not Collective, KSP returned is parallel if PC is 272 273 Input Parameters: 274 + pc - the multigrid context 275 - l - the level (0 is coarsest) to supply 276 277 Ouput Parameters: 278 . ksp - the smoother 279 280 Level: advanced 281 282 .keywords: MG, multigrid, get, smoother, down, pre-smoother, level 283 284 .seealso: PCMGGetSmootherUp(), PCMGGetSmoother() 285 @*/ 286 PetscErrorCode PETSCKSP_DLLEXPORT PCMGGetSmootherDown(PC pc,PetscInt l,KSP *ksp) 287 { 288 PetscErrorCode ierr; 289 PC_MG *mg = (PC_MG*)pc->data; 290 PC_MG_Levels **mglevels = mg->levels; 291 292 PetscFunctionBegin; 293 /* make sure smoother up and down are different */ 294 if (l != 0) { 295 ierr = PCMGGetSmootherUp(pc,l,PETSC_NULL);CHKERRQ(ierr); 296 } 297 *ksp = mglevels[l]->smoothd; 298 PetscFunctionReturn(0); 299 } 300 301 #undef __FUNCT__ 302 #define __FUNCT__ "PCMGSetCyclesOnLevel" 303 /*@ 304 PCMGSetCyclesOnLevel - Sets the number of cycles to run on this level. 305 306 Collective on PC 307 308 Input Parameters: 309 + pc - the multigrid context 310 . l - the level (0 is coarsest) this is to be used for 311 - n - the number of cycles 312 313 Level: advanced 314 315 .keywords: MG, multigrid, set, cycles, V-cycle, W-cycle, level 316 317 .seealso: PCMGSetCycles() 318 @*/ 319 PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetCyclesOnLevel(PC pc,PetscInt l,PetscInt c) 320 { 321 PC_MG *mg = (PC_MG*)pc->data; 322 PC_MG_Levels **mglevels = mg->levels; 323 324 PetscFunctionBegin; 325 if (!mglevels) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling"); 326 mglevels[l]->cycles = c; 327 PetscFunctionReturn(0); 328 } 329 330 #undef __FUNCT__ 331 #define __FUNCT__ "PCMGSetRhs" 332 /*@ 333 PCMGSetRhs - Sets the vector space to be used to store the right-hand side 334 on a particular level. 335 336 Collective on PC and Vec 337 338 Input Parameters: 339 + pc - the multigrid context 340 . l - the level (0 is coarsest) this is to be used for 341 - c - the space 342 343 Level: advanced 344 345 Notes: If this is not provided PETSc will automatically generate one. 346 347 You do not need to keep a reference to this vector if you do 348 not need it PCDestroy() will properly free it. 349 350 .keywords: MG, multigrid, set, right-hand-side, rhs, level 351 352 .seealso: PCMGSetX(), PCMGSetR() 353 @*/ 354 PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetRhs(PC pc,PetscInt l,Vec c) 355 { 356 PetscErrorCode ierr; 357 PC_MG *mg = (PC_MG*)pc->data; 358 PC_MG_Levels **mglevels = mg->levels; 359 360 PetscFunctionBegin; 361 if (!mglevels) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling"); 362 if (l == mglevels[0]->levels-1) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_INCOMP,"Do not set rhs for finest level"); 363 ierr = PetscObjectReference((PetscObject)c);CHKERRQ(ierr); 364 if (mglevels[l]->b) {ierr = VecDestroy(mglevels[l]->b);CHKERRQ(ierr);} 365 mglevels[l]->b = c; 366 PetscFunctionReturn(0); 367 } 368 369 #undef __FUNCT__ 370 #define __FUNCT__ "PCMGSetX" 371 /*@ 372 PCMGSetX - Sets the vector space to be used to store the solution on a 373 particular level. 374 375 Collective on PC and Vec 376 377 Input Parameters: 378 + pc - the multigrid context 379 . l - the level (0 is coarsest) this is to be used for 380 - c - the space 381 382 Level: advanced 383 384 Notes: If this is not provided PETSc will automatically generate one. 385 386 You do not need to keep a reference to this vector if you do 387 not need it PCDestroy() will properly free it. 388 389 .keywords: MG, multigrid, set, solution, level 390 391 .seealso: PCMGSetRhs(), PCMGSetR() 392 @*/ 393 PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetX(PC pc,PetscInt l,Vec c) 394 { 395 PetscErrorCode ierr; 396 PC_MG *mg = (PC_MG*)pc->data; 397 PC_MG_Levels **mglevels = mg->levels; 398 399 PetscFunctionBegin; 400 if (!mglevels) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling"); 401 if (l == mglevels[0]->levels-1) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_INCOMP,"Do not set x for finest level"); 402 ierr = PetscObjectReference((PetscObject)c);CHKERRQ(ierr); 403 if (mglevels[l]->x) {ierr = VecDestroy(mglevels[l]->x);CHKERRQ(ierr);} 404 mglevels[l]->x = c; 405 PetscFunctionReturn(0); 406 } 407 408 #undef __FUNCT__ 409 #define __FUNCT__ "PCMGSetR" 410 /*@ 411 PCMGSetR - Sets the vector space to be used to store the residual on a 412 particular level. 413 414 Collective on PC and Vec 415 416 Input Parameters: 417 + pc - the multigrid context 418 . l - the level (0 is coarsest) this is to be used for 419 - c - the space 420 421 Level: advanced 422 423 Notes: If this is not provided PETSc will automatically generate one. 424 425 You do not need to keep a reference to this vector if you do 426 not need it PCDestroy() will properly free it. 427 428 .keywords: MG, multigrid, set, residual, level 429 @*/ 430 PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetR(PC pc,PetscInt l,Vec c) 431 { 432 PetscErrorCode ierr; 433 PC_MG *mg = (PC_MG*)pc->data; 434 PC_MG_Levels **mglevels = mg->levels; 435 436 PetscFunctionBegin; 437 if (!mglevels) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling"); 438 if (!l) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Need not set residual vector for coarse grid"); 439 ierr = PetscObjectReference((PetscObject)c);CHKERRQ(ierr); 440 if (mglevels[l]->r) {ierr = VecDestroy(mglevels[l]->r);CHKERRQ(ierr);} 441 mglevels[l]->r = c; 442 PetscFunctionReturn(0); 443 } 444