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