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