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

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

190:     Level: developer

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

194: @*/
195: PetscErrorCode  DMGetColoring(DM dm,ISColoringType ctype,ISColoring *coloring)
196: {

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

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

210:     Collective on DM

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

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

220:     Level: developer

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

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

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

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

239:     Collective on DM

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

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

248:     Level: developer

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

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

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

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

267:     Collective on DM

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


276:     Level: beginner

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

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

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

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

295:     Collective on DM

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


304:     Level: beginner

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

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

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

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

323:     Collective on DM

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


332:     Level: beginner

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

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

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

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

351:     Collective on DM

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

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

360:     Level: developer

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

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

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

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

379:     Collective on DM

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

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

388:     Level: developer

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

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

398:   (*dm->ops->refinehierarchy)(dm,nlevels,dmf);
399:   return(0);
400: }

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

407:     Collective on DM

409:     Input Parameter:
410: +   dm - the DM object
411: -   nlevels - the number of levels of coarsening

413:     Output Parameter:
414: .   dmc - the coarsened DM hierarchy

416:     Level: developer

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

420: @*/
421: PetscErrorCode  DMCoarsenHierarchy(DM dm, PetscInt nlevels, DM **dmc)
422: {

426:   (*dm->ops->coarsenhierarchy)(dm, nlevels, dmc);
427:   return(0);
428: }

432: /*@
433:    DMGetAggregates - Gets the aggregates that map between 
434:    grids associated with two DMs.

436:    Collective on DM

438:    Input Parameters:
439: +  dmc - the coarse grid DM
440: -  dmf - the fine grid DM

442:    Output Parameters:
443: .  rest - the restriction matrix (transpose of the projection matrix)

445:    Level: intermediate

447: .keywords: interpolation, restriction, multigrid 

449: .seealso: DMRefine(), DMGetInjection(), DMGetInterpolation()
450: @*/
451: PetscErrorCode  DMGetAggregates(DM dmc, DM dmf, Mat *rest)
452: {

456:   (*dmc->ops->getaggregates)(dmc, dmf, rest);
457:   return(0);
458: }