1 #define PETSCDM_DLL 2 3 /* 4 Code for manipulating distributed regular arrays in parallel. 5 */ 6 7 #include "private/daimpl.h" /*I "petscdm.h" I*/ 8 9 #undef __FUNCT__ 10 #define __FUNCT__ "DMDASetCoordinates" 11 /*@ 12 DMDASetCoordinates - Sets into the DMDA a vector that indicates the 13 coordinates of the local nodes (NOT including ghost nodes). 14 15 Collective on DMDA 16 17 Input Parameter: 18 + da - the distributed array 19 - c - coordinate vector 20 21 Note: 22 The coordinates should NOT include those for all ghost points 23 24 Level: intermediate 25 26 .keywords: distributed array, get, corners, nodes, local indices, coordinates 27 28 .seealso: DMDASetGhostCoordinates(), DMDAGetGhostCorners(), DMDAGetCoordinates(), DMDASetUniformCoordinates(). DMDAGetGhostedCoordinates(), DMDAGetCoordinateDA() 29 @*/ 30 PetscErrorCode DMDASetCoordinates(DM da,Vec c) 31 { 32 PetscErrorCode ierr; 33 DM_DA *dd = (DM_DA*)da->data; 34 35 PetscFunctionBegin; 36 PetscValidHeaderSpecific(da,DM_CLASSID,1); 37 PetscValidHeaderSpecific(c,VEC_CLASSID,2); 38 ierr = PetscObjectReference((PetscObject)c);CHKERRQ(ierr); 39 if (dd->coordinates) {ierr = VecDestroy(dd->coordinates);CHKERRQ(ierr);} 40 dd->coordinates = c; 41 ierr = VecSetBlockSize(c,dd->dim);CHKERRQ(ierr); 42 if (dd->ghosted_coordinates) { /* The ghosted coordinates are no longer valid */ 43 ierr = VecDestroy(dd->ghosted_coordinates);CHKERRQ(ierr); 44 dd->ghosted_coordinates = PETSC_NULL; 45 } 46 PetscFunctionReturn(0); 47 } 48 49 #undef __FUNCT__ 50 #define __FUNCT__ "DMDASetGhostedCoordinates" 51 /*@ 52 DMDASetGhostedCoordinates - Sets into the DMDA a vector that indicates the 53 coordinates of the local nodes, including ghost nodes. 54 55 Collective on DMDA 56 57 Input Parameter: 58 + da - the distributed array 59 - c - coordinate vector 60 61 Note: 62 The coordinates of interior ghost points can be set using DMDASetCoordinates() 63 followed by DMDAGetGhostedCoordinates(). This is intended to enable the setting 64 of ghost coordinates outside of the domain. 65 66 Non-ghosted coordinates, if set, are assumed still valid. 67 68 Level: intermediate 69 70 .keywords: distributed array, get, corners, nodes, local indices, coordinates 71 72 .seealso: DMDASetCoordinates(), DMDAGetGhostCorners(), DMDAGetCoordinates(), DMDASetUniformCoordinates(). DMDAGetGhostedCoordinates(), DMDAGetCoordinateDA() 73 @*/ 74 PetscErrorCode DMDASetGhostedCoordinates(DM da,Vec c) 75 { 76 PetscErrorCode ierr; 77 DM_DA *dd = (DM_DA*)da->data; 78 79 PetscFunctionBegin; 80 PetscValidHeaderSpecific(da,DM_CLASSID,1); 81 PetscValidHeaderSpecific(c,VEC_CLASSID,2); 82 ierr = PetscObjectReference((PetscObject)c);CHKERRQ(ierr); 83 if (dd->ghosted_coordinates) {ierr = VecDestroy(dd->ghosted_coordinates);CHKERRQ(ierr);} 84 dd->ghosted_coordinates = c; 85 ierr = VecSetBlockSize(c,dd->dim);CHKERRQ(ierr); 86 PetscFunctionReturn(0); 87 } 88 89 #undef __FUNCT__ 90 #define __FUNCT__ "DMDAGetCoordinates" 91 /*@ 92 DMDAGetCoordinates - Gets the node coordinates associated with a DMDA. 93 94 Not Collective 95 96 Input Parameter: 97 . da - the distributed array 98 99 Output Parameter: 100 . c - coordinate vector 101 102 Note: 103 Each process has only the coordinates for its local nodes (does NOT have the 104 coordinates for the ghost nodes). 105 106 For two and three dimensions coordinates are interlaced (x_0,y_0,x_1,y_1,...) 107 and (x_0,y_0,z_0,x_1,y_1,z_1...) 108 109 Level: intermediate 110 111 .keywords: distributed array, get, corners, nodes, local indices, coordinates 112 113 .seealso: DMDAGetGhostCorners(), DMDASetCoordinates(), DMDASetUniformCoordinates(), DMDAGetGhostedCoordinates(), DMDAGetCoordinateDA() 114 @*/ 115 PetscErrorCode DMDAGetCoordinates(DM da,Vec *c) 116 { 117 DM_DA *dd = (DM_DA*)da->data; 118 PetscFunctionBegin; 119 PetscValidHeaderSpecific(da,DM_CLASSID,1); 120 PetscValidPointer(c,2); 121 *c = dd->coordinates; 122 PetscFunctionReturn(0); 123 } 124 125 #undef __FUNCT__ 126 #define __FUNCT__ "DMDAGetCoordinateDA" 127 /*@ 128 DMDAGetCoordinateDA - Gets the DMDA that scatters between global and local DMDA coordinates 129 130 Collective on DMDA 131 132 Input Parameter: 133 . da - the distributed array 134 135 Output Parameter: 136 . dac - coordinate DMDA 137 138 Level: intermediate 139 140 .keywords: distributed array, get, corners, nodes, local indices, coordinates 141 142 .seealso: DMDAGetGhostCorners(), DMDASetCoordinates(), DMDASetUniformCoordinates(), DMDAGetCoordinates(), DMDAGetGhostedCoordinates() 143 @*/ 144 PetscErrorCode DMDAGetCoordinateDA(DM da,DM *cda) 145 { 146 PetscMPIInt size; 147 PetscErrorCode ierr; 148 DM_DA *dd = (DM_DA*)da->data; 149 150 PetscFunctionBegin; 151 if (!dd->da_coordinates) { 152 ierr = MPI_Comm_size(((PetscObject)da)->comm,&size);CHKERRQ(ierr); 153 if (dd->dim == 1) { 154 PetscInt s,m,*lc,l; 155 DMDABoundaryType pt; 156 ierr = DMDAGetInfo(da,0,&m,0,0,0,0,0,0,&s,&pt,0);CHKERRQ(ierr); 157 ierr = DMDAGetCorners(da,0,0,0,&l,0,0);CHKERRQ(ierr); 158 ierr = PetscMalloc(size*sizeof(PetscInt),&lc);CHKERRQ(ierr); 159 ierr = MPI_Allgather(&l,1,MPIU_INT,lc,1,MPIU_INT,((PetscObject)da)->comm);CHKERRQ(ierr); 160 ierr = DMDACreate1d(((PetscObject)da)->comm,pt,m,1,s,lc,&dd->da_coordinates);CHKERRQ(ierr); 161 ierr = PetscFree(lc);CHKERRQ(ierr); 162 } else if (dd->dim == 2) { 163 PetscInt i,s,m,*lc,*ld,l,k,n,M,N; 164 DMDABoundaryType pt; 165 ierr = DMDAGetInfo(da,0,&m,&n,0,&M,&N,0,0,&s,&pt,0);CHKERRQ(ierr); 166 ierr = DMDAGetCorners(da,0,0,0,&l,&k,0);CHKERRQ(ierr); 167 ierr = PetscMalloc2(size,PetscInt,&lc,size,PetscInt,&ld);CHKERRQ(ierr); 168 /* only first M values in lc matter */ 169 ierr = MPI_Allgather(&l,1,MPIU_INT,lc,1,MPIU_INT,((PetscObject)da)->comm);CHKERRQ(ierr); 170 /* every Mth value in ld matters */ 171 ierr = MPI_Allgather(&k,1,MPIU_INT,ld,1,MPIU_INT,((PetscObject)da)->comm);CHKERRQ(ierr); 172 for ( i=0; i<N; i++) { 173 ld[i] = ld[M*i]; 174 } 175 ierr = DMDACreate2d(((PetscObject)da)->comm,pt,DMDA_STENCIL_BOX,m,n,M,N,2,s,lc,ld,&dd->da_coordinates);CHKERRQ(ierr); 176 ierr = PetscFree2(lc,ld);CHKERRQ(ierr); 177 } else if (dd->dim == 3) { 178 PetscInt i,s,m,*lc,*ld,*le,l,k,q,n,M,N,P,p; 179 DMDABoundaryType pt; 180 ierr = DMDAGetInfo(da,0,&m,&n,&p,&M,&N,&P,0,&s,&pt,0);CHKERRQ(ierr); 181 ierr = DMDAGetCorners(da,0,0,0,&l,&k,&q);CHKERRQ(ierr); 182 ierr = PetscMalloc3(size,PetscInt,&lc,size,PetscInt,&ld,size,PetscInt,&le);CHKERRQ(ierr); 183 /* only first M values in lc matter */ 184 ierr = MPI_Allgather(&l,1,MPIU_INT,lc,1,MPIU_INT,((PetscObject)da)->comm);CHKERRQ(ierr); 185 /* every Mth value in ld matters */ 186 ierr = MPI_Allgather(&k,1,MPIU_INT,ld,1,MPIU_INT,((PetscObject)da)->comm);CHKERRQ(ierr); 187 for ( i=0; i<N; i++) { 188 ld[i] = ld[M*i]; 189 } 190 ierr = MPI_Allgather(&q,1,MPIU_INT,le,1,MPIU_INT,((PetscObject)da)->comm);CHKERRQ(ierr); 191 for ( i=0; i<P; i++) { 192 le[i] = le[M*N*i]; 193 } 194 ierr = DMDACreate3d(((PetscObject)da)->comm,pt,DMDA_STENCIL_BOX,m,n,p,M,N,P,3,s,lc,ld,le,&dd->da_coordinates);CHKERRQ(ierr); 195 ierr = PetscFree3(lc,ld,le);CHKERRQ(ierr); 196 } 197 } 198 *cda = dd->da_coordinates; 199 PetscFunctionReturn(0); 200 } 201 202 203 #undef __FUNCT__ 204 #define __FUNCT__ "DMDAGetGhostedCoordinates" 205 /*@ 206 DMDAGetGhostedCoordinates - Gets the node coordinates associated with a DMDA. 207 208 Collective on DMDA 209 210 Input Parameter: 211 . da - the distributed array 212 213 Output Parameter: 214 . c - coordinate vector 215 216 Note: 217 Each process has only the coordinates for its local AND ghost nodes 218 219 For two and three dimensions coordinates are interlaced (x_0,y_0,x_1,y_1,...) 220 and (x_0,y_0,z_0,x_1,y_1,z_1...) 221 222 Level: intermediate 223 224 .keywords: distributed array, get, corners, nodes, local indices, coordinates 225 226 .seealso: DMDAGetGhostCorners(), DMDASetCoordinates(), DMDASetUniformCoordinates(), DMDAGetCoordinates(), DMDAGetCoordinateDA() 227 @*/ 228 PetscErrorCode DMDAGetGhostedCoordinates(DM da,Vec *c) 229 { 230 PetscErrorCode ierr; 231 DM_DA *dd = (DM_DA*)da->data; 232 233 PetscFunctionBegin; 234 PetscValidHeaderSpecific(da,DM_CLASSID,1); 235 PetscValidPointer(c,2); 236 if (!dd->coordinates) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER,"You must call DMDASetCoordinates() before this call"); 237 if (!dd->ghosted_coordinates) { 238 DM dac; 239 ierr = DMDAGetCoordinateDA(da,&dac);CHKERRQ(ierr); 240 ierr = DMCreateLocalVector(dac,&dd->ghosted_coordinates);CHKERRQ(ierr); 241 ierr = DMGlobalToLocalBegin(dac,dd->coordinates,INSERT_VALUES,dd->ghosted_coordinates);CHKERRQ(ierr); 242 ierr = DMGlobalToLocalEnd(dac,dd->coordinates,INSERT_VALUES,dd->ghosted_coordinates);CHKERRQ(ierr); 243 } 244 *c = dd->ghosted_coordinates; 245 PetscFunctionReturn(0); 246 } 247 248 #undef __FUNCT__ 249 #define __FUNCT__ "DMDASetFieldName" 250 /*@C 251 DMDASetFieldName - Sets the names of individual field components in multicomponent 252 vectors associated with a DMDA. 253 254 Not Collective 255 256 Input Parameters: 257 + da - the distributed array 258 . nf - field number for the DMDA (0, 1, ... dof-1), where dof indicates the 259 number of degrees of freedom per node within the DMDA 260 - names - the name of the field (component) 261 262 Level: intermediate 263 264 .keywords: distributed array, get, component name 265 266 .seealso: DMDAGetFieldName() 267 @*/ 268 PetscErrorCode DMDASetFieldName(DM da,PetscInt nf,const char name[]) 269 { 270 PetscErrorCode ierr; 271 DM_DA *dd = (DM_DA*)da->data; 272 273 PetscFunctionBegin; 274 PetscValidHeaderSpecific(da,DM_CLASSID,1); 275 if (nf < 0 || nf >= dd->w) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Invalid field number: %D",nf); 276 ierr = PetscFree(dd->fieldname[nf]);CHKERRQ(ierr); 277 ierr = PetscStrallocpy(name,&dd->fieldname[nf]);CHKERRQ(ierr); 278 PetscFunctionReturn(0); 279 } 280 281 #undef __FUNCT__ 282 #define __FUNCT__ "DMDAGetFieldName" 283 /*@C 284 DMDAGetFieldName - Gets the names of individual field components in multicomponent 285 vectors associated with a DMDA. 286 287 Not Collective 288 289 Input Parameter: 290 + da - the distributed array 291 - nf - field number for the DMDA (0, 1, ... dof-1), where dof indicates the 292 number of degrees of freedom per node within the DMDA 293 294 Output Parameter: 295 . names - the name of the field (component) 296 297 Level: intermediate 298 299 .keywords: distributed array, get, component name 300 301 .seealso: DMDASetFieldName() 302 @*/ 303 PetscErrorCode DMDAGetFieldName(DM da,PetscInt nf,const char **name) 304 { 305 DM_DA *dd = (DM_DA*)da->data; 306 307 PetscFunctionBegin; 308 PetscValidHeaderSpecific(da,DM_CLASSID,1); 309 PetscValidPointer(name,3); 310 if (nf < 0 || nf >= dd->w) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Invalid field number: %D",nf); 311 *name = dd->fieldname[nf]; 312 PetscFunctionReturn(0); 313 } 314 315 #undef __FUNCT__ 316 #define __FUNCT__ "DMDAGetCorners" 317 /*@ 318 DMDAGetCorners - Returns the global (x,y,z) indices of the lower left 319 corner of the local region, excluding ghost points. 320 321 Not Collective 322 323 Input Parameter: 324 . da - the distributed array 325 326 Output Parameters: 327 + x,y,z - the corner indices (where y and z are optional; these are used 328 for 2D and 3D problems) 329 - m,n,p - widths in the corresponding directions (where n and p are optional; 330 these are used for 2D and 3D problems) 331 332 Note: 333 The corner information is independent of the number of degrees of 334 freedom per node set with the DMDACreateXX() routine. Thus the x, y, z, and 335 m, n, p can be thought of as coordinates on a logical grid, where each 336 grid point has (potentially) several degrees of freedom. 337 Any of y, z, n, and p can be passed in as PETSC_NULL if not needed. 338 339 Level: beginner 340 341 .keywords: distributed array, get, corners, nodes, local indices 342 343 .seealso: DMDAGetGhostCorners(), DMDAGetOwnershipRanges() 344 @*/ 345 PetscErrorCode DMDAGetCorners(DM da,PetscInt *x,PetscInt *y,PetscInt *z,PetscInt *m,PetscInt *n,PetscInt *p) 346 { 347 PetscInt w; 348 DM_DA *dd = (DM_DA*)da->data; 349 350 PetscFunctionBegin; 351 PetscValidHeaderSpecific(da,DM_CLASSID,1); 352 /* since the xs, xe ... have all been multiplied by the number of degrees 353 of freedom per cell, w = dd->w, we divide that out before returning.*/ 354 w = dd->w; 355 if (x) *x = dd->xs/w; if(m) *m = (dd->xe - dd->xs)/w; 356 /* the y and z have NOT been multiplied by w */ 357 if (y) *y = dd->ys; if (n) *n = (dd->ye - dd->ys); 358 if (z) *z = dd->zs; if (p) *p = (dd->ze - dd->zs); 359 PetscFunctionReturn(0); 360 } 361 362 #undef __FUNCT__ 363 #define __FUNCT__ "DMDAGetLocalBoundingBox" 364 /*@ 365 DMDAGetLocalBoundingBox - Returns the local bounding box for the DMDA. 366 367 Not Collective 368 369 Input Parameter: 370 . da - the distributed array 371 372 Output Parameters: 373 + lmin - local minimum coordinates (length dim, optional) 374 - lmax - local maximim coordinates (length dim, optional) 375 376 Level: beginner 377 378 .keywords: distributed array, get, coordinates 379 380 .seealso: DMDAGetCoordinateDA(), DMDAGetCoordinates(), DMDAGetBoundingBox() 381 @*/ 382 PetscErrorCode DMDAGetLocalBoundingBox(DM da,PetscReal lmin[],PetscReal lmax[]) 383 { 384 PetscErrorCode ierr; 385 Vec coords = PETSC_NULL; 386 PetscInt dim,i,j; 387 const PetscScalar *local_coords; 388 PetscReal min[3]={PETSC_MAX,PETSC_MAX,PETSC_MAX},max[3]={PETSC_MIN,PETSC_MIN,PETSC_MIN}; 389 PetscInt N,Ni; 390 DM_DA *dd = (DM_DA*)da->data; 391 392 PetscFunctionBegin; 393 PetscValidHeaderSpecific(da,DM_CLASSID,1); 394 dim = dd->dim; 395 ierr = DMDAGetCoordinates(da,&coords);CHKERRQ(ierr); 396 ierr = VecGetArrayRead(coords,&local_coords);CHKERRQ(ierr); 397 ierr = VecGetLocalSize(coords,&N);CHKERRQ(ierr); 398 Ni = N/dim; 399 for (i=0; i<Ni; i++) { 400 for (j=0; j<dim; j++) { 401 min[j] = PetscMin(min[j],PetscRealPart(local_coords[i*dim+j]));CHKERRQ(ierr); 402 max[j] = PetscMax(min[j],PetscRealPart(local_coords[i*dim+j]));CHKERRQ(ierr); 403 } 404 } 405 ierr = VecRestoreArrayRead(coords,&local_coords);CHKERRQ(ierr); 406 if (lmin) {ierr = PetscMemcpy(lmin,min,dim*sizeof(PetscReal));CHKERRQ(ierr);} 407 if (lmax) {ierr = PetscMemcpy(lmax,max,dim*sizeof(PetscReal));CHKERRQ(ierr);} 408 PetscFunctionReturn(0); 409 } 410 411 #undef __FUNCT__ 412 #define __FUNCT__ "DMDAGetBoundingBox" 413 /*@ 414 DMDAGetBoundingBox - Returns the global bounding box for the DMDA. 415 416 Collective on DMDA 417 418 Input Parameter: 419 . da - the distributed array 420 421 Output Parameters: 422 + gmin - global minimum coordinates (length dim, optional) 423 - gmax - global maximim coordinates (length dim, optional) 424 425 Level: beginner 426 427 .keywords: distributed array, get, coordinates 428 429 .seealso: DMDAGetCoordinateDA(), DMDAGetCoordinates(), DMDAGetLocalBoundingBox() 430 @*/ 431 PetscErrorCode DMDAGetBoundingBox(DM da,PetscReal gmin[],PetscReal gmax[]) 432 { 433 PetscErrorCode ierr; 434 PetscMPIInt count; 435 PetscReal lmin[3],lmax[3]; 436 DM_DA *dd = (DM_DA*)da->data; 437 438 PetscFunctionBegin; 439 PetscValidHeaderSpecific(da,DM_CLASSID,1); 440 count = PetscMPIIntCast(dd->dim); 441 ierr = DMDAGetLocalBoundingBox(da,lmin,lmax);CHKERRQ(ierr); 442 if (gmin) {ierr = MPI_Allreduce(lmin,gmin,count,MPIU_REAL,MPI_MIN,((PetscObject)da)->comm);CHKERRQ(ierr);} 443 if (gmax) {ierr = MPI_Allreduce(lmax,gmax,count,MPIU_REAL,MPI_MAX,((PetscObject)da)->comm);CHKERRQ(ierr);} 444 PetscFunctionReturn(0); 445 } 446