Actual source code: patchcreate.c

petsc-dev 2014-02-02
Report Typos and Errors
  1: #include <petsc-private/dmpatchimpl.h>   /*I      "petscdmpatch.h"   I*/
  2: #include <petscdmda.h>

  6: PetscErrorCode DMSetFromOptions_Patch(DM dm)
  7: {
  8:   /* DM_Patch      *mesh = (DM_Patch*) dm->data; */

 13:   PetscOptionsHead("DMPatch Options");
 14:   /* Handle associated vectors */
 15:   /* Handle viewing */
 16:   PetscOptionsTail();
 17:   return(0);
 18: }

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

 30: PetscErrorCode DMInitialize_Patch(DM dm)
 31: {
 33:   dm->ops->view                            = DMView_Patch;
 34:   dm->ops->setfromoptions                  = DMSetFromOptions_Patch;
 35:   dm->ops->setup                           = DMSetUp_Patch;
 36:   dm->ops->createglobalvector              = DMCreateGlobalVector_Patch;
 37:   dm->ops->createlocalvector               = DMCreateLocalVector_Patch;
 38:   dm->ops->getlocaltoglobalmapping         = NULL;
 39:   dm->ops->getlocaltoglobalmappingblock    = NULL;
 40:   dm->ops->createfieldis                   = NULL;
 41:   dm->ops->getcoloring                     = 0;
 42:   dm->ops->creatematrix                    = 0;
 43:   dm->ops->createinterpolation             = 0;
 44:   dm->ops->getaggregates                   = 0;
 45:   dm->ops->getinjection                    = 0;
 46:   dm->ops->refine                          = 0;
 47:   dm->ops->coarsen                         = 0;
 48:   dm->ops->refinehierarchy                 = 0;
 49:   dm->ops->coarsenhierarchy                = 0;
 50:   dm->ops->globaltolocalbegin              = NULL;
 51:   dm->ops->globaltolocalend                = NULL;
 52:   dm->ops->localtoglobalbegin              = NULL;
 53:   dm->ops->localtoglobalend                = NULL;
 54:   dm->ops->destroy                         = DMDestroy_Patch;
 55:   dm->ops->createsubdm                     = DMCreateSubDM_Patch;
 56:   return(0);
 57: }

 61: PETSC_EXTERN PetscErrorCode DMCreate_Patch(DM dm)
 62: {
 63:   DM_Patch       *mesh;

 68:   PetscNewLog(dm,&mesh);
 69:   dm->data = mesh;

 71:   mesh->refct       = 1;
 72:   mesh->dmCoarse    = NULL;
 73:   mesh->patchSize.i = 0;
 74:   mesh->patchSize.j = 0;
 75:   mesh->patchSize.k = 0;
 76:   mesh->patchSize.c = 0;

 78:   DMInitialize_Patch(dm);
 79:   return(0);
 80: }

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

 87:   Collective on MPI_Comm

 89:   Input Parameter:
 90: . comm - The communicator for the DMPatch object

 92:   Output Parameter:
 93: . mesh  - The DMPatch object

 95:   Level: beginner

 97: .keywords: DMPatch, create
 98: @*/
 99: PetscErrorCode DMPatchCreate(MPI_Comm comm, DM *mesh)
100: {

105:   DMCreate(comm, mesh);
106:   DMSetType(*mesh, DMPATCH);
107:   return(0);
108: }

112: PetscErrorCode DMPatchCreateGrid(MPI_Comm comm, PetscInt dim, MatStencil patchSize, MatStencil commSize, MatStencil gridSize, DM *dm)
113: {
114:   DM_Patch       *mesh;
115:   DM             da;
116:   PetscInt       dof = 1, width = 1;

120:   DMPatchCreate(comm, dm);
121:   mesh = (DM_Patch*) (*dm)->data;
122:   if (dim < 2) {
123:     gridSize.j  = 1;
124:     patchSize.j = 1;
125:   }
126:   if (dim < 3) {
127:     gridSize.k  = 1;
128:     patchSize.k = 1;
129:   }
130:   DMCreate(comm, &da);
131:   DMSetType(da, DMDA);
132:   DMDASetDim(da, dim);
133:   DMDASetSizes(da, gridSize.i, gridSize.j, gridSize.k);
134:   DMDASetBoundaryType(da, DMDA_BOUNDARY_NONE, DMDA_BOUNDARY_NONE, DMDA_BOUNDARY_NONE);
135:   DMDASetDof(da, dof);
136:   DMDASetStencilType(da, DMDA_STENCIL_BOX);
137:   DMDASetStencilWidth(da, width);

139:   mesh->dmCoarse = da;

141:   DMPatchSetPatchSize(*dm, patchSize);
142:   DMPatchSetCommSize(*dm, commSize);
143:   return(0);
144: }