Actual source code: dm.c

  1: #define PETSCDM_DLL
  2: 
 3:  #include private/dmimpl.h

  5: /*
  6:    Provides an interface for functionality needed by the DAMG routines.
  7:    Currently this interface is supported by the DA and DMComposite objects
  8:   
  9:    Note: this is actually no such thing as a DM object, rather it is 
 10:    the common set of functions shared by DA and DMComposite.

 12: */

 16: /*@
 17:     DMDestroy - Destroys a vector packer or DA.

 19:     Collective on DM

 21:     Input Parameter:
 22: .   dm - the DM object to destroy

 24:     Level: developer

 26: .seealso DMView(), DMCreateGlobalVector(), DMGetInterpolation(), DMGetColoring(), DMGetMatrix()

 28: @*/
 29: PetscErrorCode  DMDestroy(DM dm)
 30: {

 34:   (*((PetscObject)dm)->bops->destroy)((PetscObject)dm);
 35:   return(0);
 36: }

 40: /*@
 41:     DMView - Views a vector packer or DA.

 43:     Collective on DM

 45:     Input Parameter:
 46: +   dm - the DM object to view
 47: -   v - the viewer

 49:     Level: developer

 51: .seealso DMDestroy(), DMCreateGlobalVector(), DMGetInterpolation(), DMGetColoring(), DMGetMatrix()

 53: @*/
 54: PetscErrorCode  DMView(DM dm,PetscViewer v)
 55: {

 59:   if (((PetscObject)dm)->bops->view) {
 60:     (*((PetscObject)dm)->bops->view)((PetscObject)dm,v);
 61:   }
 62:   return(0);
 63: }

 67: /*@
 68:     DMCreateGlobalVector - Creates a global vector from a DA or DMComposite object

 70:     Collective on DM

 72:     Input Parameter:
 73: .   dm - the DM object

 75:     Output Parameter:
 76: .   vec - the global vector

 78:     Level: developer

 80: .seealso DMDestroy(), DMView(), DMGetInterpolation(), DMGetColoring(), DMGetMatrix()

 82: @*/
 83: PetscErrorCode  DMCreateGlobalVector(DM dm,Vec *vec)
 84: {

 88:   (*dm->ops->createglobalvector)(dm,vec);
 89:   return(0);
 90: }

 94: /*@
 95:     DMCreateLocalVector - Creates a local vector from a DA or DMComposite object

 97:     Collective on DM

 99:     Input Parameter:
100: .   dm - the DM object

102:     Output Parameter:
103: .   vec - the local vector

105:     Level: developer

107: .seealso DMDestroy(), DMView(), DMGetInterpolation(), DMGetColoring(), DMGetMatrix()

109: @*/
110: PetscErrorCode  DMCreateLocalVector(DM dm,Vec *vec)
111: {

115:   (*dm->ops->createlocalvector)(dm,vec);
116:   return(0);
117: }

121: /*@
122:     DMGetInterpolation - Gets interpolation matrix between two DA or DMComposite objects

124:     Collective on DM

126:     Input Parameter:
127: +   dm1 - the DM object
128: -   dm2 - the second, finer DM object

130:     Output Parameter:
131: +  mat - the interpolation
132: -  vec - the scaling (optional)

134:     Level: developer

136: .seealso DMDestroy(), DMView(), DMCreateGlobalVector(), DMGetColoring(), DMGetMatrix()

138: @*/
139: PetscErrorCode  DMGetInterpolation(DM dm1,DM dm2,Mat *mat,Vec *vec)
140: {

144:   (*dm1->ops->getinterpolation)(dm1,dm2,mat,vec);
145:   return(0);
146: }

150: /*@
151:     DMGetInjection - Gets injection matrix between two DA or DMComposite objects

153:     Collective on DM

155:     Input Parameter:
156: +   dm1 - the DM object
157: -   dm2 - the second, finer DM object

159:     Output Parameter:
160: .   ctx - the injection

162:     Level: developer

164: .seealso DMDestroy(), DMView(), DMCreateGlobalVector(), DMGetColoring(), DMGetMatrix(), DMGetInterpolation()

166: @*/
167: PetscErrorCode  DMGetInjection(DM dm1,DM dm2,VecScatter *ctx)
168: {

172:   (*dm1->ops->getinjection)(dm1,dm2,ctx);
173:   return(0);
174: }

178: /*@
179:     DMGetColoring - Gets coloring for a DA or DMComposite

181:     Collective on DM

183:     Input Parameter:
184: +   dm - the DM object
185: .   ctype - IS_COLORING_GHOSTED or IS_COLORING_GLOBAL
186: -   matype - either MATAIJ or MATBAIJ

188:     Output Parameter:
189: .   coloring - the coloring

191:     Level: developer

193: .seealso DMDestroy(), DMView(), DMCreateGlobalVector(), DMGetInterpolation(), DMGetMatrix()

195: @*/
196: PetscErrorCode  DMGetColoring(DM dm,ISColoringType ctype,const MatType mtype,ISColoring *coloring)
197: {

201:   if (!dm->ops->getcoloring) SETERRQ(PETSC_ERR_SUP,"No coloring for this type of DM yet");
202:   (*dm->ops->getcoloring)(dm,ctype,mtype,coloring);
203:   return(0);
204: }

208: /*@C
209:     DMGetMatrix - Gets empty Jacobian for a DA or DMComposite

211:     Collective on DM

213:     Input Parameter:
214: +   dm - the DM object
215: -   mtype - Supported types are MATSEQAIJ, MATMPIAIJ, MATSEQBAIJ, MATMPIBAIJ, or
216:             any type which inherits from one of these (such as MATAIJ, MATLUSOL, etc.).

218:     Output Parameter:
219: .   mat - the empty Jacobian 

221:     Level: developer

223: .seealso DMDestroy(), DMView(), DMCreateGlobalVector(), DMGetInterpolation(), DMGetMatrix()

225: @*/
226: PetscErrorCode  DMGetMatrix(DM dm, const MatType mtype,Mat *mat)
227: {

231:   (*dm->ops->getmatrix)(dm,mtype,mat);
232:   return(0);
233: }

237: /*@
238:     DMRefine - Refines a DM object

240:     Collective on DM

242:     Input Parameter:
243: +   dm - the DM object
244: -   comm - the communicator to contain the new DM object (or PETSC_NULL)

246:     Output Parameter:
247: .   dmf - the refined DM

249:     Level: developer

251: .seealso DMCoarsen(), DMDestroy(), DMView(), DMCreateGlobalVector(), DMGetInterpolation()

253: @*/
254: PetscErrorCode  DMRefine(DM dm,MPI_Comm comm,DM *dmf)
255: {

259:   (*dm->ops->refine)(dm,comm,dmf);
260:   return(0);
261: }

265: /*@
266:     DMGlobalToLocalBegin - Begins updating local vectors from local vectors

268:     Collective on DM

270:     Input Parameters:
271: +   dm - the DM object
272: .   g - the global vector
273: .   mode - INSERT_VALUES or ADD_VALUES
274: -   l - the local vector


277:     Level: beginner

279: .seealso DMCoarsen(), DMDestroy(), DMView(), DMCreateGlobalVector(), DMGetInterpolation(), DMGlobalToLocalEnd(), DMLocalToGlobal()

281: @*/
282: PetscErrorCode  DMGlobalToLocalBegin(DM dm,Vec g,InsertMode mode,Vec l)
283: {

287:   (*dm->ops->globaltolocalbegin)(dm,g,mode,l);
288:   return(0);
289: }

293: /*@
294:     DMGlobalToLocalEnd - Ends updating local vectors from local vectors

296:     Collective on DM

298:     Input Parameters:
299: +   dm - the DM object
300: .   g - the global vector
301: .   mode - INSERT_VALUES or ADD_VALUES
302: -   l - the local vector


305:     Level: beginner

307: .seealso DMCoarsen(), DMDestroy(), DMView(), DMCreateGlobalVector(), DMGetInterpolation(), DMGlobalToLocalEnd(), DMLocalToGlobal()

309: @*/
310: PetscErrorCode  DMGlobalToLocalEnd(DM dm,Vec g,InsertMode mode,Vec l)
311: {

315:   (*dm->ops->globaltolocalend)(dm,g,mode,l);
316:   return(0);
317: }

321: /*@
322:     DMLocalToGlobal - updates global vectors from local vectors

324:     Collective on DM

326:     Input Parameters:
327: +   dm - the DM object
328: .   g - the global vector
329: .   mode - INSERT_VALUES or ADD_VALUES
330: -   l - the local vector


333:     Level: beginner

335: .seealso DMCoarsen(), DMDestroy(), DMView(), DMCreateGlobalVector(), DMGetInterpolation(), DMGlobalToLocalEnd(), DMGlobalToLocalBegin()

337: @*/
338: PetscErrorCode  DMLocalToGlobal(DM dm,Vec g,InsertMode mode,Vec l)
339: {

343:   (*dm->ops->localtoglobal)(dm,g,mode,l);
344:   return(0);
345: }

349: /*@
350:     DMCoarsen - Coarsens a DM object

352:     Collective on DM

354:     Input Parameter:
355: +   dm - the DM object
356: -   comm - the communicator to contain the new DM object (or PETSC_NULL)

358:     Output Parameter:
359: .   dmc - the coarsened DM

361:     Level: developer

363: .seealso DMRefine(), DMDestroy(), DMView(), DMCreateGlobalVector(), DMGetInterpolation()

365: @*/
366: PetscErrorCode  DMCoarsen(DM dm, MPI_Comm comm, DM *dmc)
367: {

371:   (*dm->ops->coarsen)(dm, comm, dmc);
372:   return(0);
373: }

377: /*@C
378:     DMRefineHierarchy - Refines a DM object, all levels at once

380:     Collective on DM

382:     Input Parameter:
383: +   dm - the DM object
384: -   nlevels - the number of levels of refinement

386:     Output Parameter:
387: .   dmf - the refined DM hierarchy

389:     Level: developer

391: .seealso DMCoarsenHierarchy(), DMDestroy(), DMView(), DMCreateGlobalVector(), DMGetInterpolation()

393: @*/
394: PetscErrorCode  DMRefineHierarchy(DM dm,PetscInt nlevels,DM dmf[])
395: {

399:   if (nlevels < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"nlevels cannot be negative");
400:   if (nlevels == 0) return(0);
401:   if (dm->ops->refinehierarchy) {
402:     (*dm->ops->refinehierarchy)(dm,nlevels,dmf);
403:   } else if (dm->ops->refine) {
404:     PetscInt i;

406:     DMRefine(dm,((PetscObject)dm)->comm,&dmf[0]);
407:     for (i=1; i<nlevels; i++) {
408:       DMRefine(dmf[i-1],((PetscObject)dm)->comm,&dmf[i]);
409:     }
410:   } else {
411:     SETERRQ(PETSC_ERR_SUP,"No RefineHierarchy for this DM yet");
412:   }
413:   return(0);
414: }

418: /*@C
419:     DMCoarsenHierarchy - Coarsens a DM object, all levels at once

421:     Collective on DM

423:     Input Parameter:
424: +   dm - the DM object
425: -   nlevels - the number of levels of coarsening

427:     Output Parameter:
428: .   dmc - the coarsened DM hierarchy

430:     Level: developer

432: .seealso DMRefineHierarchy(), DMDestroy(), DMView(), DMCreateGlobalVector(), DMGetInterpolation()

434: @*/
435: PetscErrorCode  DMCoarsenHierarchy(DM dm, PetscInt nlevels, DM dmc[])
436: {

440:   if (nlevels < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"nlevels cannot be negative");
441:   if (nlevels == 0) return(0);
443:   if (dm->ops->coarsenhierarchy) {
444:     (*dm->ops->coarsenhierarchy)(dm, nlevels, dmc);
445:   } else if (dm->ops->coarsen) {
446:     PetscInt i;

448:     DMCoarsen(dm,((PetscObject)dm)->comm,&dmc[0]);
449:     for (i=1; i<nlevels; i++) {
450:       DMCoarsen(dmc[i-1],((PetscObject)dm)->comm,&dmc[i]);
451:     }
452:   } else {
453:     SETERRQ(PETSC_ERR_SUP,"No CoarsenHierarchy for this DM yet");
454:   }
455:   return(0);
456: }

460: /*@
461:    DMGetAggregates - Gets the aggregates that map between 
462:    grids associated with two DMs.

464:    Collective on DM

466:    Input Parameters:
467: +  dmc - the coarse grid DM
468: -  dmf - the fine grid DM

470:    Output Parameters:
471: .  rest - the restriction matrix (transpose of the projection matrix)

473:    Level: intermediate

475: .keywords: interpolation, restriction, multigrid 

477: .seealso: DMRefine(), DMGetInjection(), DMGetInterpolation()
478: @*/
479: PetscErrorCode  DMGetAggregates(DM dmc, DM dmf, Mat *rest)
480: {

484:   (*dmc->ops->getaggregates)(dmc, dmf, rest);
485:   return(0);
486: }