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