Actual source code: daview.c
1: /*
2: Code for manipulating distributed regular arrays in parallel.
3: */
5: #include <petsc/private/dmdaimpl.h>
7: #if defined(PETSC_HAVE_MATLAB)
8: #include <mat.h> /* MATLAB include file */
10: PetscErrorCode DMView_DA_Matlab(DM da, PetscViewer viewer)
11: {
12: PetscMPIInt rank;
13: PetscInt dim, m, n, p, dof, swidth;
14: DMDAStencilType stencil;
15: DMBoundaryType bx, by, bz;
16: mxArray *mx;
17: const char *fnames[] = {"dimension", "m", "n", "p", "dof", "stencil_width", "bx", "by", "bz", "stencil_type"};
19: PetscFunctionBegin;
20: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)da), &rank));
21: if (rank == 0) {
22: PetscCall(DMDAGetInfo(da, &dim, &m, &n, &p, 0, 0, 0, &dof, &swidth, &bx, &by, &bz, &stencil));
23: mx = mxCreateStructMatrix(1, 1, 8, (const char **)fnames);
24: PetscCheck(mx, PETSC_COMM_SELF, PETSC_ERR_LIB, "Unable to generate MATLAB struct array to hold DMDA information");
25: mxSetFieldByNumber(mx, 0, 0, mxCreateDoubleScalar((double)dim));
26: mxSetFieldByNumber(mx, 0, 1, mxCreateDoubleScalar((double)m));
27: mxSetFieldByNumber(mx, 0, 2, mxCreateDoubleScalar((double)n));
28: mxSetFieldByNumber(mx, 0, 3, mxCreateDoubleScalar((double)p));
29: mxSetFieldByNumber(mx, 0, 4, mxCreateDoubleScalar((double)dof));
30: mxSetFieldByNumber(mx, 0, 5, mxCreateDoubleScalar((double)swidth));
31: mxSetFieldByNumber(mx, 0, 6, mxCreateDoubleScalar((double)bx));
32: mxSetFieldByNumber(mx, 0, 7, mxCreateDoubleScalar((double)by));
33: mxSetFieldByNumber(mx, 0, 8, mxCreateDoubleScalar((double)bz));
34: mxSetFieldByNumber(mx, 0, 9, mxCreateDoubleScalar((double)stencil));
35: PetscCall(PetscObjectName((PetscObject)da));
36: PetscCall(PetscViewerMatlabPutVariable(viewer, ((PetscObject)da)->name, mx));
37: }
38: PetscFunctionReturn(PETSC_SUCCESS);
39: }
40: #endif
42: PetscErrorCode DMView_DA_Binary(DM da, PetscViewer viewer)
43: {
44: PetscMPIInt rank;
45: PetscInt dim, m, n, p, dof, swidth, M, N, P;
46: DMDAStencilType stencil;
47: DMBoundaryType bx, by, bz;
48: MPI_Comm comm;
49: PetscBool coors = PETSC_FALSE;
50: Vec coordinates;
52: PetscFunctionBegin;
53: PetscCall(PetscObjectGetComm((PetscObject)da, &comm));
55: PetscCall(DMDAGetInfo(da, &dim, &m, &n, &p, &M, &N, &P, &dof, &swidth, &bx, &by, &bz, &stencil));
56: PetscCallMPI(MPI_Comm_rank(comm, &rank));
57: PetscCall(DMGetCoordinates(da, &coordinates));
58: if (rank == 0) {
59: PetscCall(PetscViewerBinaryWrite(viewer, &dim, 1, PETSC_INT));
60: PetscCall(PetscViewerBinaryWrite(viewer, &m, 1, PETSC_INT));
61: PetscCall(PetscViewerBinaryWrite(viewer, &n, 1, PETSC_INT));
62: PetscCall(PetscViewerBinaryWrite(viewer, &p, 1, PETSC_INT));
63: PetscCall(PetscViewerBinaryWrite(viewer, &dof, 1, PETSC_INT));
64: PetscCall(PetscViewerBinaryWrite(viewer, &swidth, 1, PETSC_INT));
65: PetscCall(PetscViewerBinaryWrite(viewer, &bx, 1, PETSC_ENUM));
66: PetscCall(PetscViewerBinaryWrite(viewer, &by, 1, PETSC_ENUM));
67: PetscCall(PetscViewerBinaryWrite(viewer, &bz, 1, PETSC_ENUM));
68: PetscCall(PetscViewerBinaryWrite(viewer, &stencil, 1, PETSC_ENUM));
69: if (coordinates) coors = PETSC_TRUE;
70: PetscCall(PetscViewerBinaryWrite(viewer, &coors, 1, PETSC_BOOL));
71: }
73: /* save the coordinates if they exist to disk (in the natural ordering) */
74: if (coordinates) PetscCall(VecView(coordinates, viewer));
75: PetscFunctionReturn(PETSC_SUCCESS);
76: }
78: /*@
79: DMDAGetInfo - Gets information about a given distributed array.
81: Not Collective
83: Input Parameter:
84: . da - the `DMDA`
86: Output Parameters:
87: + dim - dimension of the `DMDA` (1, 2, or 3)
88: . M - global dimension in first direction of the array
89: . N - global dimension in second direction of the array
90: . P - global dimension in third direction of the array
91: . m - corresponding number of MPI processes in first dimension
92: . n - corresponding number of MPI processes in second dimension
93: . p - corresponding number of MPI processes in third dimension
94: . dof - number of degrees of freedom per node
95: . s - stencil width
96: . bx - type of ghost nodes at boundary in first dimension
97: . by - type of ghost nodes at boundary in second dimension
98: . bz - type of ghost nodes at boundary in third dimension
99: - st - stencil type, either `DMDA_STENCIL_STAR` or `DMDA_STENCIL_BOX`
101: Level: beginner
103: Note:
104: Use `NULL` (`PETSC_NULL_INTEGER` in Fortran) in place of any output parameter that is not of interest.
106: .seealso: [](sec_struct), `DM`, `DMDA`, `DMView()`, `DMDAGetCorners()`, `DMDAGetLocalInfo()`
107: @*/
108: PetscErrorCode DMDAGetInfo(DM da, PeOp PetscInt *dim, PeOp PetscInt *M, PeOp PetscInt *N, PeOp PetscInt *P, PeOp PetscInt *m, PeOp PetscInt *n, PeOp PetscInt *p, PeOp PetscInt *dof, PeOp PetscInt *s, PeOp DMBoundaryType *bx, PeOp DMBoundaryType *by, PeOp DMBoundaryType *bz, PeOp DMDAStencilType *st)
109: {
110: DM_DA *dd = (DM_DA *)da->data;
112: PetscFunctionBegin;
114: if (dim) *dim = da->dim;
115: if (M) {
116: if (dd->Mo < 0) *M = dd->M;
117: else *M = dd->Mo;
118: }
119: if (N) {
120: if (dd->No < 0) *N = dd->N;
121: else *N = dd->No;
122: }
123: if (P) {
124: if (dd->Po < 0) *P = dd->P;
125: else *P = dd->Po;
126: }
127: if (m) *m = dd->m;
128: if (n) *n = dd->n;
129: if (p) *p = dd->p;
130: if (dof) *dof = dd->w;
131: if (s) *s = dd->s;
132: if (bx) *bx = dd->bx;
133: if (by) *by = dd->by;
134: if (bz) *bz = dd->bz;
135: if (st) *st = dd->stencil_type;
136: PetscFunctionReturn(PETSC_SUCCESS);
137: }
139: /*@C
140: DMDAGetLocalInfo - Gets information about a given `DMDA` and this MPI process's location in it
142: Not Collective
144: Input Parameter:
145: . da - the `DMDA`
147: Output Parameter:
148: . info - structure containing the information
150: Level: beginner
152: Note:
153: See `DMDALocalInfo` for the information that is returned
155: .seealso: [](sec_struct), `DM`, `DMDA`, `DMDAGetInfo()`, `DMDAGetCorners()`, `DMDALocalInfo`
156: @*/
157: PetscErrorCode DMDAGetLocalInfo(DM da, DMDALocalInfo *info)
158: {
159: PetscInt w;
160: DM_DA *dd = (DM_DA *)da->data;
162: PetscFunctionBegin;
164: PetscAssertPointer(info, 2);
165: info->da = da;
166: info->dim = da->dim;
167: if (dd->Mo < 0) info->mx = dd->M;
168: else info->mx = dd->Mo;
169: if (dd->No < 0) info->my = dd->N;
170: else info->my = dd->No;
171: if (dd->Po < 0) info->mz = dd->P;
172: else info->mz = dd->Po;
173: info->dof = dd->w;
174: info->sw = dd->s;
175: info->bx = dd->bx;
176: info->by = dd->by;
177: info->bz = dd->bz;
178: info->st = dd->stencil_type;
180: /* since the xs, xe ... have all been multiplied by the number of degrees
181: of freedom per cell, w = dd->w, we divide that out before returning.*/
182: w = dd->w;
183: info->xs = dd->xs / w + dd->xo;
184: info->xm = (dd->xe - dd->xs) / w;
185: /* the y and z have NOT been multiplied by w */
186: info->ys = dd->ys + dd->yo;
187: info->ym = (dd->ye - dd->ys);
188: info->zs = dd->zs + dd->zo;
189: info->zm = (dd->ze - dd->zs);
191: info->gxs = dd->Xs / w + dd->xo;
192: info->gxm = (dd->Xe - dd->Xs) / w;
193: /* the y and z have NOT been multiplied by w */
194: info->gys = dd->Ys + dd->yo;
195: info->gym = (dd->Ye - dd->Ys);
196: info->gzs = dd->Zs + dd->zo;
197: info->gzm = (dd->Ze - dd->Zs);
198: PetscFunctionReturn(PETSC_SUCCESS);
199: }