Actual source code: daindex.c
2: /*
3: Code for manipulating distributed regular arrays in parallel.
4: */
6: #include <petsc/private/dmdaimpl.h>
8: /*
9: Gets the natural number for each global number on the process.
11: Used by DMDAGetAO() and DMDAGlobalToNatural_Create()
12: */
13: PetscErrorCode DMDAGetNatural_Private(DM da, PetscInt *outNlocal, IS *isnatural)
14: {
15: PetscInt Nlocal, i, j, k, *lidx, lict = 0, dim = da->dim;
16: DM_DA *dd = (DM_DA *)da->data;
18: Nlocal = (dd->xe - dd->xs);
19: if (dim > 1) Nlocal *= (dd->ye - dd->ys);
20: if (dim > 2) Nlocal *= (dd->ze - dd->zs);
22: PetscMalloc1(Nlocal, &lidx);
24: if (dim == 1) {
25: for (i = dd->xs; i < dd->xe; i++) {
26: /* global number in natural ordering */
27: lidx[lict++] = i;
28: }
29: } else if (dim == 2) {
30: for (j = dd->ys; j < dd->ye; j++) {
31: for (i = dd->xs; i < dd->xe; i++) {
32: /* global number in natural ordering */
33: lidx[lict++] = i + j * dd->M * dd->w;
34: }
35: }
36: } else if (dim == 3) {
37: for (k = dd->zs; k < dd->ze; k++) {
38: for (j = dd->ys; j < dd->ye; j++) {
39: for (i = dd->xs; i < dd->xe; i++) lidx[lict++] = i + j * dd->M * dd->w + k * dd->M * dd->N * dd->w;
40: }
41: }
42: }
43: *outNlocal = Nlocal;
44: ISCreateGeneral(PetscObjectComm((PetscObject)da), Nlocal, lidx, PETSC_OWN_POINTER, isnatural);
45: return 0;
46: }
48: /*@C
49: DMDASetAOType - Sets the type of application ordering for a distributed array.
51: Collective on da
53: Input Parameters:
54: + da - the distributed array
55: - aotype - type of AO
57: Output Parameters:
59: Level: intermediate
61: Notes:
62: It will generate and error if an AO has already been obtained with a call to DMDAGetAO and the user sets a different AOType
64: .seealso: `DMDACreate2d()`, `DMDAGetAO()`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `DMLocalToGlobal()`
65: `DMGlobalToLocalBegin()`, `DMGlobalToLocalEnd()`, `DMLocalToLocalBegin()`, `DMLocalToLocalEnd()`, `DMDAGetGlobalIndices()`, `DMDAGetOwnershipRanges()`,
66: `AO`, `AOPetscToApplication()`, `AOApplicationToPetsc()`
67: @*/
68: PetscErrorCode DMDASetAOType(DM da, AOType aotype)
69: {
70: DM_DA *dd;
71: PetscBool isdmda;
74: PetscObjectTypeCompare((PetscObject)da, DMDA, &isdmda);
76: /* now we can safely dereference */
77: dd = (DM_DA *)da->data;
78: if (dd->ao) { /* check if the already computed AO has the same type as requested */
79: PetscBool match;
80: PetscObjectTypeCompare((PetscObject)dd->ao, aotype, &match);
82: return 0;
83: }
84: PetscFree(dd->aotype);
85: PetscStrallocpy(aotype, (char **)&dd->aotype);
86: return 0;
87: }
89: /*@
90: DMDAGetAO - Gets the application ordering context for a distributed array.
92: Collective on da
94: Input Parameter:
95: . da - the distributed array
97: Output Parameters:
98: . ao - the application ordering context for DMDAs
100: Level: intermediate
102: Notes:
103: In this case, the AO maps to the natural grid ordering that would be used
104: for the DMDA if only 1 processor were employed (ordering most rapidly in the
105: x-direction, then y, then z). Multiple degrees of freedom are numbered
106: for each node (rather than 1 component for the whole grid, then the next
107: component, etc.)
109: Do NOT call AODestroy() on the ao returned by this function.
111: .seealso: `DMDACreate2d()`, `DMDASetAOType()`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `DMLocalToGlobal()`
112: `DMGlobalToLocalBegin()`, `DMGlobalToLocalEnd()`, `DMLocalToLocalBegin()`, `DMLocalToLocalEnd()`, `DMDAGetOwnershipRanges()`,
113: `AO`, `AOPetscToApplication()`, `AOApplicationToPetsc()`
114: @*/
115: PetscErrorCode DMDAGetAO(DM da, AO *ao)
116: {
117: DM_DA *dd;
118: PetscBool isdmda;
122: PetscObjectTypeCompare((PetscObject)da, DMDA, &isdmda);
124: /* now we can safely dereference */
125: dd = (DM_DA *)da->data;
127: /*
128: Build the natural ordering to PETSc ordering mappings.
129: */
130: if (!dd->ao) {
131: IS ispetsc, isnatural;
132: PetscInt Nlocal;
134: DMDAGetNatural_Private(da, &Nlocal, &isnatural);
135: ISCreateStride(PetscObjectComm((PetscObject)da), Nlocal, dd->base, 1, &ispetsc);
136: AOCreate(PetscObjectComm((PetscObject)da), &dd->ao);
137: AOSetIS(dd->ao, isnatural, ispetsc);
138: AOSetType(dd->ao, dd->aotype);
139: ISDestroy(&ispetsc);
140: ISDestroy(&isnatural);
141: }
142: *ao = dd->ao;
143: return 0;
144: }