Actual source code: mgfunc.c
1: /*$Id: mgfunc.c,v 1.41 2001/08/07 03:03:36 balay Exp $*/
3: #include src/sles/pc/impls/mg/mgimpl.h
4: /*I "petscmg.h" I*/
6: #undef __FUNCT__
8: /*@C
9: MGDefaultResidual - Default routine to calculate the residual.
11: Collective on Mat and Vec
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
21: Level: advanced
23: .keywords: MG, default, multigrid, residual
25: .seealso: MGSetResidual()
26: @*/
27: int MGDefaultResidual(Mat mat,Vec b,Vec x,Vec r)
28: {
29: int ierr;
30: PetscScalar mone = -1.0;
33: MatMult(mat,x,r);
34: VecAYPX(&mone,b,r);
35: return(0);
36: }
38: /* ---------------------------------------------------------------------------*/
40: #undef __FUNCT__
42: /*@C
43: MGGetCoarseSolve - Gets the solver context to be used on the coarse grid.
45: Not Collective
47: Input Parameter:
48: . pc - the multigrid context
50: Output Parameter:
51: . sles - the coarse grid solver context
53: Level: advanced
55: .keywords: MG, multigrid, get, coarse grid
56: @*/
57: int MGGetCoarseSolve(PC pc,SLES *sles)
58: {
59: MG *mg = (MG*)pc->data;
62: *sles = mg[0]->smoothd;
63: return(0);
64: }
66: #undef __FUNCT__
68: /*@C
69: MGSetResidual - Sets the function to be used to calculate the residual
70: on the lth level.
72: Collective on PC and Mat
74: Input Parameters:
75: + pc - the multigrid context
76: . l - the level to supply
77: . residual - function used to form residual (usually MGDefaultResidual)
78: - mat - matrix associated with residual
80: Level: advanced
82: .keywords: MG, set, multigrid, residual, level
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;
91: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
93: mg[l]->residual = residual;
94: mg[l]->A = mat;
95: return(0);
96: }
98: #undef __FUNCT__
100: /*@
101: MGSetInterpolate - Sets the function to be used to calculate the
102: interpolation on the lth level.
104: Collective on PC and Mat
106: Input Parameters:
107: + pc - the multigrid context
108: . mat - the interpolation operator
109: - l - the level to supply
111: Level: advanced
113: Notes:
114: Usually this is the same matrix used also to set the restriction
115: for the same level.
117: One can pass in the interpolation matrix or its transpose; PETSc figures
118: out from the matrix size which one it is.
120: .keywords: multigrid, set, interpolate, level
122: .seealso: MGSetRestriction()
123: @*/
124: int MGSetInterpolate(PC pc,int l,Mat mat)
125: {
126: MG *mg = (MG*)pc->data;
129: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
130: mg[l]->interpolate = mat;
131: return(0);
132: }
134: #undef __FUNCT__
136: /*@
137: MGSetRestriction - Sets the function to be used to restrict vector
138: from level l to l-1.
140: Collective on PC and Mat
142: Input Parameters:
143: + pc - the multigrid context
144: . mat - the restriction matrix
145: - l - the level to supply
147: Level: advanced
149: Notes:
150: Usually this is the same matrix used also to set the interpolation
151: for the same level.
153: One can pass in the interpolation matrix or its transpose; PETSc figures
154: out from the matrix size which one it is.
156: .keywords: MG, set, multigrid, restriction, level
158: .seealso: MGSetInterpolate()
159: @*/
160: int MGSetRestriction(PC pc,int l,Mat mat)
161: {
162: MG *mg = (MG*)pc->data;
165: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
166: mg[l]->restrct = mat;
167: return(0);
168: }
170: #undef __FUNCT__
172: /*@C
173: MGGetSmoother - Gets the SLES 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.
178: Not Collective, SLES returned is parallel if PC is
180: Input Parameters:
181: + pc - the multigrid context
182: - l - the level to supply
184: Ouput Parameters:
185: . sles - the smoother
187: Level: advanced
189: .keywords: MG, get, multigrid, level, smoother, pre-smoother, post-smoother
191: .seealso: MGGetSmootherUp(), MGGetSmootherDown()
192: @*/
193: int MGGetSmoother(PC pc,int l,SLES *sles)
194: {
195: MG *mg = (MG*)pc->data;
198: *sles = mg[l]->smoothd;
199: return(0);
200: }
202: #undef __FUNCT__
204: /*@C
205: MGGetSmootherUp - Gets the SLES context to be used as smoother after
206: coarse grid correction (post-smoother).
208: Not Collective, SLES returned is parallel if PC is
210: Input Parameters:
211: + pc - the multigrid context
212: - l - the level to supply
214: Ouput Parameters:
215: . sles - the smoother
217: Level: advanced
219: .keywords: MG, multigrid, get, smoother, up, post-smoother, level
221: .seealso: MGGetSmootherUp(), MGGetSmootherDown()
222: @*/
223: int MGGetSmootherUp(PC pc,int l,SLES *sles)
224: {
225: MG *mg = (MG*)pc->data;
226: int ierr;
227: char *prefix;
228: KSP ksp;
229: MPI_Comm comm;
232: /*
233: This is called only if user wants a different pre-smoother from post.
234: Thus we check if a different one has already been allocated,
235: if not we allocate it.
236: */
237: PCGetOptionsPrefix(pc,&prefix);
239: if (mg[l]->smoothu == mg[l]->smoothd) {
240: PetscObjectGetComm((PetscObject)mg[l]->smoothd,&comm);
241: SLESCreate(comm,&mg[l]->smoothu);
242: SLESGetKSP(mg[l]->smoothd,&ksp);
243: KSPSetTolerances(ksp,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,1);
244: SLESSetOptionsPrefix(mg[l]->smoothu,prefix);
245: SLESAppendOptionsPrefix(mg[l]->smoothd,"mg_levels_");
246: PetscLogObjectParent(pc,mg[l]->smoothu);
247: }
248: *sles = mg[l]->smoothu;
249: return(0);
250: }
252: #undef __FUNCT__
254: /*@C
255: MGGetSmootherDown - Gets the SLES context to be used as smoother before
256: coarse grid correction (pre-smoother).
258: Not Collective, SLES returned is parallel if PC is
260: Input Parameters:
261: + pc - the multigrid context
262: - l - the level to supply
264: Ouput Parameters:
265: . sles - the smoother
267: Level: advanced
269: .keywords: MG, multigrid, get, smoother, down, pre-smoother, level
271: .seealso: MGGetSmootherUp(), MGGetSmoother()
272: @*/
273: int MGGetSmootherDown(PC pc,int l,SLES *sles)
274: {
275: MG *mg = (MG*)pc->data;
278: *sles = mg[l]->smoothd;
279: return(0);
280: }
282: #undef __FUNCT__
284: /*@
285: MGSetCyclesOnLevel - Sets the number of cycles to run on this level.
287: Collective on PC
289: Input Parameters:
290: + pc - the multigrid context
291: . l - the level this is to be used for
292: - n - the number of cycles
294: Level: advanced
296: .keywords: MG, multigrid, set, cycles, V-cycle, W-cycle, level
298: .seealso: MGSetCycles()
299: @*/
300: int MGSetCyclesOnLevel(PC pc,int l,int c)
301: {
302: MG *mg = (MG*)pc->data;
305: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
306: mg[l]->cycles = c;
307: return(0);
308: }
310: #undef __FUNCT__
312: /*@
313: MGSetRhs - Sets the vector space to be used to store the right-hand side
314: on a particular level. The user should free this space at the conclusion
315: of multigrid use.
317: Collective on PC and Vec
319: Input Parameters:
320: + pc - the multigrid context
321: . l - the level this is to be used for
322: - c - the space
324: Level: advanced
326: .keywords: MG, multigrid, set, right-hand-side, rhs, level
328: .seealso: MGSetX(), MGSetR()
329: @*/
330: int MGSetRhs(PC pc,int l,Vec c)
331: {
332: MG *mg = (MG*)pc->data;
335: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
336: mg[l]->b = c;
337: return(0);
338: }
340: #undef __FUNCT__
342: /*@
343: MGSetX - Sets the vector space to be used to store the solution on a
344: particular level. The user should free this space at the conclusion
345: of multigrid use.
347: Collective on PC and Vec
349: Input Parameters:
350: + pc - the multigrid context
351: . l - the level this is to be used for
352: - c - the space
354: Level: advanced
356: .keywords: MG, multigrid, set, solution, level
358: .seealso: MGSetRhs(), MGSetR()
359: @*/
360: int MGSetX(PC pc,int l,Vec c)
361: {
362: MG *mg = (MG*)pc->data;
365: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
366: mg[l]->x = c;
367: return(0);
368: }
370: #undef __FUNCT__
372: /*@
373: MGSetR - Sets the vector space to be used to store the residual on a
374: particular level. The user should free this space at the conclusion of
375: multigrid use.
377: Collective on PC and Vec
379: Input Parameters:
380: + pc - the multigrid context
381: . l - the level this is to be used for
382: - c - the space
384: Level: advanced
386: .keywords: MG, multigrid, set, residual, level
387: @*/
388: int MGSetR(PC pc,int l,Vec c)
389: {
390: MG *mg = (MG*)pc->data;
393: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
394: mg[l]->r = c;
395: return(0);
396: }