1 2 #include <petsc/private/dmdaimpl.h> /*I "petscdmda.h" I*/ 3 4 PetscErrorCode DMSetFromOptions_DA(PetscOptionItems *PetscOptionsObject,DM da) 5 { 6 PetscErrorCode ierr; 7 DM_DA *dd = (DM_DA*)da->data; 8 PetscInt refine = 0,dim = da->dim,maxnlevels = 100,refx[100],refy[100],refz[100],n,i; 9 PetscBool flg; 10 11 PetscFunctionBegin; 12 PetscValidHeaderSpecific(da,DM_CLASSID,1); 13 14 if (dd->M < 0) SETERRQ(PetscObjectComm((PetscObject)da),PETSC_ERR_ARG_OUTOFRANGE,"Dimension must be non-negative, call DMSetFromOptions() if you want to change the value at runtime"); 15 if (dd->N < 0) SETERRQ(PetscObjectComm((PetscObject)da),PETSC_ERR_ARG_OUTOFRANGE,"Dimension must be non-negative, call DMSetFromOptions() if you want to change the value at runtime"); 16 if (dd->P < 0) SETERRQ(PetscObjectComm((PetscObject)da),PETSC_ERR_ARG_OUTOFRANGE,"Dimension must be non-negative, call DMSetFromOptions() if you want to change the value at runtime"); 17 18 ierr = PetscOptionsHead(PetscOptionsObject,"DMDA Options");CHKERRQ(ierr); 19 ierr = PetscOptionsBoundedInt("-da_grid_x","Number of grid points in x direction","DMDASetSizes",dd->M,&dd->M,NULL,1);CHKERRQ(ierr); 20 ierr = PetscOptionsBoundedInt("-da_grid_y","Number of grid points in y direction","DMDASetSizes",dd->N,&dd->N,NULL,1);CHKERRQ(ierr); 21 ierr = PetscOptionsBoundedInt("-da_grid_z","Number of grid points in z direction","DMDASetSizes",dd->P,&dd->P,NULL,1);CHKERRQ(ierr); 22 23 ierr = PetscOptionsBoundedInt("-da_overlap","Decomposition overlap in all directions","DMDASetOverlap",dd->xol,&dd->xol,&flg,0);CHKERRQ(ierr); 24 if (flg) {ierr = DMDASetOverlap(da,dd->xol,dd->xol,dd->xol);CHKERRQ(ierr);} 25 ierr = PetscOptionsBoundedInt("-da_overlap_x","Decomposition overlap in x direction","DMDASetOverlap",dd->xol,&dd->xol,NULL,0);CHKERRQ(ierr); 26 if (dim > 1) {ierr = PetscOptionsBoundedInt("-da_overlap_y","Decomposition overlap in y direction","DMDASetOverlap",dd->yol,&dd->yol,NULL,0);CHKERRQ(ierr);} 27 if (dim > 2) {ierr = PetscOptionsBoundedInt("-da_overlap_z","Decomposition overlap in z direction","DMDASetOverlap",dd->zol,&dd->zol,NULL,0);CHKERRQ(ierr);} 28 29 ierr = PetscOptionsBoundedInt("-da_local_subdomains","","DMDASetNumLocalSubdomains",dd->Nsub,&dd->Nsub,&flg,PETSC_DECIDE);CHKERRQ(ierr); 30 if (flg) {ierr = DMDASetNumLocalSubDomains(da,dd->Nsub);CHKERRQ(ierr);} 31 32 /* Handle DMDA parallel distribution */ 33 ierr = PetscOptionsBoundedInt("-da_processors_x","Number of processors in x direction","DMDASetNumProcs",dd->m,&dd->m,NULL,PETSC_DECIDE);CHKERRQ(ierr); 34 if (dim > 1) {ierr = PetscOptionsBoundedInt("-da_processors_y","Number of processors in y direction","DMDASetNumProcs",dd->n,&dd->n,NULL,PETSC_DECIDE);CHKERRQ(ierr);} 35 if (dim > 2) {ierr = PetscOptionsBoundedInt("-da_processors_z","Number of processors in z direction","DMDASetNumProcs",dd->p,&dd->p,NULL,PETSC_DECIDE);CHKERRQ(ierr);} 36 /* Handle DMDA refinement */ 37 ierr = PetscOptionsBoundedInt("-da_refine_x","Refinement ratio in x direction","DMDASetRefinementFactor",dd->refine_x,&dd->refine_x,NULL,1);CHKERRQ(ierr); 38 if (dim > 1) {ierr = PetscOptionsBoundedInt("-da_refine_y","Refinement ratio in y direction","DMDASetRefinementFactor",dd->refine_y,&dd->refine_y,NULL,1);CHKERRQ(ierr);} 39 if (dim > 2) {ierr = PetscOptionsBoundedInt("-da_refine_z","Refinement ratio in z direction","DMDASetRefinementFactor",dd->refine_z,&dd->refine_z,NULL,1);CHKERRQ(ierr);} 40 dd->coarsen_x = dd->refine_x; dd->coarsen_y = dd->refine_y; dd->coarsen_z = dd->refine_z; 41 42 /* Get refinement factors, defaults taken from the coarse DMDA */ 43 ierr = DMDAGetRefinementFactor(da,&refx[0],&refy[0],&refz[0]);CHKERRQ(ierr); 44 for (i=1; i<maxnlevels; i++) { 45 refx[i] = refx[0]; 46 refy[i] = refy[0]; 47 refz[i] = refz[0]; 48 } 49 n = maxnlevels; 50 ierr = PetscOptionsIntArray("-da_refine_hierarchy_x","Refinement factor for each level","None",refx,&n,&flg);CHKERRQ(ierr); 51 if (flg) { 52 dd->refine_x = refx[0]; 53 dd->refine_x_hier_n = n; 54 ierr = PetscMalloc1(n,&dd->refine_x_hier);CHKERRQ(ierr); 55 ierr = PetscArraycpy(dd->refine_x_hier,refx,n);CHKERRQ(ierr); 56 } 57 if (dim > 1) { 58 n = maxnlevels; 59 ierr = PetscOptionsIntArray("-da_refine_hierarchy_y","Refinement factor for each level","None",refy,&n,&flg);CHKERRQ(ierr); 60 if (flg) { 61 dd->refine_y = refy[0]; 62 dd->refine_y_hier_n = n; 63 ierr = PetscMalloc1(n,&dd->refine_y_hier);CHKERRQ(ierr); 64 ierr = PetscArraycpy(dd->refine_y_hier,refy,n);CHKERRQ(ierr); 65 } 66 } 67 if (dim > 2) { 68 n = maxnlevels; 69 ierr = PetscOptionsIntArray("-da_refine_hierarchy_z","Refinement factor for each level","None",refz,&n,&flg);CHKERRQ(ierr); 70 if (flg) { 71 dd->refine_z = refz[0]; 72 dd->refine_z_hier_n = n; 73 ierr = PetscMalloc1(n,&dd->refine_z_hier);CHKERRQ(ierr); 74 ierr = PetscArraycpy(dd->refine_z_hier,refz,n);CHKERRQ(ierr); 75 } 76 } 77 78 ierr = PetscOptionsBoundedInt("-da_refine","Uniformly refine DA one or more times","None",refine,&refine,NULL,0);CHKERRQ(ierr); 79 ierr = PetscOptionsTail();CHKERRQ(ierr); 80 81 while (refine--) { 82 if (dd->bx == DM_BOUNDARY_PERIODIC || dd->interptype == DMDA_Q0) { 83 dd->M = dd->refine_x*dd->M; 84 } else { 85 dd->M = 1 + dd->refine_x*(dd->M - 1); 86 } 87 if (dd->by == DM_BOUNDARY_PERIODIC || dd->interptype == DMDA_Q0) { 88 dd->N = dd->refine_y*dd->N; 89 } else { 90 dd->N = 1 + dd->refine_y*(dd->N - 1); 91 } 92 if (dd->bz == DM_BOUNDARY_PERIODIC || dd->interptype == DMDA_Q0) { 93 dd->P = dd->refine_z*dd->P; 94 } else { 95 dd->P = 1 + dd->refine_z*(dd->P - 1); 96 } 97 da->levelup++; 98 if (da->levelup - da->leveldown >= 0) { 99 dd->refine_x = refx[da->levelup - da->leveldown]; 100 dd->refine_y = refy[da->levelup - da->leveldown]; 101 dd->refine_z = refz[da->levelup - da->leveldown]; 102 } 103 if (da->levelup - da->leveldown >= 1) { 104 dd->coarsen_x = refx[da->levelup - da->leveldown - 1]; 105 dd->coarsen_y = refy[da->levelup - da->leveldown - 1]; 106 dd->coarsen_z = refz[da->levelup - da->leveldown - 1]; 107 } 108 } 109 PetscFunctionReturn(0); 110 } 111 112 extern PetscErrorCode DMCreateGlobalVector_DA(DM,Vec*); 113 extern PetscErrorCode DMCreateLocalVector_DA(DM,Vec*); 114 extern PetscErrorCode DMGlobalToLocalBegin_DA(DM,Vec,InsertMode,Vec); 115 extern PetscErrorCode DMGlobalToLocalEnd_DA(DM,Vec,InsertMode,Vec); 116 extern PetscErrorCode DMLocalToGlobalBegin_DA(DM,Vec,InsertMode,Vec); 117 extern PetscErrorCode DMLocalToGlobalEnd_DA(DM,Vec,InsertMode,Vec); 118 extern PetscErrorCode DMLocalToLocalBegin_DA(DM,Vec,InsertMode,Vec); 119 extern PetscErrorCode DMLocalToLocalEnd_DA(DM,Vec,InsertMode,Vec); 120 extern PetscErrorCode DMCreateInterpolation_DA(DM,DM,Mat*,Vec*); 121 extern PetscErrorCode DMCreateColoring_DA(DM,ISColoringType,ISColoring*); 122 extern PetscErrorCode DMCreateMatrix_DA(DM,Mat*); 123 extern PetscErrorCode DMCreateCoordinateDM_DA(DM,DM*); 124 extern PetscErrorCode DMRefine_DA(DM,MPI_Comm,DM*); 125 extern PetscErrorCode DMCoarsen_DA(DM,MPI_Comm,DM*); 126 extern PetscErrorCode DMRefineHierarchy_DA(DM,PetscInt,DM[]); 127 extern PetscErrorCode DMCoarsenHierarchy_DA(DM,PetscInt,DM[]); 128 extern PetscErrorCode DMCreateInjection_DA(DM,DM,Mat*); 129 extern PetscErrorCode DMCreateAggregates_DA(DM,DM,Mat*); 130 extern PetscErrorCode DMView_DA(DM,PetscViewer); 131 extern PetscErrorCode DMSetUp_DA(DM); 132 extern PetscErrorCode DMDestroy_DA(DM); 133 extern PetscErrorCode DMCreateDomainDecomposition_DA(DM,PetscInt*,char***,IS**,IS**,DM**); 134 extern PetscErrorCode DMCreateDomainDecompositionScatters_DA(DM,PetscInt,DM*,VecScatter**,VecScatter**,VecScatter**); 135 PETSC_INTERN PetscErrorCode DMGetCompatibility_DA(DM,DM,PetscBool*,PetscBool*); 136 137 PetscErrorCode DMLoad_DA(DM da,PetscViewer viewer) 138 { 139 PetscErrorCode ierr; 140 PetscInt dim,m,n,p,dof,swidth; 141 DMDAStencilType stencil; 142 DMBoundaryType bx,by,bz; 143 PetscBool coors; 144 DM dac; 145 Vec c; 146 147 PetscFunctionBegin; 148 ierr = PetscViewerBinaryRead(viewer,&dim,1,NULL,PETSC_INT);CHKERRQ(ierr); 149 ierr = PetscViewerBinaryRead(viewer,&m,1,NULL,PETSC_INT);CHKERRQ(ierr); 150 ierr = PetscViewerBinaryRead(viewer,&n,1,NULL,PETSC_INT);CHKERRQ(ierr); 151 ierr = PetscViewerBinaryRead(viewer,&p,1,NULL,PETSC_INT);CHKERRQ(ierr); 152 ierr = PetscViewerBinaryRead(viewer,&dof,1,NULL,PETSC_INT);CHKERRQ(ierr); 153 ierr = PetscViewerBinaryRead(viewer,&swidth,1,NULL,PETSC_INT);CHKERRQ(ierr); 154 ierr = PetscViewerBinaryRead(viewer,&bx,1,NULL,PETSC_ENUM);CHKERRQ(ierr); 155 ierr = PetscViewerBinaryRead(viewer,&by,1,NULL,PETSC_ENUM);CHKERRQ(ierr); 156 ierr = PetscViewerBinaryRead(viewer,&bz,1,NULL,PETSC_ENUM);CHKERRQ(ierr); 157 ierr = PetscViewerBinaryRead(viewer,&stencil,1,NULL,PETSC_ENUM);CHKERRQ(ierr); 158 159 ierr = DMSetDimension(da, dim);CHKERRQ(ierr); 160 ierr = DMDASetSizes(da, m,n,p);CHKERRQ(ierr); 161 ierr = DMDASetBoundaryType(da, bx, by, bz);CHKERRQ(ierr); 162 ierr = DMDASetDof(da, dof);CHKERRQ(ierr); 163 ierr = DMDASetStencilType(da, stencil);CHKERRQ(ierr); 164 ierr = DMDASetStencilWidth(da, swidth);CHKERRQ(ierr); 165 ierr = DMSetUp(da);CHKERRQ(ierr); 166 ierr = PetscViewerBinaryRead(viewer,&coors,1,NULL,PETSC_ENUM);CHKERRQ(ierr); 167 if (coors) { 168 ierr = DMGetCoordinateDM(da,&dac);CHKERRQ(ierr); 169 ierr = DMCreateGlobalVector(dac,&c);CHKERRQ(ierr); 170 ierr = VecLoad(c,viewer);CHKERRQ(ierr); 171 ierr = DMSetCoordinates(da,c);CHKERRQ(ierr); 172 ierr = VecDestroy(&c);CHKERRQ(ierr); 173 } 174 PetscFunctionReturn(0); 175 } 176 177 PetscErrorCode DMCreateSubDM_DA(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm) 178 { 179 DM_DA *da = (DM_DA*) dm->data; 180 PetscErrorCode ierr; 181 182 PetscFunctionBegin; 183 if (subdm) { 184 PetscSF sf; 185 Vec coords; 186 void *ctx; 187 /* Cannot use DMClone since the dof stuff is mixed in. Ugh 188 ierr = DMClone(dm, subdm);CHKERRQ(ierr); */ 189 ierr = DMCreate(PetscObjectComm((PetscObject)dm), subdm);CHKERRQ(ierr); 190 ierr = DMGetPointSF(dm, &sf);CHKERRQ(ierr); 191 ierr = DMSetPointSF(*subdm, sf);CHKERRQ(ierr); 192 ierr = DMGetApplicationContext(dm, &ctx);CHKERRQ(ierr); 193 ierr = DMSetApplicationContext(*subdm, ctx);CHKERRQ(ierr); 194 ierr = DMGetCoordinatesLocal(dm, &coords);CHKERRQ(ierr); 195 if (coords) { 196 ierr = DMSetCoordinatesLocal(*subdm, coords);CHKERRQ(ierr); 197 } else { 198 ierr = DMGetCoordinates(dm, &coords);CHKERRQ(ierr); 199 if (coords) {ierr = DMSetCoordinates(*subdm, coords);CHKERRQ(ierr);} 200 } 201 202 ierr = DMSetType(*subdm, DMDA);CHKERRQ(ierr); 203 ierr = DMSetDimension(*subdm, dm->dim);CHKERRQ(ierr); 204 ierr = DMDASetSizes(*subdm, da->M, da->N, da->P);CHKERRQ(ierr); 205 ierr = DMDASetNumProcs(*subdm, da->m, da->n, da->p);CHKERRQ(ierr); 206 ierr = DMDASetBoundaryType(*subdm, da->bx, da->by, da->bz);CHKERRQ(ierr); 207 ierr = DMDASetDof(*subdm, numFields);CHKERRQ(ierr); 208 ierr = DMDASetStencilType(*subdm, da->stencil_type);CHKERRQ(ierr); 209 ierr = DMDASetStencilWidth(*subdm, da->s);CHKERRQ(ierr); 210 ierr = DMDASetOwnershipRanges(*subdm, da->lx, da->ly, da->lz);CHKERRQ(ierr); 211 } 212 if (is) { 213 PetscInt *indices, cnt = 0, dof = da->w, i, j; 214 215 ierr = PetscMalloc1(da->Nlocal*numFields/dof, &indices);CHKERRQ(ierr); 216 for (i = da->base/dof; i < (da->base+da->Nlocal)/dof; ++i) { 217 for (j = 0; j < numFields; ++j) { 218 indices[cnt++] = dof*i + fields[j]; 219 } 220 } 221 if (cnt != da->Nlocal*numFields/dof) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Count %D does not equal expected value %D", cnt, da->Nlocal*numFields/dof); 222 ierr = ISCreateGeneral(PetscObjectComm((PetscObject) dm), cnt, indices, PETSC_OWN_POINTER, is);CHKERRQ(ierr); 223 } 224 PetscFunctionReturn(0); 225 } 226 227 PetscErrorCode DMCreateFieldDecomposition_DA(DM dm, PetscInt *len,char ***namelist, IS **islist, DM **dmlist) 228 { 229 PetscInt i; 230 PetscErrorCode ierr; 231 DM_DA *dd = (DM_DA*)dm->data; 232 PetscInt dof = dd->w; 233 234 PetscFunctionBegin; 235 if (len) *len = dof; 236 if (islist) { 237 Vec v; 238 PetscInt rstart,n; 239 240 ierr = DMGetGlobalVector(dm,&v);CHKERRQ(ierr); 241 ierr = VecGetOwnershipRange(v,&rstart,NULL);CHKERRQ(ierr); 242 ierr = VecGetLocalSize(v,&n);CHKERRQ(ierr); 243 ierr = DMRestoreGlobalVector(dm,&v);CHKERRQ(ierr); 244 ierr = PetscMalloc1(dof,islist);CHKERRQ(ierr); 245 for (i=0; i<dof; i++) { 246 ierr = ISCreateStride(PetscObjectComm((PetscObject)dm),n/dof,rstart+i,dof,&(*islist)[i]);CHKERRQ(ierr); 247 } 248 } 249 if (namelist) { 250 ierr = PetscMalloc1(dof, namelist);CHKERRQ(ierr); 251 if (dd->fieldname) { 252 for (i=0; i<dof; i++) { 253 ierr = PetscStrallocpy(dd->fieldname[i],&(*namelist)[i]);CHKERRQ(ierr); 254 } 255 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Currently DMDA must have fieldnames"); 256 } 257 if (dmlist) { 258 DM da; 259 260 ierr = DMDACreate(PetscObjectComm((PetscObject)dm), &da);CHKERRQ(ierr); 261 ierr = DMSetDimension(da, dm->dim);CHKERRQ(ierr); 262 ierr = DMDASetSizes(da, dd->M, dd->N, dd->P);CHKERRQ(ierr); 263 ierr = DMDASetNumProcs(da, dd->m, dd->n, dd->p);CHKERRQ(ierr); 264 ierr = DMDASetBoundaryType(da, dd->bx, dd->by, dd->bz);CHKERRQ(ierr); 265 ierr = DMDASetDof(da, 1);CHKERRQ(ierr); 266 ierr = DMDASetStencilType(da, dd->stencil_type);CHKERRQ(ierr); 267 ierr = DMDASetStencilWidth(da, dd->s);CHKERRQ(ierr); 268 ierr = DMSetUp(da);CHKERRQ(ierr); 269 ierr = PetscMalloc1(dof,dmlist);CHKERRQ(ierr); 270 for (i=0; i<dof-1; i++) {ierr = PetscObjectReference((PetscObject)da);CHKERRQ(ierr);} 271 for (i=0; i<dof; i++) (*dmlist)[i] = da; 272 } 273 PetscFunctionReturn(0); 274 } 275 276 PetscErrorCode DMClone_DA(DM dm, DM *newdm) 277 { 278 DM_DA *da = (DM_DA *) dm->data; 279 PetscErrorCode ierr; 280 281 PetscFunctionBegin; 282 ierr = DMSetType(*newdm, DMDA);CHKERRQ(ierr); 283 ierr = DMSetDimension(*newdm, dm->dim);CHKERRQ(ierr); 284 ierr = DMDASetSizes(*newdm, da->M, da->N, da->P);CHKERRQ(ierr); 285 ierr = DMDASetNumProcs(*newdm, da->m, da->n, da->p);CHKERRQ(ierr); 286 ierr = DMDASetBoundaryType(*newdm, da->bx, da->by, da->bz);CHKERRQ(ierr); 287 ierr = DMDASetDof(*newdm, da->w);CHKERRQ(ierr); 288 ierr = DMDASetStencilType(*newdm, da->stencil_type);CHKERRQ(ierr); 289 ierr = DMDASetStencilWidth(*newdm, da->s);CHKERRQ(ierr); 290 ierr = DMDASetOwnershipRanges(*newdm, da->lx, da->ly, da->lz);CHKERRQ(ierr); 291 ierr = DMSetUp(*newdm);CHKERRQ(ierr); 292 PetscFunctionReturn(0); 293 } 294 295 static PetscErrorCode DMHasCreateInjection_DA(DM dm, PetscBool *flg) 296 { 297 DM_DA *da = (DM_DA *)dm->data; 298 299 PetscFunctionBegin; 300 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 301 PetscValidPointer(flg,2); 302 *flg = da->interptype == DMDA_Q1 ? PETSC_TRUE : PETSC_FALSE; 303 PetscFunctionReturn(0); 304 } 305 306 static PetscErrorCode DMGetDimPoints_DA(DM dm, PetscInt dim, PetscInt *pStart, PetscInt *pEnd) 307 { 308 PetscErrorCode ierr; 309 310 PetscFunctionBegin; 311 ierr = DMDAGetDepthStratum(dm, dim, pStart, pEnd);CHKERRQ(ierr); 312 PetscFunctionReturn(0); 313 } 314 315 static PetscErrorCode DMGetNeighbors_DA(DM dm, PetscInt *nranks, const PetscMPIInt *ranks[]) 316 { 317 PetscErrorCode ierr; 318 PetscInt dim; 319 DMDAStencilType st; 320 321 PetscFunctionBegin; 322 ierr = DMDAGetNeighbors(dm,ranks);CHKERRQ(ierr); 323 ierr = DMDAGetInfo(dm,&dim,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,&st);CHKERRQ(ierr); 324 325 switch (dim) { 326 case 1: 327 *nranks = 3; 328 /* if (st == DMDA_STENCIL_STAR) { *nranks = 3; } */ 329 break; 330 case 2: 331 *nranks = 9; 332 /* if (st == DMDA_STENCIL_STAR) { *nranks = 5; } */ 333 break; 334 case 3: 335 *nranks = 27; 336 /* if (st == DMDA_STENCIL_STAR) { *nranks = 7; } */ 337 break; 338 default: 339 break; 340 } 341 PetscFunctionReturn(0); 342 } 343 344 /*MC 345 DMDA = "da" - A DM object that is used to manage data for a structured grid in 1, 2, or 3 dimensions. 346 In the global representation of the vector each process stores a non-overlapping rectangular (or slab in 3d) portion of the grid points. 347 In the local representation these rectangular regions (slabs) are extended in all directions by a stencil width. 348 349 The vectors can be thought of as either cell centered or vertex centered on the mesh. But some variables cannot be cell centered and others 350 vertex centered. 351 352 Level: intermediate 353 354 .seealso: DMType, DMCOMPOSITE, DMDACreate(), DMCreate(), DMSetType() 355 M*/ 356 357 extern PetscErrorCode DMLocatePoints_DA_Regular(DM,Vec,DMPointLocationType,PetscSF); 358 PETSC_INTERN PetscErrorCode DMSetUpGLVisViewer_DMDA(PetscObject,PetscViewer); 359 360 PETSC_EXTERN PetscErrorCode DMCreate_DA(DM da) 361 { 362 PetscErrorCode ierr; 363 DM_DA *dd; 364 365 PetscFunctionBegin; 366 PetscValidPointer(da,1); 367 ierr = PetscNewLog(da,&dd);CHKERRQ(ierr); 368 da->data = dd; 369 370 da->dim = -1; 371 dd->interptype = DMDA_Q1; 372 dd->refine_x = 2; 373 dd->refine_y = 2; 374 dd->refine_z = 2; 375 dd->coarsen_x = 2; 376 dd->coarsen_y = 2; 377 dd->coarsen_z = 2; 378 dd->fieldname = NULL; 379 dd->nlocal = -1; 380 dd->Nlocal = -1; 381 dd->M = -1; 382 dd->N = -1; 383 dd->P = -1; 384 dd->m = -1; 385 dd->n = -1; 386 dd->p = -1; 387 dd->w = -1; 388 dd->s = -1; 389 390 dd->xs = -1; dd->xe = -1; dd->ys = -1; dd->ye = -1; dd->zs = -1; dd->ze = -1; 391 dd->Xs = -1; dd->Xe = -1; dd->Ys = -1; dd->Ye = -1; dd->Zs = -1; dd->Ze = -1; 392 393 dd->Nsub = 1; 394 dd->xol = 0; 395 dd->yol = 0; 396 dd->zol = 0; 397 dd->xo = 0; 398 dd->yo = 0; 399 dd->zo = 0; 400 dd->Mo = -1; 401 dd->No = -1; 402 dd->Po = -1; 403 404 dd->gtol = NULL; 405 dd->ltol = NULL; 406 dd->ao = NULL; 407 PetscStrallocpy(AOBASIC,(char**)&dd->aotype); 408 dd->base = -1; 409 dd->bx = DM_BOUNDARY_NONE; 410 dd->by = DM_BOUNDARY_NONE; 411 dd->bz = DM_BOUNDARY_NONE; 412 dd->stencil_type = DMDA_STENCIL_BOX; 413 dd->interptype = DMDA_Q1; 414 dd->lx = NULL; 415 dd->ly = NULL; 416 dd->lz = NULL; 417 418 dd->elementtype = DMDA_ELEMENT_Q1; 419 420 da->ops->globaltolocalbegin = DMGlobalToLocalBegin_DA; 421 da->ops->globaltolocalend = DMGlobalToLocalEnd_DA; 422 da->ops->localtoglobalbegin = DMLocalToGlobalBegin_DA; 423 da->ops->localtoglobalend = DMLocalToGlobalEnd_DA; 424 da->ops->localtolocalbegin = DMLocalToLocalBegin_DA; 425 da->ops->localtolocalend = DMLocalToLocalEnd_DA; 426 da->ops->createglobalvector = DMCreateGlobalVector_DA; 427 da->ops->createlocalvector = DMCreateLocalVector_DA; 428 da->ops->createinterpolation = DMCreateInterpolation_DA; 429 da->ops->getcoloring = DMCreateColoring_DA; 430 da->ops->creatematrix = DMCreateMatrix_DA; 431 da->ops->refine = DMRefine_DA; 432 da->ops->coarsen = DMCoarsen_DA; 433 da->ops->refinehierarchy = DMRefineHierarchy_DA; 434 da->ops->coarsenhierarchy = DMCoarsenHierarchy_DA; 435 da->ops->getinjection = DMCreateInjection_DA; 436 da->ops->hascreateinjection = DMHasCreateInjection_DA; 437 da->ops->getaggregates = DMCreateAggregates_DA; 438 da->ops->destroy = DMDestroy_DA; 439 da->ops->view = 0; 440 da->ops->setfromoptions = DMSetFromOptions_DA; 441 da->ops->setup = DMSetUp_DA; 442 da->ops->clone = DMClone_DA; 443 da->ops->load = DMLoad_DA; 444 da->ops->createcoordinatedm = DMCreateCoordinateDM_DA; 445 da->ops->createsubdm = DMCreateSubDM_DA; 446 da->ops->createfielddecomposition = DMCreateFieldDecomposition_DA; 447 da->ops->createdomaindecomposition = DMCreateDomainDecomposition_DA; 448 da->ops->createddscatters = DMCreateDomainDecompositionScatters_DA; 449 da->ops->getdimpoints = DMGetDimPoints_DA; 450 da->ops->getneighbors = DMGetNeighbors_DA; 451 da->ops->locatepoints = DMLocatePoints_DA_Regular; 452 da->ops->getcompatibility = DMGetCompatibility_DA; 453 ierr = PetscObjectComposeFunction((PetscObject)da,"DMSetUpGLVisViewer_C",DMSetUpGLVisViewer_DMDA);CHKERRQ(ierr); 454 PetscFunctionReturn(0); 455 } 456 457 /*@ 458 DMDACreate - Creates a DMDA object. 459 460 Collective 461 462 Input Parameter: 463 . comm - The communicator for the DMDA object 464 465 Output Parameter: 466 . da - The DMDA object 467 468 Level: advanced 469 470 Developers Note: 471 Since there exists DMDACreate1/2/3d() should this routine even exist? 472 473 .seealso: DMDASetSizes(), DMClone(), DMDACreate1d(), DMDACreate2d(), DMDACreate3d() 474 @*/ 475 PetscErrorCode DMDACreate(MPI_Comm comm, DM *da) 476 { 477 PetscErrorCode ierr; 478 479 PetscFunctionBegin; 480 PetscValidPointer(da,2); 481 ierr = DMCreate(comm,da);CHKERRQ(ierr); 482 ierr = DMSetType(*da,DMDA);CHKERRQ(ierr); 483 PetscFunctionReturn(0); 484 } 485 486 487