1*47c6ae99SBarry Smith #define PETSCDM_DLL 2*47c6ae99SBarry Smith 3*47c6ae99SBarry Smith /* 4*47c6ae99SBarry Smith Code for manipulating distributed regular arrays in parallel. 5*47c6ae99SBarry Smith */ 6*47c6ae99SBarry Smith 7*47c6ae99SBarry Smith #include "private/daimpl.h" /*I "petscda.h" I*/ 8*47c6ae99SBarry Smith 9*47c6ae99SBarry Smith #undef __FUNCT__ 10*47c6ae99SBarry Smith #define __FUNCT__ "DAGetGhostCorners" 11*47c6ae99SBarry Smith /*@ 12*47c6ae99SBarry Smith DAGetGhostCorners - Returns the global (x,y,z) indices of the lower left 13*47c6ae99SBarry Smith corner of the local region, including ghost points. 14*47c6ae99SBarry Smith 15*47c6ae99SBarry Smith Not Collective 16*47c6ae99SBarry Smith 17*47c6ae99SBarry Smith Input Parameter: 18*47c6ae99SBarry Smith . da - the distributed array 19*47c6ae99SBarry Smith 20*47c6ae99SBarry Smith Output Parameters: 21*47c6ae99SBarry Smith + x,y,z - the corner indices (where y and z are optional; these are used 22*47c6ae99SBarry Smith for 2D and 3D problems) 23*47c6ae99SBarry Smith - m,n,p - widths in the corresponding directions (where n and p are optional; 24*47c6ae99SBarry Smith these are used for 2D and 3D problems) 25*47c6ae99SBarry Smith 26*47c6ae99SBarry Smith Level: beginner 27*47c6ae99SBarry Smith 28*47c6ae99SBarry Smith Note: 29*47c6ae99SBarry Smith The corner information is independent of the number of degrees of 30*47c6ae99SBarry Smith freedom per node set with the DACreateXX() routine. Thus the x, y, z, and 31*47c6ae99SBarry Smith m, n, p can be thought of as coordinates on a logical grid, where each 32*47c6ae99SBarry Smith grid point has (potentially) several degrees of freedom. 33*47c6ae99SBarry Smith Any of y, z, n, and p can be passed in as PETSC_NULL if not needed. 34*47c6ae99SBarry Smith 35*47c6ae99SBarry Smith .keywords: distributed array, get, ghost, corners, nodes, local indices 36*47c6ae99SBarry Smith 37*47c6ae99SBarry Smith .seealso: DAGetCorners(), DACreate1d(), DACreate2d(), DACreate3d(), DAGetOwnershipRanges() 38*47c6ae99SBarry Smith 39*47c6ae99SBarry Smith @*/ 40*47c6ae99SBarry Smith PetscErrorCode PETSCDM_DLLEXPORT DAGetGhostCorners(DA da,PetscInt *x,PetscInt *y,PetscInt *z,PetscInt *m,PetscInt *n,PetscInt *p) 41*47c6ae99SBarry Smith { 42*47c6ae99SBarry Smith PetscInt w; 43*47c6ae99SBarry Smith DM_DA *dd = (DM_DA*)da->data; 44*47c6ae99SBarry Smith 45*47c6ae99SBarry Smith PetscFunctionBegin; 46*47c6ae99SBarry Smith PetscValidHeaderSpecific(da,DM_CLASSID,1); 47*47c6ae99SBarry Smith /* since the xs, xe ... have all been multiplied by the number of degrees 48*47c6ae99SBarry Smith of freedom per cell, w = dd->w, we divide that out before returning.*/ 49*47c6ae99SBarry Smith w = dd->w; 50*47c6ae99SBarry Smith if (x) *x = dd->Xs/w; if (m) *m = (dd->Xe - dd->Xs)/w; 51*47c6ae99SBarry Smith if (y) *y = dd->Ys; if (n) *n = (dd->Ye - dd->Ys); 52*47c6ae99SBarry Smith if (z) *z = dd->Zs; if (p) *p = (dd->Ze - dd->Zs); 53*47c6ae99SBarry Smith PetscFunctionReturn(0); 54*47c6ae99SBarry Smith } 55*47c6ae99SBarry Smith 56