1
2 #include <petsc/private/pcmgimpl.h> /*I "petscksp.h" I*/
3
4 /* ---------------------------------------------------------------------------*/
5 /*@C
6 PCMGResidualDefault - Default routine to calculate the residual.
7
8 Collective on Mat
9
10 Input Parameters:
11 + mat - the matrix
12 . b - the right-hand-side
13 - x - the approximate solution
14
15 Output Parameter:
16 . r - location to store the residual
17
18 Level: developer
19
20 .seealso: PCMGSetResidual()
21 @*/
PCMGResidualDefault(Mat mat,Vec b,Vec x,Vec r)22 PetscErrorCode PCMGResidualDefault(Mat mat,Vec b,Vec x,Vec r)
23 {
24 PetscErrorCode ierr;
25
26 PetscFunctionBegin;
27 ierr = MatResidual(mat,b,x,r);CHKERRQ(ierr);
28 PetscFunctionReturn(0);
29 }
30
31 /*@
32 PCMGGetCoarseSolve - Gets the solver context to be used on the coarse grid.
33
34 Not Collective
35
36 Input Parameter:
37 . pc - the multigrid context
38
39 Output Parameter:
40 . ksp - the coarse grid solver context
41
42 Level: advanced
43
44 .seealso: PCMGGetSmootherUp(), PCMGGetSmootherDown(), PCMGGetSmoother()
45 @*/
PCMGGetCoarseSolve(PC pc,KSP * ksp)46 PetscErrorCode PCMGGetCoarseSolve(PC pc,KSP *ksp)
47 {
48 PC_MG *mg = (PC_MG*)pc->data;
49 PC_MG_Levels **mglevels = mg->levels;
50
51 PetscFunctionBegin;
52 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
53 *ksp = mglevels[0]->smoothd;
54 PetscFunctionReturn(0);
55 }
56
57 /*@C
58 PCMGSetResidual - Sets the function to be used to calculate the residual
59 on the lth level.
60
61 Logically Collective on PC
62
63 Input Parameters:
64 + pc - the multigrid context
65 . l - the level (0 is coarsest) to supply
66 . residual - function used to form residual, if none is provided the previously provide one is used, if no
67 previous one were provided then a default is used
68 - mat - matrix associated with residual
69
70 Level: advanced
71
72 .seealso: PCMGResidualDefault()
73 @*/
PCMGSetResidual(PC pc,PetscInt l,PetscErrorCode (* residual)(Mat,Vec,Vec,Vec),Mat mat)74 PetscErrorCode PCMGSetResidual(PC pc,PetscInt l,PetscErrorCode (*residual)(Mat,Vec,Vec,Vec),Mat mat)
75 {
76 PC_MG *mg = (PC_MG*)pc->data;
77 PC_MG_Levels **mglevels = mg->levels;
78 PetscErrorCode ierr;
79
80 PetscFunctionBegin;
81 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
82 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
83 if (residual) mglevels[l]->residual = residual;
84 if (!mglevels[l]->residual) mglevels[l]->residual = PCMGResidualDefault;
85 if (mat) {ierr = PetscObjectReference((PetscObject)mat);CHKERRQ(ierr);}
86 ierr = MatDestroy(&mglevels[l]->A);CHKERRQ(ierr);
87 mglevels[l]->A = mat;
88 PetscFunctionReturn(0);
89 }
90
91 /*@
92 PCMGSetInterpolation - Sets the function to be used to calculate the
93 interpolation from l-1 to the lth level
94
95 Logically Collective on PC
96
97 Input Parameters:
98 + pc - the multigrid context
99 . mat - the interpolation operator
100 - l - the level (0 is coarsest) to supply [do not supply 0]
101
102 Level: advanced
103
104 Notes:
105 Usually this is the same matrix used also to set the restriction
106 for the same level.
107
108 One can pass in the interpolation matrix or its transpose; PETSc figures
109 out from the matrix size which one it is.
110
111 .seealso: PCMGSetRestriction()
112 @*/
PCMGSetInterpolation(PC pc,PetscInt l,Mat mat)113 PetscErrorCode PCMGSetInterpolation(PC pc,PetscInt l,Mat mat)
114 {
115 PC_MG *mg = (PC_MG*)pc->data;
116 PC_MG_Levels **mglevels = mg->levels;
117 PetscErrorCode ierr;
118
119 PetscFunctionBegin;
120 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
121 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
122 if (!l) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Do not set interpolation routine for coarsest level");
123 ierr = PetscObjectReference((PetscObject)mat);CHKERRQ(ierr);
124 ierr = MatDestroy(&mglevels[l]->interpolate);CHKERRQ(ierr);
125
126 mglevels[l]->interpolate = mat;
127 PetscFunctionReturn(0);
128 }
129
130 /*@
131 PCMGSetOperators - Sets operator and preconditioning matrix for lth level
132
133 Logically Collective on PC
134
135 Input Parameters:
136 + pc - the multigrid context
137 . Amat - the operator
138 . pmat - the preconditioning operator
139 - l - the level (0 is the coarsest) to supply
140
141 Level: advanced
142
143 .keywords: multigrid, set, interpolate, level
144
145 .seealso: PCMGSetRestriction(), PCMGSetInterpolation()
146 @*/
PCMGSetOperators(PC pc,PetscInt l,Mat Amat,Mat Pmat)147 PetscErrorCode PCMGSetOperators(PC pc,PetscInt l,Mat Amat,Mat Pmat)
148 {
149 PC_MG *mg = (PC_MG*)pc->data;
150 PC_MG_Levels **mglevels = mg->levels;
151 PetscErrorCode ierr;
152
153 PetscFunctionBegin;
154 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
155 PetscValidHeaderSpecific(Amat,MAT_CLASSID,3);
156 PetscValidHeaderSpecific(Pmat,MAT_CLASSID,4);
157 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
158 ierr = KSPSetOperators(mglevels[l]->smoothd,Amat,Pmat);CHKERRQ(ierr);
159 PetscFunctionReturn(0);
160 }
161
162 /*@
163 PCMGGetInterpolation - Gets the function to be used to calculate the
164 interpolation from l-1 to the lth level
165
166 Logically Collective on PC
167
168 Input Parameters:
169 + pc - the multigrid context
170 - l - the level (0 is coarsest) to supply [Do not supply 0]
171
172 Output Parameter:
173 . mat - the interpolation matrix, can be NULL
174
175 Level: advanced
176
177 .seealso: PCMGGetRestriction(), PCMGSetInterpolation(), PCMGGetRScale()
178 @*/
PCMGGetInterpolation(PC pc,PetscInt l,Mat * mat)179 PetscErrorCode PCMGGetInterpolation(PC pc,PetscInt l,Mat *mat)
180 {
181 PC_MG *mg = (PC_MG*)pc->data;
182 PC_MG_Levels **mglevels = mg->levels;
183 PetscErrorCode ierr;
184
185 PetscFunctionBegin;
186 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
187 if (mat) PetscValidPointer(mat,3);
188 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
189 if (l <= 0 || mg->nlevels <= l) SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Level %D must be in range {1,...,%D}",l,mg->nlevels-1);
190 if (!mglevels[l]->interpolate) {
191 if (!mglevels[l]->restrct) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must call PCMGSetInterpolation() or PCMGSetRestriction()");
192 ierr = PCMGSetInterpolation(pc,l,mglevels[l]->restrct);CHKERRQ(ierr);
193 }
194 if (mat) *mat = mglevels[l]->interpolate;
195 PetscFunctionReturn(0);
196 }
197
198 /*@
199 PCMGSetRestriction - Sets the function to be used to restrict dual vectors
200 from level l to l-1.
201
202 Logically Collective on PC
203
204 Input Parameters:
205 + pc - the multigrid context
206 . l - the level (0 is coarsest) to supply [Do not supply 0]
207 - mat - the restriction matrix
208
209 Level: advanced
210
211 Notes:
212 Usually this is the same matrix used also to set the interpolation
213 for the same level.
214
215 One can pass in the interpolation matrix or its transpose; PETSc figures
216 out from the matrix size which one it is.
217
218 If you do not set this, the transpose of the Mat set with PCMGSetInterpolation()
219 is used.
220
221 .seealso: PCMGSetInterpolation()
222 @*/
PCMGSetRestriction(PC pc,PetscInt l,Mat mat)223 PetscErrorCode PCMGSetRestriction(PC pc,PetscInt l,Mat mat)
224 {
225 PetscErrorCode ierr;
226 PC_MG *mg = (PC_MG*)pc->data;
227 PC_MG_Levels **mglevels = mg->levels;
228
229 PetscFunctionBegin;
230 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
231 PetscValidHeaderSpecific(mat,MAT_CLASSID,3);
232 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
233 if (!l) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Do not set restriction routine for coarsest level");
234 ierr = PetscObjectReference((PetscObject)mat);CHKERRQ(ierr);
235 ierr = MatDestroy(&mglevels[l]->restrct);CHKERRQ(ierr);
236
237 mglevels[l]->restrct = mat;
238 PetscFunctionReturn(0);
239 }
240
241 /*@
242 PCMGGetRestriction - Gets the function to be used to restrict dual vectors
243 from level l to l-1.
244
245 Logically Collective on PC
246
247 Input Parameters:
248 + pc - the multigrid context
249 - l - the level (0 is coarsest) to supply [Do not supply 0]
250
251 Output Parameter:
252 . mat - the restriction matrix
253
254 Level: advanced
255
256 .seealso: PCMGGetInterpolation(), PCMGSetRestriction(), PCMGGetRScale(), PCMGGetInjection()
257 @*/
PCMGGetRestriction(PC pc,PetscInt l,Mat * mat)258 PetscErrorCode PCMGGetRestriction(PC pc,PetscInt l,Mat *mat)
259 {
260 PC_MG *mg = (PC_MG*)pc->data;
261 PC_MG_Levels **mglevels = mg->levels;
262 PetscErrorCode ierr;
263
264 PetscFunctionBegin;
265 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
266 if (mat) PetscValidPointer(mat,3);
267 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
268 if (l <= 0 || mg->nlevels <= l) SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Level %D must be in range {1,...,%D}",l,mg->nlevels-1);
269 if (!mglevels[l]->restrct) {
270 if (!mglevels[l]->interpolate) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must call PCMGSetRestriction() or PCMGSetInterpolation()");
271 ierr = PCMGSetRestriction(pc,l,mglevels[l]->interpolate);CHKERRQ(ierr);
272 }
273 if (mat) *mat = mglevels[l]->restrct;
274 PetscFunctionReturn(0);
275 }
276
277 /*@
278 PCMGSetRScale - Sets the pointwise scaling for the restriction operator from level l to l-1.
279
280 Logically Collective on PC
281
282 Input Parameters:
283 + pc - the multigrid context
284 - l - the level (0 is coarsest) to supply [Do not supply 0]
285 . rscale - the scaling
286
287 Level: advanced
288
289 Notes:
290 When evaluating a function on a coarse level one does not want to do F(R * x) one does F(rscale * R * x) where rscale is 1 over the row sums of R. It is preferable to use PCMGSetInjection() to control moving primal vectors.
291
292 .seealso: PCMGSetInterpolation(), PCMGSetRestriction(), PCMGGetRScale(), PCMGSetInjection()
293 @*/
PCMGSetRScale(PC pc,PetscInt l,Vec rscale)294 PetscErrorCode PCMGSetRScale(PC pc,PetscInt l,Vec rscale)
295 {
296 PetscErrorCode ierr;
297 PC_MG *mg = (PC_MG*)pc->data;
298 PC_MG_Levels **mglevels = mg->levels;
299
300 PetscFunctionBegin;
301 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
302 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
303 if (l <= 0 || mg->nlevels <= l) SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Level %D must be in range {1,...,%D}",l,mg->nlevels-1);
304 ierr = PetscObjectReference((PetscObject)rscale);CHKERRQ(ierr);
305 ierr = VecDestroy(&mglevels[l]->rscale);CHKERRQ(ierr);
306
307 mglevels[l]->rscale = rscale;
308 PetscFunctionReturn(0);
309 }
310
311 /*@
312 PCMGGetRScale - Gets the pointwise scaling for the restriction operator from level l to l-1.
313
314 Collective on PC
315
316 Input Parameters:
317 + pc - the multigrid context
318 . rscale - the scaling
319 - l - the level (0 is coarsest) to supply [Do not supply 0]
320
321 Level: advanced
322
323 Notes:
324 When evaluating a function on a coarse level one does not want to do F(R * x) one does F(rscale * R * x) where rscale is 1 over the row sums of R. It is preferable to use PCMGGetInjection() to control moving primal vectors.
325
326 .seealso: PCMGSetInterpolation(), PCMGGetRestriction(), PCMGGetInjection()
327 @*/
PCMGGetRScale(PC pc,PetscInt l,Vec * rscale)328 PetscErrorCode PCMGGetRScale(PC pc,PetscInt l,Vec *rscale)
329 {
330 PetscErrorCode ierr;
331 PC_MG *mg = (PC_MG*)pc->data;
332 PC_MG_Levels **mglevels = mg->levels;
333
334 PetscFunctionBegin;
335 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
336 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
337 if (l <= 0 || mg->nlevels <= l) SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Level %D must be in range {1,...,%D}",l,mg->nlevels-1);
338 if (!mglevels[l]->rscale) {
339 Mat R;
340 Vec X,Y,coarse,fine;
341 PetscInt M,N;
342 ierr = PCMGGetRestriction(pc,l,&R);CHKERRQ(ierr);
343 ierr = MatCreateVecs(R,&X,&Y);CHKERRQ(ierr);
344 ierr = MatGetSize(R,&M,&N);CHKERRQ(ierr);
345 if (M < N) {
346 fine = X;
347 coarse = Y;
348 } else if (N < M) {
349 fine = Y; coarse = X;
350 } else SETERRQ(PetscObjectComm((PetscObject)R),PETSC_ERR_SUP,"Restriction matrix is square, cannot determine which Vec is coarser");
351 ierr = VecSet(fine,1.);CHKERRQ(ierr);
352 ierr = MatRestrict(R,fine,coarse);CHKERRQ(ierr);
353 ierr = VecDestroy(&fine);CHKERRQ(ierr);
354 ierr = VecReciprocal(coarse);CHKERRQ(ierr);
355 mglevels[l]->rscale = coarse;
356 }
357 *rscale = mglevels[l]->rscale;
358 PetscFunctionReturn(0);
359 }
360
361 /*@
362 PCMGSetInjection - Sets the function to be used to inject primal vectors
363 from level l to l-1.
364
365 Logically Collective on PC
366
367 Input Parameters:
368 + pc - the multigrid context
369 . l - the level (0 is coarsest) to supply [Do not supply 0]
370 - mat - the injection matrix
371
372 Level: advanced
373
374 .seealso: PCMGSetRestriction()
375 @*/
PCMGSetInjection(PC pc,PetscInt l,Mat mat)376 PetscErrorCode PCMGSetInjection(PC pc,PetscInt l,Mat mat)
377 {
378 PetscErrorCode ierr;
379 PC_MG *mg = (PC_MG*)pc->data;
380 PC_MG_Levels **mglevels = mg->levels;
381
382 PetscFunctionBegin;
383 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
384 PetscValidHeaderSpecific(mat,MAT_CLASSID,3);
385 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
386 if (!l) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Do not set restriction routine for coarsest level");
387 ierr = PetscObjectReference((PetscObject)mat);CHKERRQ(ierr);
388 ierr = MatDestroy(&mglevels[l]->inject);CHKERRQ(ierr);
389
390 mglevels[l]->inject = mat;
391 PetscFunctionReturn(0);
392 }
393
394 /*@
395 PCMGGetInjection - Gets the function to be used to inject primal vectors
396 from level l to l-1.
397
398 Logically Collective on PC
399
400 Input Parameters:
401 + pc - the multigrid context
402 - l - the level (0 is coarsest) to supply [Do not supply 0]
403
404 Output Parameter:
405 . mat - the restriction matrix (may be NULL if no injection is available).
406
407 Level: advanced
408
409 .seealso: PCMGSetInjection(), PCMGetGetRestriction()
410 @*/
PCMGGetInjection(PC pc,PetscInt l,Mat * mat)411 PetscErrorCode PCMGGetInjection(PC pc,PetscInt l,Mat *mat)
412 {
413 PC_MG *mg = (PC_MG*)pc->data;
414 PC_MG_Levels **mglevels = mg->levels;
415
416 PetscFunctionBegin;
417 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
418 if (mat) PetscValidPointer(mat,3);
419 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
420 if (l <= 0 || mg->nlevels <= l) SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Level %D must be in range {1,...,%D}",l,mg->nlevels-1);
421 if (mat) *mat = mglevels[l]->inject;
422 PetscFunctionReturn(0);
423 }
424
425 /*@
426 PCMGGetSmoother - Gets the KSP context to be used as smoother for
427 both pre- and post-smoothing. Call both PCMGGetSmootherUp() and
428 PCMGGetSmootherDown() to use different functions for pre- and
429 post-smoothing.
430
431 Not Collective, KSP returned is parallel if PC is
432
433 Input Parameters:
434 + pc - the multigrid context
435 - l - the level (0 is coarsest) to supply
436
437 Ouput Parameters:
438 . ksp - the smoother
439
440 Notes:
441 Once you have called this routine, you can call KSPSetOperators(ksp,...) on the resulting ksp to provide the operators for the smoother for this level.
442 You can also modify smoother options by calling the various KSPSetXXX() options on this ksp. In addition you can call KSPGetPC(ksp,&pc)
443 and modify PC options for the smoother; for example PCSetType(pc,PCSOR); to use SOR smoothing.
444
445 Level: advanced
446
447 .seealso: PCMGGetSmootherUp(), PCMGGetSmootherDown(), PCMGGetCoarseSolve()
448 @*/
PCMGGetSmoother(PC pc,PetscInt l,KSP * ksp)449 PetscErrorCode PCMGGetSmoother(PC pc,PetscInt l,KSP *ksp)
450 {
451 PC_MG *mg = (PC_MG*)pc->data;
452 PC_MG_Levels **mglevels = mg->levels;
453
454 PetscFunctionBegin;
455 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
456 *ksp = mglevels[l]->smoothd;
457 PetscFunctionReturn(0);
458 }
459
460 /*@
461 PCMGGetSmootherUp - Gets the KSP context to be used as smoother after
462 coarse grid correction (post-smoother).
463
464 Not Collective, KSP returned is parallel if PC is
465
466 Input Parameters:
467 + pc - the multigrid context
468 - l - the level (0 is coarsest) to supply
469
470 Ouput Parameters:
471 . ksp - the smoother
472
473 Level: advanced
474
475 Notes:
476 calling this will result in a different pre and post smoother so you may need to
477 set options on the pre smoother also
478
479 .seealso: PCMGGetSmootherUp(), PCMGGetSmootherDown()
480 @*/
PCMGGetSmootherUp(PC pc,PetscInt l,KSP * ksp)481 PetscErrorCode PCMGGetSmootherUp(PC pc,PetscInt l,KSP *ksp)
482 {
483 PC_MG *mg = (PC_MG*)pc->data;
484 PC_MG_Levels **mglevels = mg->levels;
485 PetscErrorCode ierr;
486 const char *prefix;
487 MPI_Comm comm;
488
489 PetscFunctionBegin;
490 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
491 /*
492 This is called only if user wants a different pre-smoother from post.
493 Thus we check if a different one has already been allocated,
494 if not we allocate it.
495 */
496 if (!l) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"There is no such thing as a up smoother on the coarse grid");
497 if (mglevels[l]->smoothu == mglevels[l]->smoothd) {
498 KSPType ksptype;
499 PCType pctype;
500 PC ipc;
501 PetscReal rtol,abstol,dtol;
502 PetscInt maxits;
503 KSPNormType normtype;
504 ierr = PetscObjectGetComm((PetscObject)mglevels[l]->smoothd,&comm);CHKERRQ(ierr);
505 ierr = KSPGetOptionsPrefix(mglevels[l]->smoothd,&prefix);CHKERRQ(ierr);
506 ierr = KSPGetTolerances(mglevels[l]->smoothd,&rtol,&abstol,&dtol,&maxits);CHKERRQ(ierr);
507 ierr = KSPGetType(mglevels[l]->smoothd,&ksptype);CHKERRQ(ierr);
508 ierr = KSPGetNormType(mglevels[l]->smoothd,&normtype);CHKERRQ(ierr);
509 ierr = KSPGetPC(mglevels[l]->smoothd,&ipc);CHKERRQ(ierr);
510 ierr = PCGetType(ipc,&pctype);CHKERRQ(ierr);
511
512 ierr = KSPCreate(comm,&mglevels[l]->smoothu);CHKERRQ(ierr);
513 ierr = KSPSetErrorIfNotConverged(mglevels[l]->smoothu,pc->erroriffailure);CHKERRQ(ierr);
514 ierr = PetscObjectIncrementTabLevel((PetscObject)mglevels[l]->smoothu,(PetscObject)pc,mglevels[0]->levels-l);CHKERRQ(ierr);
515 ierr = KSPSetOptionsPrefix(mglevels[l]->smoothu,prefix);CHKERRQ(ierr);
516 ierr = KSPSetTolerances(mglevels[l]->smoothu,rtol,abstol,dtol,maxits);CHKERRQ(ierr);
517 ierr = KSPSetType(mglevels[l]->smoothu,ksptype);CHKERRQ(ierr);
518 ierr = KSPSetNormType(mglevels[l]->smoothu,normtype);CHKERRQ(ierr);
519 ierr = KSPSetConvergenceTest(mglevels[l]->smoothu,KSPConvergedSkip,NULL,NULL);CHKERRQ(ierr);
520 ierr = KSPGetPC(mglevels[l]->smoothu,&ipc);CHKERRQ(ierr);
521 ierr = PCSetType(ipc,pctype);CHKERRQ(ierr);
522 ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)mglevels[l]->smoothu);CHKERRQ(ierr);
523 ierr = PetscObjectComposedDataSetInt((PetscObject) mglevels[l]->smoothu, PetscMGLevelId, mglevels[l]->level);CHKERRQ(ierr);
524 }
525 if (ksp) *ksp = mglevels[l]->smoothu;
526 PetscFunctionReturn(0);
527 }
528
529 /*@
530 PCMGGetSmootherDown - Gets the KSP context to be used as smoother before
531 coarse grid correction (pre-smoother).
532
533 Not Collective, KSP returned is parallel if PC is
534
535 Input Parameters:
536 + pc - the multigrid context
537 - l - the level (0 is coarsest) to supply
538
539 Ouput Parameters:
540 . ksp - the smoother
541
542 Level: advanced
543
544 Notes:
545 calling this will result in a different pre and post smoother so you may need to
546 set options on the post smoother also
547
548 .seealso: PCMGGetSmootherUp(), PCMGGetSmoother()
549 @*/
PCMGGetSmootherDown(PC pc,PetscInt l,KSP * ksp)550 PetscErrorCode PCMGGetSmootherDown(PC pc,PetscInt l,KSP *ksp)
551 {
552 PetscErrorCode ierr;
553 PC_MG *mg = (PC_MG*)pc->data;
554 PC_MG_Levels **mglevels = mg->levels;
555
556 PetscFunctionBegin;
557 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
558 /* make sure smoother up and down are different */
559 if (l) {
560 ierr = PCMGGetSmootherUp(pc,l,NULL);CHKERRQ(ierr);
561 }
562 *ksp = mglevels[l]->smoothd;
563 PetscFunctionReturn(0);
564 }
565
566 /*@
567 PCMGSetCycleTypeOnLevel - Sets the type of cycle (aka cycle index) to run on the specified level.
568
569 Logically Collective on PC
570
571 Input Parameters:
572 + pc - the multigrid context
573 . l - the level (0 is coarsest)
574 - c - either PC_MG_CYCLE_V or PC_MG_CYCLE_W
575
576 Level: advanced
577
578 .seealso: PCMGSetCycleType()
579 @*/
PCMGSetCycleTypeOnLevel(PC pc,PetscInt l,PCMGCycleType c)580 PetscErrorCode PCMGSetCycleTypeOnLevel(PC pc,PetscInt l,PCMGCycleType c)
581 {
582 PC_MG *mg = (PC_MG*)pc->data;
583 PC_MG_Levels **mglevels = mg->levels;
584
585 PetscFunctionBegin;
586 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
587 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
588 PetscValidLogicalCollectiveInt(pc,l,2);
589 PetscValidLogicalCollectiveEnum(pc,c,3);
590 mglevels[l]->cycles = c;
591 PetscFunctionReturn(0);
592 }
593
594 /*@
595 PCMGSetRhs - Sets the vector to be used to store the right-hand side on a particular level.
596
597 Logically Collective on PC
598
599 Input Parameters:
600 + pc - the multigrid context
601 . l - the level (0 is coarsest) this is to be used for
602 - c - the Vec
603
604 Level: advanced
605
606 Notes:
607 If this is not provided PETSc will automatically generate one. You do not need to keep a reference to this vector if you do not need it. PCDestroy() will properly free it.
608
609 .keywords: MG, multigrid, set, right-hand-side, rhs, level
610 .seealso: PCMGSetX(), PCMGSetR()
611 @*/
PCMGSetRhs(PC pc,PetscInt l,Vec c)612 PetscErrorCode PCMGSetRhs(PC pc,PetscInt l,Vec c)
613 {
614 PetscErrorCode ierr;
615 PC_MG *mg = (PC_MG*)pc->data;
616 PC_MG_Levels **mglevels = mg->levels;
617
618 PetscFunctionBegin;
619 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
620 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
621 if (l == mglevels[0]->levels-1) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_INCOMP,"Do not set rhs for finest level");
622 ierr = PetscObjectReference((PetscObject)c);CHKERRQ(ierr);
623 ierr = VecDestroy(&mglevels[l]->b);CHKERRQ(ierr);
624
625 mglevels[l]->b = c;
626 PetscFunctionReturn(0);
627 }
628
629 /*@
630 PCMGSetX - Sets the vector to be used to store the solution on a particular level.
631
632 Logically Collective on PC
633
634 Input Parameters:
635 + pc - the multigrid context
636 . l - the level (0 is coarsest) this is to be used for (do not supply the finest level)
637 - c - the Vec
638
639 Level: advanced
640
641 Notes:
642 If this is not provided PETSc will automatically generate one. You do not need to keep a reference to this vector if you do not need it. PCDestroy() will properly free it.
643
644 .keywords: MG, multigrid, set, solution, level
645 .seealso: PCMGSetRhs(), PCMGSetR()
646 @*/
PCMGSetX(PC pc,PetscInt l,Vec c)647 PetscErrorCode PCMGSetX(PC pc,PetscInt l,Vec c)
648 {
649 PetscErrorCode ierr;
650 PC_MG *mg = (PC_MG*)pc->data;
651 PC_MG_Levels **mglevels = mg->levels;
652
653 PetscFunctionBegin;
654 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
655 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
656 if (l == mglevels[0]->levels-1) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_INCOMP,"Do not set x for finest level");
657 ierr = PetscObjectReference((PetscObject)c);CHKERRQ(ierr);
658 ierr = VecDestroy(&mglevels[l]->x);CHKERRQ(ierr);
659
660 mglevels[l]->x = c;
661 PetscFunctionReturn(0);
662 }
663
664 /*@
665 PCMGSetR - Sets the vector to be used to store the residual on a particular level.
666
667 Logically Collective on PC
668
669 Input Parameters:
670 + pc - the multigrid context
671 . l - the level (0 is coarsest) this is to be used for
672 - c - the Vec
673
674 Level: advanced
675
676 Notes:
677 If this is not provided PETSc will automatically generate one. You do not need to keep a reference to this vector if you do not need it. PCDestroy() will properly free it.
678
679 .keywords: MG, multigrid, set, residual, level
680 .seealso: PCMGSetRhs(), PCMGSetX()
681 @*/
PCMGSetR(PC pc,PetscInt l,Vec c)682 PetscErrorCode PCMGSetR(PC pc,PetscInt l,Vec c)
683 {
684 PetscErrorCode ierr;
685 PC_MG *mg = (PC_MG*)pc->data;
686 PC_MG_Levels **mglevels = mg->levels;
687
688 PetscFunctionBegin;
689 PetscValidHeaderSpecific(pc,PC_CLASSID,1);
690 if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
691 if (!l) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Need not set residual vector for coarse grid");
692 ierr = PetscObjectReference((PetscObject)c);CHKERRQ(ierr);
693 ierr = VecDestroy(&mglevels[l]->r);CHKERRQ(ierr);
694
695 mglevels[l]->r = c;
696 PetscFunctionReturn(0);
697 }
698