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