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