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