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 = PetscOptionsInt("-da_grid_x","Number of grid points in x direction","DMDASetSizes",dd->M,&dd->M,NULL);CHKERRQ(ierr); 20 ierr = PetscOptionsInt("-da_grid_y","Number of grid points in y direction","DMDASetSizes",dd->N,&dd->N,NULL);CHKERRQ(ierr); 21 ierr = PetscOptionsInt("-da_grid_z","Number of grid points in z direction","DMDASetSizes",dd->P,&dd->P,NULL);CHKERRQ(ierr); 22 23 ierr = PetscOptionsInt("-da_overlap","Decomposition overlap in all directions","DMDASetOverlap",dd->xol,&dd->xol,&flg);CHKERRQ(ierr); 24 if (flg) {ierr = DMDASetOverlap(da,dd->xol,dd->xol,dd->xol);CHKERRQ(ierr);} 25 ierr = PetscOptionsInt("-da_overlap_x","Decomposition overlap in x direction","DMDASetOverlap",dd->xol,&dd->xol,NULL);CHKERRQ(ierr); 26 if (dim > 1) {ierr = PetscOptionsInt("-da_overlap_y","Decomposition overlap in y direction","DMDASetOverlap",dd->yol,&dd->yol,NULL);CHKERRQ(ierr);} 27 if (dim > 2) {ierr = PetscOptionsInt("-da_overlap_z","Decomposition overlap in z direction","DMDASetOverlap",dd->zol,&dd->zol,NULL);CHKERRQ(ierr);} 28 29 ierr = PetscOptionsInt("-da_local_subdomains","","DMDASetNumLocalSubdomains",dd->Nsub,&dd->Nsub,&flg);CHKERRQ(ierr); 30 if (flg) {ierr = DMDASetNumLocalSubDomains(da,dd->Nsub);CHKERRQ(ierr);} 31 32 /* Handle DMDA parallel distribution */ 33 ierr = PetscOptionsInt("-da_processors_x","Number of processors in x direction","DMDASetNumProcs",dd->m,&dd->m,NULL);CHKERRQ(ierr); 34 if (dim > 1) {ierr = PetscOptionsInt("-da_processors_y","Number of processors in y direction","DMDASetNumProcs",dd->n,&dd->n,NULL);CHKERRQ(ierr);} 35 if (dim > 2) {ierr = PetscOptionsInt("-da_processors_z","Number of processors in z direction","DMDASetNumProcs",dd->p,&dd->p,NULL);CHKERRQ(ierr);} 36 /* Handle DMDA refinement */ 37 ierr = PetscOptionsInt("-da_refine_x","Refinement ratio in x direction","DMDASetRefinementFactor",dd->refine_x,&dd->refine_x,NULL);CHKERRQ(ierr); 38 if (dim > 1) {ierr = PetscOptionsInt("-da_refine_y","Refinement ratio in y direction","DMDASetRefinementFactor",dd->refine_y,&dd->refine_y,NULL);CHKERRQ(ierr);} 39 if (dim > 2) {ierr = PetscOptionsInt("-da_refine_z","Refinement ratio in z direction","DMDASetRefinementFactor",dd->refine_z,&dd->refine_z,NULL);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 = PetscMemcpy(dd->refine_x_hier,refx,n*sizeof(PetscInt));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 = PetscMemcpy(dd->refine_y_hier,refy,n*sizeof(PetscInt));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 = PetscMemcpy(dd->refine_z_hier,refz,n*sizeof(PetscInt));CHKERRQ(ierr); 75 } 76 } 77 78 ierr = PetscOptionsInt("-da_refine","Uniformly refine DA one or more times","None",refine,&refine,NULL);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 PetscSection section; 181 PetscErrorCode ierr; 182 183 PetscFunctionBegin; 184 if (subdm) { 185 PetscSF sf; 186 Vec coords; 187 void *ctx; 188 /* Cannot use DMClone since the dof stuff is mixed in. Ugh 189 ierr = DMClone(dm, subdm);CHKERRQ(ierr); */ 190 ierr = DMCreate(PetscObjectComm((PetscObject)dm), subdm);CHKERRQ(ierr); 191 ierr = DMGetPointSF(dm, &sf);CHKERRQ(ierr); 192 ierr = DMSetPointSF(*subdm, sf);CHKERRQ(ierr); 193 ierr = DMGetApplicationContext(dm, &ctx);CHKERRQ(ierr); 194 ierr = DMSetApplicationContext(*subdm, ctx);CHKERRQ(ierr); 195 ierr = DMGetCoordinatesLocal(dm, &coords);CHKERRQ(ierr); 196 if (coords) { 197 ierr = DMSetCoordinatesLocal(*subdm, coords);CHKERRQ(ierr); 198 } else { 199 ierr = DMGetCoordinates(dm, &coords);CHKERRQ(ierr); 200 if (coords) {ierr = DMSetCoordinates(*subdm, coords);CHKERRQ(ierr);} 201 } 202 203 ierr = DMSetType(*subdm, DMDA);CHKERRQ(ierr); 204 ierr = DMSetDimension(*subdm, dm->dim);CHKERRQ(ierr); 205 ierr = DMDASetSizes(*subdm, da->M, da->N, da->P);CHKERRQ(ierr); 206 ierr = DMDASetNumProcs(*subdm, da->m, da->n, da->p);CHKERRQ(ierr); 207 ierr = DMDASetBoundaryType(*subdm, da->bx, da->by, da->bz);CHKERRQ(ierr); 208 ierr = DMDASetDof(*subdm, numFields);CHKERRQ(ierr); 209 ierr = DMDASetStencilType(*subdm, da->stencil_type);CHKERRQ(ierr); 210 ierr = DMDASetStencilWidth(*subdm, da->s);CHKERRQ(ierr); 211 ierr = DMDASetOwnershipRanges(*subdm, da->lx, da->ly, da->lz);CHKERRQ(ierr); 212 } 213 ierr = DMGetSection(dm, §ion);CHKERRQ(ierr); 214 if (section) { 215 ierr = DMCreateSectionSubDM(dm, numFields, fields, is, subdm);CHKERRQ(ierr); 216 } else { 217 if (is) { 218 PetscInt *indices, cnt = 0, dof = da->w, i, j; 219 220 ierr = PetscMalloc1(da->Nlocal*numFields/dof, &indices);CHKERRQ(ierr); 221 for (i = da->base/dof; i < (da->base+da->Nlocal)/dof; ++i) { 222 for (j = 0; j < numFields; ++j) { 223 indices[cnt++] = dof*i + fields[j]; 224 } 225 } 226 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); 227 ierr = ISCreateGeneral(PetscObjectComm((PetscObject) dm), cnt, indices, PETSC_OWN_POINTER, is);CHKERRQ(ierr); 228 } 229 } 230 PetscFunctionReturn(0); 231 } 232 233 PetscErrorCode DMCreateFieldDecomposition_DA(DM dm, PetscInt *len,char ***namelist, IS **islist, DM **dmlist) 234 { 235 PetscInt i; 236 PetscErrorCode ierr; 237 DM_DA *dd = (DM_DA*)dm->data; 238 PetscInt dof = dd->w; 239 240 PetscFunctionBegin; 241 if (len) *len = dof; 242 if (islist) { 243 Vec v; 244 PetscInt rstart,n; 245 246 ierr = DMGetGlobalVector(dm,&v);CHKERRQ(ierr); 247 ierr = VecGetOwnershipRange(v,&rstart,NULL);CHKERRQ(ierr); 248 ierr = VecGetLocalSize(v,&n);CHKERRQ(ierr); 249 ierr = DMRestoreGlobalVector(dm,&v);CHKERRQ(ierr); 250 ierr = PetscMalloc1(dof,islist);CHKERRQ(ierr); 251 for (i=0; i<dof; i++) { 252 ierr = ISCreateStride(PetscObjectComm((PetscObject)dm),n/dof,rstart+i,dof,&(*islist)[i]);CHKERRQ(ierr); 253 } 254 } 255 if (namelist) { 256 ierr = PetscMalloc1(dof, namelist);CHKERRQ(ierr); 257 if (dd->fieldname) { 258 for (i=0; i<dof; i++) { 259 ierr = PetscStrallocpy(dd->fieldname[i],&(*namelist)[i]);CHKERRQ(ierr); 260 } 261 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Currently DMDA must have fieldnames"); 262 } 263 if (dmlist) { 264 DM da; 265 266 ierr = DMDACreate(PetscObjectComm((PetscObject)dm), &da);CHKERRQ(ierr); 267 ierr = DMSetDimension(da, dm->dim);CHKERRQ(ierr); 268 ierr = DMDASetSizes(da, dd->M, dd->N, dd->P);CHKERRQ(ierr); 269 ierr = DMDASetNumProcs(da, dd->m, dd->n, dd->p);CHKERRQ(ierr); 270 ierr = DMDASetBoundaryType(da, dd->bx, dd->by, dd->bz);CHKERRQ(ierr); 271 ierr = DMDASetDof(da, 1);CHKERRQ(ierr); 272 ierr = DMDASetStencilType(da, dd->stencil_type);CHKERRQ(ierr); 273 ierr = DMDASetStencilWidth(da, dd->s);CHKERRQ(ierr); 274 ierr = DMSetUp(da);CHKERRQ(ierr); 275 ierr = PetscMalloc1(dof,dmlist);CHKERRQ(ierr); 276 for (i=0; i<dof-1; i++) {ierr = PetscObjectReference((PetscObject)da);CHKERRQ(ierr);} 277 for (i=0; i<dof; i++) (*dmlist)[i] = da; 278 } 279 PetscFunctionReturn(0); 280 } 281 282 PetscErrorCode DMClone_DA(DM dm, DM *newdm) 283 { 284 DM_DA *da = (DM_DA *) dm->data; 285 PetscErrorCode ierr; 286 287 PetscFunctionBegin; 288 ierr = DMSetType(*newdm, DMDA);CHKERRQ(ierr); 289 ierr = DMSetDimension(*newdm, dm->dim);CHKERRQ(ierr); 290 ierr = DMDASetSizes(*newdm, da->M, da->N, da->P);CHKERRQ(ierr); 291 ierr = DMDASetNumProcs(*newdm, da->m, da->n, da->p);CHKERRQ(ierr); 292 ierr = DMDASetBoundaryType(*newdm, da->bx, da->by, da->bz);CHKERRQ(ierr); 293 ierr = DMDASetDof(*newdm, da->w);CHKERRQ(ierr); 294 ierr = DMDASetStencilType(*newdm, da->stencil_type);CHKERRQ(ierr); 295 ierr = DMDASetStencilWidth(*newdm, da->s);CHKERRQ(ierr); 296 ierr = DMDASetOwnershipRanges(*newdm, da->lx, da->ly, da->lz);CHKERRQ(ierr); 297 ierr = DMSetUp(*newdm);CHKERRQ(ierr); 298 PetscFunctionReturn(0); 299 } 300 301 static PetscErrorCode DMHasCreateInjection_DA(DM dm, PetscBool *flg) 302 { 303 DM_DA *da = (DM_DA *)dm->data; 304 305 PetscFunctionBegin; 306 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 307 PetscValidPointer(flg,2); 308 *flg = da->interptype == DMDA_Q1 ? PETSC_TRUE : PETSC_FALSE; 309 PetscFunctionReturn(0); 310 } 311 312 static PetscErrorCode DMGetDimPoints_DA(DM dm, PetscInt dim, PetscInt *pStart, PetscInt *pEnd) 313 { 314 PetscErrorCode ierr; 315 316 PetscFunctionBegin; 317 ierr = DMDAGetDepthStratum(dm, dim, pStart, pEnd);CHKERRQ(ierr); 318 PetscFunctionReturn(0); 319 } 320 321 static PetscErrorCode DMGetNeighbors_DA(DM dm, PetscInt *nranks, const PetscMPIInt *ranks[]) 322 { 323 PetscErrorCode ierr; 324 PetscInt dim; 325 DMDAStencilType st; 326 327 PetscFunctionBegin; 328 ierr = DMDAGetNeighbors(dm,ranks);CHKERRQ(ierr); 329 ierr = DMDAGetInfo(dm,&dim,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,&st);CHKERRQ(ierr); 330 331 switch (dim) { 332 case 1: 333 *nranks = 3; 334 /* if (st == DMDA_STENCIL_STAR) { *nranks = 3; } */ 335 break; 336 case 2: 337 *nranks = 9; 338 /* if (st == DMDA_STENCIL_STAR) { *nranks = 5; } */ 339 break; 340 case 3: 341 *nranks = 27; 342 /* if (st == DMDA_STENCIL_STAR) { *nranks = 7; } */ 343 break; 344 default: 345 break; 346 } 347 PetscFunctionReturn(0); 348 } 349 350 /*MC 351 DMDA = "da" - A DM object that is used to manage data for a structured grid in 1, 2, or 3 dimensions. 352 In the global representation of the vector each process stores a non-overlapping rectangular (or slab in 3d) portion of the grid points. 353 In the local representation these rectangular regions (slabs) are extended in all directions by a stencil width. 354 355 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 356 vertex centered. 357 358 Level: intermediate 359 360 .seealso: DMType, DMCOMPOSITE, DMDACreate(), DMCreate(), DMSetType() 361 M*/ 362 363 extern PetscErrorCode DMProjectFunctionLocal_DA(DM, PetscReal, PetscErrorCode (**)(PetscInt, PetscReal, const PetscReal [], PetscInt, PetscScalar *, void *), void **, InsertMode, Vec); 364 extern PetscErrorCode DMComputeL2Diff_DA(DM, PetscReal, PetscErrorCode (**)(PetscInt, PetscReal, const PetscReal [], PetscInt, PetscScalar *, void *), void **, Vec, PetscReal *); 365 extern PetscErrorCode DMComputeL2GradientDiff_DA(DM, PetscReal, PetscErrorCode (**)(PetscInt, PetscReal, const PetscReal [], const PetscReal [],PetscInt, PetscScalar *, void *), void **, Vec,const PetscReal [], PetscReal *); 366 extern PetscErrorCode DMLocatePoints_DA_Regular(DM,Vec,DMPointLocationType,PetscSF); 367 PETSC_INTERN PetscErrorCode DMSetUpGLVisViewer_DMDA(PetscObject,PetscViewer); 368 369 PETSC_EXTERN PetscErrorCode DMCreate_DA(DM da) 370 { 371 PetscErrorCode ierr; 372 DM_DA *dd; 373 374 PetscFunctionBegin; 375 PetscValidPointer(da,1); 376 ierr = PetscNewLog(da,&dd);CHKERRQ(ierr); 377 da->data = dd; 378 379 da->dim = -1; 380 dd->interptype = DMDA_Q1; 381 dd->refine_x = 2; 382 dd->refine_y = 2; 383 dd->refine_z = 2; 384 dd->coarsen_x = 2; 385 dd->coarsen_y = 2; 386 dd->coarsen_z = 2; 387 dd->fieldname = NULL; 388 dd->nlocal = -1; 389 dd->Nlocal = -1; 390 dd->M = -1; 391 dd->N = -1; 392 dd->P = -1; 393 dd->m = -1; 394 dd->n = -1; 395 dd->p = -1; 396 dd->w = -1; 397 dd->s = -1; 398 399 dd->xs = -1; dd->xe = -1; dd->ys = -1; dd->ye = -1; dd->zs = -1; dd->ze = -1; 400 dd->Xs = -1; dd->Xe = -1; dd->Ys = -1; dd->Ye = -1; dd->Zs = -1; dd->Ze = -1; 401 402 dd->Nsub = 1; 403 dd->xol = 0; 404 dd->yol = 0; 405 dd->zol = 0; 406 dd->xo = 0; 407 dd->yo = 0; 408 dd->zo = 0; 409 dd->Mo = -1; 410 dd->No = -1; 411 dd->Po = -1; 412 413 dd->gtol = NULL; 414 dd->ltol = NULL; 415 dd->ao = NULL; 416 PetscStrallocpy(AOBASIC,(char**)&dd->aotype); 417 dd->base = -1; 418 dd->bx = DM_BOUNDARY_NONE; 419 dd->by = DM_BOUNDARY_NONE; 420 dd->bz = DM_BOUNDARY_NONE; 421 dd->stencil_type = DMDA_STENCIL_BOX; 422 dd->interptype = DMDA_Q1; 423 dd->lx = NULL; 424 dd->ly = NULL; 425 dd->lz = NULL; 426 427 dd->elementtype = DMDA_ELEMENT_Q1; 428 429 da->ops->globaltolocalbegin = DMGlobalToLocalBegin_DA; 430 da->ops->globaltolocalend = DMGlobalToLocalEnd_DA; 431 da->ops->localtoglobalbegin = DMLocalToGlobalBegin_DA; 432 da->ops->localtoglobalend = DMLocalToGlobalEnd_DA; 433 da->ops->localtolocalbegin = DMLocalToLocalBegin_DA; 434 da->ops->localtolocalend = DMLocalToLocalEnd_DA; 435 da->ops->createglobalvector = DMCreateGlobalVector_DA; 436 da->ops->createlocalvector = DMCreateLocalVector_DA; 437 da->ops->createinterpolation = DMCreateInterpolation_DA; 438 da->ops->getcoloring = DMCreateColoring_DA; 439 da->ops->creatematrix = DMCreateMatrix_DA; 440 da->ops->refine = DMRefine_DA; 441 da->ops->coarsen = DMCoarsen_DA; 442 da->ops->refinehierarchy = DMRefineHierarchy_DA; 443 da->ops->coarsenhierarchy = DMCoarsenHierarchy_DA; 444 da->ops->getinjection = DMCreateInjection_DA; 445 da->ops->hascreateinjection = DMHasCreateInjection_DA; 446 da->ops->getaggregates = DMCreateAggregates_DA; 447 da->ops->destroy = DMDestroy_DA; 448 da->ops->view = 0; 449 da->ops->setfromoptions = DMSetFromOptions_DA; 450 da->ops->setup = DMSetUp_DA; 451 da->ops->clone = DMClone_DA; 452 da->ops->load = DMLoad_DA; 453 da->ops->createcoordinatedm = DMCreateCoordinateDM_DA; 454 da->ops->createsubdm = DMCreateSubDM_DA; 455 da->ops->createfielddecomposition = DMCreateFieldDecomposition_DA; 456 da->ops->createdomaindecomposition = DMCreateDomainDecomposition_DA; 457 da->ops->createddscatters = DMCreateDomainDecompositionScatters_DA; 458 da->ops->getdimpoints = DMGetDimPoints_DA; 459 da->ops->projectfunctionlocal = DMProjectFunctionLocal_DA; 460 da->ops->computel2diff = DMComputeL2Diff_DA; 461 da->ops->computel2gradientdiff = DMComputeL2GradientDiff_DA; 462 da->ops->getneighbors = DMGetNeighbors_DA; 463 da->ops->locatepoints = DMLocatePoints_DA_Regular; 464 da->ops->getcompatibility = DMGetCompatibility_DA; 465 ierr = PetscObjectComposeFunction((PetscObject)da,"DMSetUpGLVisViewer_C",DMSetUpGLVisViewer_DMDA);CHKERRQ(ierr); 466 PetscFunctionReturn(0); 467 } 468 469 /*@ 470 DMDACreate - Creates a DMDA object. 471 472 Collective 473 474 Input Parameter: 475 . comm - The communicator for the DMDA object 476 477 Output Parameter: 478 . da - The DMDA object 479 480 Level: advanced 481 482 Developers Note: 483 Since there exists DMDACreate1/2/3d() should this routine even exist? 484 485 .keywords: DMDA, create 486 .seealso: DMDASetSizes(), DMClone(), DMDACreate1d(), DMDACreate2d(), DMDACreate3d() 487 @*/ 488 PetscErrorCode DMDACreate(MPI_Comm comm, DM *da) 489 { 490 PetscErrorCode ierr; 491 492 PetscFunctionBegin; 493 PetscValidPointer(da,2); 494 ierr = DMCreate(comm,da);CHKERRQ(ierr); 495 ierr = DMSetType(*da,DMDA);CHKERRQ(ierr); 496 PetscFunctionReturn(0); 497 } 498 499 500