Actual source code: daindex.c
1: /*
2: Code for manipulating distributed regular arrays in parallel.
3: */
5: #include src/dm/da/daimpl.h
9: /*@C
10: DAGetGlobalIndices - Returns the global node number of all local nodes,
11: including ghost nodes.
13: Not Collective
15: Input Parameter:
16: . da - the distributed array
18: Output Parameters:
19: + n - the number of local elements, including ghost nodes (or PETSC_NULL)
20: - idx - the global indices
22: Level: intermediate
24: Note:
25: For DA_STENCIL_STAR stencils the inactive corner ghost nodes are also included
26: in the list of local indices (even though those nodes are not updated
27: during calls to DAXXXToXXX().
29: Essentially the same data is returned in the form of a local-to-global mapping
30: with the routine DAGetISLocalToGlobalMapping();
32: Fortran Note:
33: This routine is used differently from Fortran
34: .vb
35: DA da
36: integer n,da_array(1)
37: PetscOffset i_da
38: integer ierr
39: call DAGetGlobalIndices(da,n,da_array,i_da,ierr)
41: C Access first local entry in list
42: value = da_array(i_da + 1)
43: .ve
45: See the Fortran chapter of the users manual for details
47: .keywords: distributed array, get, global, indices, local-to-global
49: .seealso: DACreate2d(), DAGetGhostCorners(), DAGetCorners(), DALocalToGlobal()
50: DAGlobalToLocalBegin(), DAGlobalToLocalEnd(), DALocalToLocalBegin(), DAGetAO(), DAGetGlobalIndicesF90()
51: DAGetISLocalToGlobalMapping(), DACreate3d(), DACreate1d(), DALocalToLocalEnd()
52: @*/
53: PetscErrorCode DAGetGlobalIndices(DA da,PetscInt *n,PetscInt **idx)
54: {
57: if (n) *n = da->Nl;
58: if (idx) *idx = da->idx;
59: return(0);
60: }
64: /*
65: Gets the natural number for each global number on the process.
67: Used by DAGetAO() and DAGlobalToNatural_Create()
68: */
69: PetscErrorCode DAGetNatural_Private(DA da,PetscInt *outNlocal,IS *isnatural)
70: {
72: PetscInt Nlocal,i,j,k,*lidx,lict = 0;
75: Nlocal = (da->xe-da->xs);
76: if (da->dim > 1) {
77: Nlocal *= (da->ye-da->ys);
78: }
79: if (da->dim > 2) {
80: Nlocal *= (da->ze-da->zs);
81: }
82:
83: PetscMalloc(Nlocal*sizeof(PetscInt),&lidx);
84:
85: if (da->dim == 1) {
86: for (i=da->xs; i<da->xe; i++) {
87: /* global number in natural ordering */
88: lidx[lict++] = i;
89: }
90: } else if (da->dim == 2) {
91: for (j=da->ys; j<da->ye; j++) {
92: for (i=da->xs; i<da->xe; i++) {
93: /* global number in natural ordering */
94: lidx[lict++] = i + j*da->M*da->w;
95: }
96: }
97: } else if (da->dim == 3) {
98: for (k=da->zs; k<da->ze; k++) {
99: for (j=da->ys; j<da->ye; j++) {
100: for (i=da->xs; i<da->xe; i++) {
101: lidx[lict++] = i + j*da->M*da->w + k*da->M*da->N*da->w;
102: }
103: }
104: }
105: }
106: *outNlocal = Nlocal;
107: ISCreateGeneral(da->comm,Nlocal,lidx,isnatural);
108: PetscFree(lidx);
109: return(0);
110: }
114: /*@C
115: DAGetAO - Gets the application ordering context for a distributed array.
117: Collective on DA
119: Input Parameter:
120: . da - the distributed array
122: Output Parameters:
123: . ao - the application ordering context for DAs
125: Level: intermediate
127: Notes:
128: In this case, the AO maps to the natural grid ordering that would be used
129: for the DA if only 1 processor were employed (ordering most rapidly in the
130: x-direction, then y, then z). Multiple degrees of freedom are numbered
131: for each node (rather than 1 component for the whole grid, then the next
132: component, etc.)
134: .keywords: distributed array, get, global, indices, local-to-global
136: .seealso: DACreate2d(), DAGetGhostCorners(), DAGetCorners(), DALocalToGlocal()
137: DAGlobalToLocalBegin(), DAGlobalToLocalEnd(), DALocalToLocalBegin(), DALocalToLocalEnd(), DAGetGlobalIndices()
138: @*/
139: PetscErrorCode DAGetAO(DA da,AO *ao)
140: {
145: /*
146: Build the natural ordering to PETSc ordering mappings.
147: */
148: if (!da->ao) {
149: IS ispetsc,isnatural;
151: PetscInt Nlocal;
153: DAGetNatural_Private(da,&Nlocal,&isnatural);
154: ISCreateStride(da->comm,Nlocal,da->base,1,&ispetsc);
155: AOCreateBasicIS(isnatural,ispetsc,&da->ao);
156: PetscLogObjectParent(da,da->ao);
157: ISDestroy(ispetsc);
158: ISDestroy(isnatural);
159: }
160: *ao = da->ao;
161: return(0);
162: }
164: /*MC
165: DAGetGlobalIndicesF90 - Returns a Fortran90 pointer to the list of
166: global indices (global node number of all local nodes, including
167: ghost nodes).
169: Synopsis:
170: DAGetGlobalIndicesF90(DA da,integer n,{integer, pointer :: idx(:)},integer ierr)
172: Input Parameter:
173: . da - the distributed array
175: Output Parameters:
176: + n - the number of local elements, including ghost nodes (or PETSC_NULL)
177: . idx - the Fortran90 pointer to the global indices
178: - ierr - error code
180: Level: intermediate
182: Notes:
183: Not yet supported for all F90 compilers
185: .keywords: distributed array, get, global, indices, local-to-global, f90
187: .seealso: DAGetGlobalIndices()
188: M*/