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