Actual source code: dacreate.c

  2: #include <private/daimpl.h>    /*I   "petscdmda.h"   I*/

  6: PetscErrorCode  DMSetFromOptions_DA(DM da)
  7: {
  9:   PetscBool      flg;
 10:   char           typeName[256];
 11:   DM_DA          *dd = (DM_DA*)da->data;
 12:   PetscInt       refine = 0;
 13:   PetscBool      negativeMNP = PETSC_FALSE,bM = PETSC_FALSE,bN = PETSC_FALSE, bP = PETSC_FALSE;


 18:   if (dd->M < 0) {
 19:     dd->M       = -dd->M;
 20:     bM          = PETSC_TRUE;
 21:     negativeMNP = PETSC_TRUE;
 22:   }
 23:   if (dd->dim > 1 && dd->N < 0) {
 24:     dd->N       = -dd->N;
 25:     bN          = PETSC_TRUE;
 26:     negativeMNP = PETSC_TRUE;
 27:   }
 28:   if (dd->dim > 2 && dd->P < 0) {
 29:     dd->P       = -dd->P;
 30:     bP          = PETSC_TRUE;
 31:     negativeMNP = PETSC_TRUE;
 32:   }

 34:   PetscOptionsBegin(((PetscObject)da)->comm,((PetscObject)da)->prefix,"DMDA Options","DMDA");
 35:     if (bM) {PetscOptionsInt("-da_grid_x","Number of grid points in x direction","DMDASetSizes",dd->M,&dd->M,PETSC_NULL);}
 36:     if (bN) {PetscOptionsInt("-da_grid_y","Number of grid points in y direction","DMDASetSizes",dd->N,&dd->N,PETSC_NULL);}
 37:     if (bP) {PetscOptionsInt("-da_grid_z","Number of grid points in z direction","DMDASetSizes",dd->P,&dd->P,PETSC_NULL);}
 38:     /* Handle DMDA parallel distibution */
 39:     PetscOptionsInt("-da_processors_x","Number of processors in x direction","DMDASetNumProcs",dd->m,&dd->m,PETSC_NULL);
 40:     if (dd->dim > 1) {PetscOptionsInt("-da_processors_y","Number of processors in y direction","DMDASetNumProcs",dd->n,&dd->n,PETSC_NULL);}
 41:     if (dd->dim > 2) {PetscOptionsInt("-da_processors_z","Number of processors in z direction","DMDASetNumProcs",dd->p,&dd->p,PETSC_NULL);}
 42:     /* Handle DMDA refinement */
 43:     PetscOptionsInt("-da_refine_x","Refinement ratio in x direction","DMDASetRefinementFactor",dd->refine_x,&dd->refine_x,PETSC_NULL);
 44:     if (dd->dim > 1) {PetscOptionsInt("-da_refine_y","Refinement ratio in y direction","DMDASetRefinementFactor",dd->refine_y,&dd->refine_y,PETSC_NULL);}
 45:     if (dd->dim > 2) {PetscOptionsInt("-da_refine_z","Refinement ratio in z direction","DMDASetRefinementFactor",dd->refine_z,&dd->refine_z,PETSC_NULL);}

 47:     if (!VecRegisterAllCalled) {VecRegisterAll(PETSC_NULL);}
 48:     PetscOptionsList("-da_vec_type","Vector type used for created vectors","DMSetVecType",VecList,da->vectype,typeName,256,&flg);
 49:     if (flg) {
 50:       DMSetVecType(da,typeName);
 51:     }
 52:     if (negativeMNP) {PetscOptionsInt("-da_refine","Uniformly refine DA one or more times","None",refine,&refine,PETSC_NULL);}

 54:     /* process any options handlers added with PetscObjectAddOptionsHandler() */
 55:     PetscObjectProcessOptionsHandlers((PetscObject)da);
 56:   PetscOptionsEnd();

 58:   while (refine--) {
 59:     if (dd->bx == DMDA_BOUNDARY_PERIODIC || dd->interptype == DMDA_Q0){
 60:       dd->M = dd->refine_x*dd->M;
 61:     } else {
 62:       dd->M = 1 + dd->refine_x*(dd->M - 1);
 63:     }
 64:     if (dd->by == DMDA_BOUNDARY_PERIODIC || dd->interptype == DMDA_Q0){
 65:       dd->N = dd->refine_y*dd->N;
 66:     } else {
 67:       dd->N = 1 + dd->refine_y*(dd->N - 1);
 68:     }
 69:     if (dd->bz == DMDA_BOUNDARY_PERIODIC || dd->interptype == DMDA_Q0){
 70:       dd->P = dd->refine_z*dd->P;
 71:     } else {
 72:       dd->P = 1 + dd->refine_z*(dd->P - 1);
 73:     }
 74:     da->levelup++;
 75:   }
 76:   return(0);
 77: }


100: PetscErrorCode DMLoad_DA(DM da,PetscViewer viewer)
101: {
102:   PetscErrorCode   ierr;
103:   PetscInt         dim,m,n,p,dof,swidth;
104:   DMDAStencilType  stencil;
105:   DMDABoundaryType bx,by,bz;
106:   PetscInt         classid = DM_FILE_CLASSID,subclassid = DMDA_FILE_CLASSID;
107:   PetscBool        coors;
108:   DM               dac;
109:   Vec              c;

112:   PetscViewerBinaryRead(viewer,&classid,1,PETSC_INT);
113:   if (classid != DM_FILE_CLASSID) SETERRQ(((PetscObject)da)->comm,PETSC_ERR_ARG_WRONG,"Not DM next in file");
114:   PetscViewerBinaryRead(viewer,&subclassid,1,PETSC_INT);
115:   if (subclassid != DMDA_FILE_CLASSID) SETERRQ(((PetscObject)da)->comm,PETSC_ERR_ARG_WRONG,"Not DM DA next in file");
116:   PetscViewerBinaryRead(viewer,&dim,1,PETSC_INT);
117:   PetscViewerBinaryRead(viewer,&m,1,PETSC_INT);
118:   PetscViewerBinaryRead(viewer,&n,1,PETSC_INT);
119:   PetscViewerBinaryRead(viewer,&p,1,PETSC_INT);
120:   PetscViewerBinaryRead(viewer,&dof,1,PETSC_INT);
121:   PetscViewerBinaryRead(viewer,&swidth,1,PETSC_INT);
122:   PetscViewerBinaryRead(viewer,&bx,1,PETSC_ENUM);
123:   PetscViewerBinaryRead(viewer,&by,1,PETSC_ENUM);
124:   PetscViewerBinaryRead(viewer,&bz,1,PETSC_ENUM);
125:   PetscViewerBinaryRead(viewer,&stencil,1,PETSC_ENUM);

127:   DMDASetDim(da, dim);
128:   DMDASetSizes(da, m,n,p);
129:   DMDASetBoundaryType(da, bx, by, bz);
130:   DMDASetDof(da, dof);
131:   DMDASetStencilType(da, stencil);
132:   DMDASetStencilWidth(da, swidth);
133:   DMSetUp(da);
134:   PetscViewerBinaryRead(viewer,&coors,1,PETSC_ENUM);
135:   if (coors) {
136:     DMDAGetCoordinateDA(da,&dac);
137:     DMCreateGlobalVector(dac,&c);
138:     VecLoad(c,viewer);
139:     DMDASetCoordinates(da,c);
140:     VecDestroy(&c);
141:   }
142:   return(0);
143: }

148: PetscErrorCode  DMCreate_DA(DM da)
149: {
151:   DM_DA          *dd;

155:   PetscNewLog(da,DM_DA,&dd);
156:   da->data = dd;

158:   dd->dim        = -1;
159:   dd->interptype = DMDA_Q1;
160:   dd->refine_x   = 2;
161:   dd->refine_y   = 2;
162:   dd->refine_z   = 2;
163:   dd->fieldname  = PETSC_NULL;
164:   dd->nlocal     = -1;
165:   dd->Nlocal     = -1;
166:   dd->M          = -1;
167:   dd->N          = -1;
168:   dd->P          = -1;
169:   dd->m          = -1;
170:   dd->n          = -1;
171:   dd->p          = -1;
172:   dd->w          = -1;
173:   dd->s          = -1;
174:   dd->xs = -1; dd->xe = -1; dd->ys = -1; dd->ye = -1; dd->zs = -1; dd->ze = -1;
175:   dd->Xs = -1; dd->Xe = -1; dd->Ys = -1; dd->Ye = -1; dd->Zs = -1; dd->Ze = -1;

177:   dd->gtol         = PETSC_NULL;
178:   dd->ltog         = PETSC_NULL;
179:   dd->ltol         = PETSC_NULL;
180:   dd->ao           = PETSC_NULL;
181:   dd->base         = -1;
182:   dd->bx         = DMDA_BOUNDARY_NONE;
183:   dd->by         = DMDA_BOUNDARY_NONE;
184:   dd->bz         = DMDA_BOUNDARY_NONE;
185:   dd->stencil_type = DMDA_STENCIL_BOX;
186:   dd->interptype   = DMDA_Q1;
187:   dd->idx          = PETSC_NULL;
188:   dd->Nl           = -1;
189:   dd->lx           = PETSC_NULL;
190:   dd->ly           = PETSC_NULL;
191:   dd->lz           = PETSC_NULL;

193:   dd->elementtype  = DMDA_ELEMENT_Q1;

195:   PetscStrallocpy(VECSTANDARD,&da->vectype);
196:   da->ops->globaltolocalbegin = DMGlobalToLocalBegin_DA;
197:   da->ops->globaltolocalend   = DMGlobalToLocalEnd_DA;
198:   da->ops->localtoglobalbegin = DMLocalToGlobalBegin_DA;
199:   da->ops->localtoglobalend   = DMLocalToGlobalEnd_DA;
200:   da->ops->createglobalvector = DMCreateGlobalVector_DA;
201:   da->ops->createlocalvector  = DMCreateLocalVector_DA;
202:   da->ops->getinterpolation   = DMGetInterpolation_DA;
203:   da->ops->getcoloring        = DMGetColoring_DA;
204:   da->ops->getmatrix          = DMGetMatrix_DA;
205:   da->ops->refine             = DMRefine_DA;
206:   da->ops->coarsen            = DMCoarsen_DA;
207:   da->ops->refinehierarchy    = DMRefineHierarchy_DA;
208:   da->ops->coarsenhierarchy   = DMCoarsenHierarchy_DA;
209:   da->ops->getinjection       = DMGetInjection_DA;
210:   da->ops->getaggregates      = DMGetAggregates_DA;
211:   da->ops->destroy            = DMDestroy_DA;
212:   da->ops->view               = 0;
213:   da->ops->setfromoptions     = DMSetFromOptions_DA;
214:   da->ops->setup              = DMSetUp_DA;
215:   da->ops->load               = DMLoad_DA;
216:   return(0);
217: }

222: /*@
223:   DMDACreate - Creates a DMDA object. 

225:   Collective on MPI_Comm

227:   Input Parameter:
228: . comm - The communicator for the DMDA object

230:   Output Parameter:
231: . da  - The DMDA object

233:   Level: advanced

235:   Developers Note: Since there exists DMDACreate1/2/3d() should this routine even exist?

237: .keywords: DMDA, create
238: .seealso:  DMDASetSizes(), DMDADuplicate(),  DMDACreate1d(), DMDACreate2d(), DMDACreate3d()
239: @*/
240: PetscErrorCode  DMDACreate(MPI_Comm comm, DM *da)
241: {

246:   DMCreate(comm,da);
247:   DMSetType(*da,DMDA);
248:   return(0);
249: }