Actual source code: dm.c

  1: #define PETSCDM_DLL
  2: 
 3:  #include src/dm/da/daimpl.h

  5: /*
  6:    Provides an interface for functionality needed by the DAMG routines.
  7:    Currently this interface is supported by the DA and VecPack 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 VecPack.

 12: */

 16: /*@C
 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:   (*dm->bops->destroy)((PetscObject)dm);
 35:   return(0);
 36: }

 40: /*@C
 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 (dm->bops->view) {
 60:     (*dm->bops->view)((PetscObject)dm,v);
 61:   }
 62:   return(0);
 63: }

 67: /*@C
 68:     DMCreateGlobalVector - Creates a global vector from a DA or VecPack 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: /*@C
 95:     DMGetInterpolation - Gets interpolation matrix between two DA or VecPack objects

 97:     Collective on DM

 99:     Input Parameter:
100: +   dm1 - the DM object
101: -   dm2 - the second, coarser DM object

103:     Output Parameter:
104: +  mat - the interpolation
105: -  vec - the scaling (optional)

107:     Level: developer

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

111: @*/
112: PetscErrorCode  DMGetInterpolation(DM dm1,DM dm2,Mat *mat,Vec *vec)
113: {

117:   (*dm1->ops->getinterpolation)(dm1,dm2,mat,vec);
118:   return(0);
119: }

123: /*@C
124:     DMGetInjection - Gets injection matrix between two DA or VecPack objects

126:     Collective on DM

128:     Input Parameter:
129: +   dm1 - the DM object
130: -   dm2 - the second, coarser DM object

132:     Output Parameter:
133: .   ctx - the injection

135:     Level: developer

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

139: @*/
140: PetscErrorCode  DMGetInjection(DM dm1,DM dm2,VecScatter *ctx)
141: {

145:   (*dm1->ops->getinjection)(dm1,dm2,ctx);
146:   return(0);
147: }

151: /*@C
152:     DMGetColoring - Gets coloring and empty Jacobian for a DA or VecPack

154:     Collective on DM

156:     Input Parameter:
157: +   dm - the DM object
158: -   ctype - IS_COLORING_GHOSTED or IS_COLORING_GLOBAL

160:     Output Parameter:
161: .   coloring - the coloring

163:     Level: developer

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

167: @*/
168: PetscErrorCode  DMGetColoring(DM dm,ISColoringType ctype,ISColoring *coloring)
169: {

173:   if (!dm->ops->getcoloring) SETERRQ(PETSC_ERR_SUP,"No coloring for this type of DM yet");
174:   (*dm->ops->getcoloring)(dm,ctype,coloring);
175:   return(0);
176: }

180: /*@C
181:     DMGetMatrix - Gets empty Jacobian for a DA or VecPack

183:     Collective on DM

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

190:     Output Parameter:
191: .   mat - the empty Jacobian 

193:     Level: developer

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

197: @*/
198: PetscErrorCode  DMGetMatrix(DM dm, MatType mtype,Mat *mat)
199: {

203:   (*dm->ops->getmatrix)(dm,mtype,mat);
204:   return(0);
205: }

209: /*@C
210:     DMRefine - Refines a DM object

212:     Collective on DM

214:     Input Parameter:
215: +   dm - the DM object
216: -   comm - the communicator to contain the new DM object (or PETSC_NULL)

218:     Output Parameter:
219: .   dmf - the refined DM

221:     Level: developer

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

225: @*/
226: PetscErrorCode  DMRefine(DM dm,MPI_Comm comm,DM *dmf)
227: {

231:   (*dm->ops->refine)(dm,comm,dmf);
232:   return(0);
233: }

237: /*@C
238:     DMCoarsen - Coarsens 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: .   dmc - the coarsened DM

249:     Level: developer

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

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

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

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

268:     Collective on DM

270:     Input Parameter:
271: +   dm - the DM object
272: -   nlevels - the number of levels of refinement

274:     Output Parameter:
275: .   dmf - the refined DM hierarchy

277:     Level: developer

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

281: @*/
282: PetscErrorCode  DMRefineHierarchy(DM dm,int nlevels,DM **dmf)
283: {

287:   (*dm->ops->refinehierarchy)(dm,nlevels,dmf);
288:   return(0);
289: }

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

296:     Collective on DM

298:     Input Parameter:
299: +   dm - the DM object
300: -   nlevels - the number of levels of coarsening

302:     Output Parameter:
303: .   dmc - the coarsened DM hierarchy

305:     Level: developer

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

309: @*/
310: PetscErrorCode  DMCoarsenHierarchy(DM dm, int nlevels, DM **dmc)
311: {

315:   (*dm->ops->coarsenhierarchy)(dm, nlevels, dmc);
316:   return(0);
317: }