Actual source code: patchcreate.c

  1: #include <petsc/private/dmpatchimpl.h>
  2: #include <petscdmda.h>

  4: PetscErrorCode DMSetFromOptions_Patch(PetscOptionItems *PetscOptionsObject,DM dm)
  5: {
  6:   /* DM_Patch      *mesh = (DM_Patch*) dm->data; */

  9:   PetscOptionsHead(PetscOptionsObject,"DMPatch Options");
 10:   /* Handle associated vectors */
 11:   /* Handle viewing */
 12:   PetscOptionsTail();
 13:   return 0;
 14: }

 16: /* External function declarations here */
 17: extern PetscErrorCode DMSetUp_Patch(DM dm);
 18: extern PetscErrorCode DMView_Patch(DM dm, PetscViewer viewer);
 19: extern PetscErrorCode DMCreateGlobalVector_Patch(DM dm, Vec *g);
 20: extern PetscErrorCode DMCreateLocalVector_Patch(DM dm, Vec *l);
 21: extern PetscErrorCode DMDestroy_Patch(DM dm);
 22: extern PetscErrorCode DMCreateSubDM_Patch(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm);

 24: PetscErrorCode DMInitialize_Patch(DM dm)
 25: {
 26:   dm->ops->view                            = DMView_Patch;
 27:   dm->ops->setfromoptions                  = DMSetFromOptions_Patch;
 28:   dm->ops->setup                           = DMSetUp_Patch;
 29:   dm->ops->createglobalvector              = DMCreateGlobalVector_Patch;
 30:   dm->ops->createlocalvector               = DMCreateLocalVector_Patch;
 31:   dm->ops->getlocaltoglobalmapping         = NULL;
 32:   dm->ops->createfieldis                   = NULL;
 33:   dm->ops->getcoloring                     = NULL;
 34:   dm->ops->creatematrix                    = NULL;
 35:   dm->ops->createinterpolation             = NULL;
 36:   dm->ops->createinjection                 = NULL;
 37:   dm->ops->refine                          = NULL;
 38:   dm->ops->coarsen                         = NULL;
 39:   dm->ops->refinehierarchy                 = NULL;
 40:   dm->ops->coarsenhierarchy                = NULL;
 41:   dm->ops->globaltolocalbegin              = NULL;
 42:   dm->ops->globaltolocalend                = NULL;
 43:   dm->ops->localtoglobalbegin              = NULL;
 44:   dm->ops->localtoglobalend                = NULL;
 45:   dm->ops->destroy                         = DMDestroy_Patch;
 46:   dm->ops->createsubdm                     = DMCreateSubDM_Patch;
 47:   return 0;
 48: }

 50: PETSC_EXTERN PetscErrorCode DMCreate_Patch(DM dm)
 51: {
 52:   DM_Patch       *mesh;

 55:   PetscNewLog(dm,&mesh);
 56:   dm->data = mesh;

 58:   mesh->refct       = 1;
 59:   mesh->dmCoarse    = NULL;
 60:   mesh->patchSize.i = 0;
 61:   mesh->patchSize.j = 0;
 62:   mesh->patchSize.k = 0;
 63:   mesh->patchSize.c = 0;

 65:   DMInitialize_Patch(dm);
 66:   return 0;
 67: }

 69: /*@
 70:   DMPatchCreate - Creates a DMPatch object, which is a collections of DMs called patches.

 72:   Collective

 74:   Input Parameter:
 75: . comm - The communicator for the DMPatch object

 77:   Output Parameter:
 78: . mesh  - The DMPatch object

 80:   Notes:

 82:   This code is incomplete and not used by other parts of PETSc.

 84:   Level: beginner

 86: .seealso: DMPatchZoom()

 88: @*/
 89: PetscErrorCode DMPatchCreate(MPI_Comm comm, DM *mesh)
 90: {
 92:   DMCreate(comm, mesh);
 93:   DMSetType(*mesh, DMPATCH);
 94:   return 0;
 95: }

 97: PetscErrorCode DMPatchCreateGrid(MPI_Comm comm, PetscInt dim, MatStencil patchSize, MatStencil commSize, MatStencil gridSize, DM *dm)
 98: {
 99:   DM_Patch       *mesh;
100:   DM             da;
101:   PetscInt       dof = 1, width = 1;

103:   DMPatchCreate(comm, dm);
104:   mesh = (DM_Patch*) (*dm)->data;
105:   if (dim < 2) {
106:     gridSize.j  = 1;
107:     patchSize.j = 1;
108:   }
109:   if (dim < 3) {
110:     gridSize.k  = 1;
111:     patchSize.k = 1;
112:   }
113:   DMCreate(comm, &da);
114:   DMSetType(da, DMDA);
115:   DMSetDimension(da, dim);
116:   DMDASetSizes(da, gridSize.i, gridSize.j, gridSize.k);
117:   DMDASetBoundaryType(da, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE);
118:   DMDASetDof(da, dof);
119:   DMDASetStencilType(da, DMDA_STENCIL_BOX);
120:   DMDASetStencilWidth(da, width);

122:   mesh->dmCoarse = da;

124:   DMPatchSetPatchSize(*dm, patchSize);
125:   DMPatchSetCommSize(*dm, commSize);
126:   return 0;
127: }