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