Actual source code: grglvis.c
1: /* Routines to visualize DMDAs and fields through GLVis */
3: #include <petsc/private/dmdaimpl.h>
4: #include <petsc/private/glvisviewerimpl.h>
6: typedef struct {
7: PetscBool ll;
8: } DMDAGhostedGLVisViewerCtx;
10: static PetscErrorCode DMDAGhostedDestroyGLVisViewerCtx_Private(void **vctx)
11: {
12: PetscFunctionBegin;
13: PetscCall(PetscFree(*vctx));
14: PetscFunctionReturn(PETSC_SUCCESS);
15: }
17: typedef struct {
18: Vec xlocal;
19: } DMDAFieldGLVisViewerCtx;
21: static PetscErrorCode DMDAFieldDestroyGLVisViewerCtx_Private(void *vctx)
22: {
23: DMDAFieldGLVisViewerCtx *ctx = (DMDAFieldGLVisViewerCtx *)vctx;
25: PetscFunctionBegin;
26: PetscCall(VecDestroy(&ctx->xlocal));
27: PetscCall(PetscFree(vctx));
28: PetscFunctionReturn(PETSC_SUCCESS);
29: }
31: /*
32: dactx->ll is false -> all but the last proc per dimension claim the ghosted node on the right
33: dactx->ll is true -> all but the first proc per dimension claim the ghosted node on the left
34: */
35: static PetscErrorCode DMDAGetNumElementsGhosted(DM da, PetscInt *nex, PetscInt *ney, PetscInt *nez)
36: {
37: DMDAGhostedGLVisViewerCtx *dactx;
38: PetscInt sx, sy, sz, ien, jen, ken;
40: PetscFunctionBegin;
41: /* Appease -Wmaybe-uninitialized */
42: if (nex) *nex = -1;
43: if (ney) *ney = -1;
44: if (nez) *nez = -1;
45: PetscCall(DMDAGetCorners(da, &sx, &sy, &sz, &ien, &jen, &ken));
46: PetscCall(DMGetApplicationContext(da, &dactx));
47: if (dactx->ll) {
48: PetscInt dim;
50: PetscCall(DMGetDimension(da, &dim));
51: if (!sx) ien--;
52: if (!sy && dim > 1) jen--;
53: if (!sz && dim > 2) ken--;
54: } else {
55: PetscInt M, N, P;
57: PetscCall(DMDAGetInfo(da, NULL, &M, &N, &P, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL));
58: if (sx + ien == M) ien--;
59: if (sy + jen == N) jen--;
60: if (sz + ken == P) ken--;
61: }
62: if (nex) *nex = ien;
63: if (ney) *ney = jen;
64: if (nez) *nez = ken;
65: PetscFunctionReturn(PETSC_SUCCESS);
66: }
68: /* inherits number of vertices from DMDAGetNumElementsGhosted */
69: static PetscErrorCode DMDAGetNumVerticesGhosted(DM da, PetscInt *nvx, PetscInt *nvy, PetscInt *nvz)
70: {
71: PetscInt ien = 0, jen = 0, ken = 0, dim;
72: PetscInt tote;
74: PetscFunctionBegin;
75: PetscCall(DMGetDimension(da, &dim));
76: PetscCall(DMDAGetNumElementsGhosted(da, &ien, &jen, &ken));
77: tote = ien * (dim > 1 ? jen : 1) * (dim > 2 ? ken : 1);
78: if (tote) {
79: ien = ien + 1;
80: jen = dim > 1 ? jen + 1 : 1;
81: ken = dim > 2 ? ken + 1 : 1;
82: }
83: if (nvx) *nvx = ien;
84: if (nvy) *nvy = jen;
85: if (nvz) *nvz = ken;
86: PetscFunctionReturn(PETSC_SUCCESS);
87: }
89: static PetscErrorCode DMDASampleGLVisFields_Private(PetscObject oX, PetscInt nf, PetscObject oXf[], void *vctx)
90: {
91: DM da;
92: DMDAFieldGLVisViewerCtx *ctx = (DMDAFieldGLVisViewerCtx *)vctx;
93: DMDAGhostedGLVisViewerCtx *dactx;
94: const PetscScalar *array;
95: PetscScalar **arrayf;
96: PetscInt i, f, ii, ien, jen, ken, ie, je, ke, bs, *bss, *nfs;
97: PetscInt sx, sy, sz, gsx, gsy, gsz, ist, jst, kst, gm, gn, gp;
99: PetscFunctionBegin;
100: PetscCall(VecGetDM(ctx->xlocal, &da));
101: PetscCheck(da, PetscObjectComm(oX), PETSC_ERR_ARG_WRONG, "Vector not generated from a DMDA");
102: PetscCall(DMGetApplicationContext(da, &dactx));
103: PetscCall(VecGetBlockSize(ctx->xlocal, &bs));
104: PetscCall(DMGlobalToLocalBegin(da, (Vec)oX, INSERT_VALUES, ctx->xlocal));
105: PetscCall(DMGlobalToLocalEnd(da, (Vec)oX, INSERT_VALUES, ctx->xlocal));
106: PetscCall(DMDAGetNumVerticesGhosted(da, &ien, &jen, &ken));
107: PetscCall(DMDAGetGhostCorners(da, &gsx, &gsy, &gsz, &gm, &gn, &gp));
108: PetscCall(DMDAGetCorners(da, &sx, &sy, &sz, NULL, NULL, NULL));
109: if (dactx->ll) {
110: kst = jst = ist = 0;
111: } else {
112: kst = gsz != sz ? 1 : 0;
113: jst = gsy != sy ? 1 : 0;
114: ist = gsx != sx ? 1 : 0;
115: }
116: PetscCall(PetscMalloc3(nf, &arrayf, nf, &bss, nf, &nfs));
117: PetscCall(VecGetArrayRead(ctx->xlocal, &array));
118: for (f = 0; f < nf; f++) {
119: PetscCall(VecGetBlockSize((Vec)oXf[f], &bss[f]));
120: PetscCall(VecGetArray((Vec)oXf[f], &arrayf[f]));
121: PetscCall(VecGetLocalSize((Vec)oXf[f], &nfs[f]));
122: }
123: for (ke = kst, ii = 0; ke < kst + ken; ke++) {
124: for (je = jst; je < jst + jen; je++) {
125: for (ie = ist; ie < ist + ien; ie++) {
126: PetscInt cf, b;
127: i = ke * gm * gn + je * gm + ie;
128: for (f = 0, cf = 0; f < nf; f++) {
129: if (!nfs[f]) continue;
130: for (b = 0; b < bss[f]; b++) arrayf[f][bss[f] * ii + b] = array[i * bs + cf++];
131: }
132: ii++;
133: }
134: }
135: }
136: for (f = 0; f < nf; f++) PetscCall(VecRestoreArray((Vec)oXf[f], &arrayf[f]));
137: PetscCall(VecRestoreArrayRead(ctx->xlocal, &array));
138: PetscCall(PetscFree3(arrayf, bss, nfs));
139: PetscFunctionReturn(PETSC_SUCCESS);
140: }
142: PetscErrorCode DMSetUpGLVisViewer_DMDA(PetscObject oda, PetscViewer viewer)
143: {
144: DM da = (DM)oda, daview;
146: PetscFunctionBegin;
147: PetscCall(PetscObjectQuery(oda, "GLVisGraphicsDMDAGhosted", (PetscObject *)&daview));
148: if (!daview) {
149: DMDAGhostedGLVisViewerCtx *dactx;
150: DM dacoord = NULL;
151: Vec xcoor, xcoorl;
152: PetscBool hashocoord = PETSC_FALSE;
153: const PetscInt *lx, *ly, *lz;
154: PetscInt dim, M, N, P, m, n, p, dof, s, i;
156: PetscCall(PetscNew(&dactx));
157: dactx->ll = PETSC_TRUE; /* default to match elements layout obtained by DMDAGetElements */
158: PetscOptionsBegin(PetscObjectComm(oda), oda->prefix, "GLVis PetscViewer DMDA Options", "PetscViewer");
159: PetscCall(PetscOptionsBool("-viewer_glvis_dm_da_ll", "Left-looking subdomain view", NULL, dactx->ll, &dactx->ll, NULL));
160: PetscOptionsEnd();
161: /* Create a properly ghosted DMDA to visualize the mesh and the fields associated with */
162: PetscCall(DMDAGetInfo(da, &dim, &M, &N, &P, &m, &n, &p, &dof, &s, NULL, NULL, NULL, NULL));
163: PetscCall(DMDAGetOwnershipRanges(da, &lx, &ly, &lz));
164: PetscCall(PetscObjectQuery((PetscObject)da, "_glvis_mesh_coords", (PetscObject *)&xcoor));
165: if (!xcoor) {
166: PetscCall(DMGetCoordinates(da, &xcoor));
167: } else {
168: hashocoord = PETSC_TRUE;
169: }
170: PetscCall(PetscInfo(da, "Creating auxiliary DMDA for managing GLVis graphics\n"));
171: switch (dim) {
172: case 1:
173: PetscCall(DMDACreate1d(PetscObjectComm((PetscObject)da), DM_BOUNDARY_NONE, M, dof, 1, lx, &daview));
174: if (!hashocoord) PetscCall(DMDACreate1d(PetscObjectComm((PetscObject)da), DM_BOUNDARY_NONE, M, 1, 1, lx, &dacoord));
175: break;
176: case 2:
177: PetscCall(DMDACreate2d(PetscObjectComm((PetscObject)da), DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, M, N, m, n, dof, 1, lx, ly, &daview));
178: if (!hashocoord) PetscCall(DMDACreate2d(PetscObjectComm((PetscObject)da), DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, M, N, m, n, 2, 1, lx, ly, &dacoord));
179: break;
180: case 3:
181: PetscCall(DMDACreate3d(PetscObjectComm((PetscObject)da), DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, M, N, P, m, n, p, dof, 1, lx, ly, lz, &daview));
182: if (!hashocoord) PetscCall(DMDACreate3d(PetscObjectComm((PetscObject)da), DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, M, N, P, m, n, p, 3, 1, lx, ly, lz, &dacoord));
183: break;
184: default:
185: SETERRQ(PetscObjectComm((PetscObject)da), PETSC_ERR_SUP, "Unsupported dimension %" PetscInt_FMT, dim);
186: }
187: PetscCall(DMSetApplicationContext(daview, dactx));
188: PetscCall(DMSetApplicationContextDestroy(daview, DMDAGhostedDestroyGLVisViewerCtx_Private));
189: PetscCall(DMSetUp(daview));
190: if (!xcoor) {
191: PetscCall(DMDASetUniformCoordinates(daview, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0));
192: PetscCall(DMGetCoordinates(daview, &xcoor));
193: }
194: if (dacoord) {
195: PetscCall(DMSetUp(dacoord));
196: PetscCall(DMCreateLocalVector(dacoord, &xcoorl));
197: PetscCall(DMGlobalToLocalBegin(dacoord, xcoor, INSERT_VALUES, xcoorl));
198: PetscCall(DMGlobalToLocalEnd(dacoord, xcoor, INSERT_VALUES, xcoorl));
199: PetscCall(DMDestroy(&dacoord));
200: } else {
201: PetscInt ien, jen, ken, nc, nl, cdof, deg;
202: char fecmesh[64];
203: const char *name;
204: PetscBool flg;
206: PetscCall(DMDAGetNumElementsGhosted(daview, &ien, &jen, &ken));
207: nc = ien * (jen > 0 ? jen : 1) * (ken > 0 ? ken : 1);
209: PetscCall(VecGetLocalSize(xcoor, &nl));
210: PetscCheck(!nc || (nl % nc) == 0, PETSC_COMM_SELF, PETSC_ERR_SUP, "Incompatible local coordinate size %" PetscInt_FMT " and number of cells %" PetscInt_FMT, nl, nc);
211: PetscCall(VecDuplicate(xcoor, &xcoorl));
212: PetscCall(VecCopy(xcoor, xcoorl));
213: PetscCall(VecSetDM(xcoorl, NULL));
214: PetscCall(PetscObjectGetName((PetscObject)xcoor, &name));
215: PetscCall(PetscStrbeginswith(name, "FiniteElementCollection:", &flg));
216: if (!flg) {
217: deg = 0;
218: if (nc && nl) {
219: cdof = nl / (nc * dim);
220: deg = 1;
221: while (1) {
222: PetscInt degd = 1;
223: for (i = 0; i < dim; i++) degd *= (deg + 1);
224: PetscCheck(degd <= cdof, PETSC_COMM_SELF, PETSC_ERR_SUP, "Cell dofs %" PetscInt_FMT, cdof);
225: if (degd == cdof) break;
226: deg++;
227: }
228: }
229: PetscCall(PetscSNPrintf(fecmesh, sizeof(fecmesh), "FiniteElementCollection: L2_T1_%" PetscInt_FMT "D_P%" PetscInt_FMT, dim, deg));
230: PetscCall(PetscObjectSetName((PetscObject)xcoorl, fecmesh));
231: } else PetscCall(PetscObjectSetName((PetscObject)xcoorl, name));
232: }
234: /* xcoorl is composed with the ghosted DMDA, the ghosted coordinate DMDA (if present) is only available through this vector */
235: PetscCall(PetscObjectCompose((PetscObject)daview, "GLVisGraphicsCoordsGhosted", (PetscObject)xcoorl));
236: PetscCall(PetscObjectDereference((PetscObject)xcoorl));
238: /* daview is composed with the original DMDA */
239: PetscCall(PetscObjectCompose(oda, "GLVisGraphicsDMDAGhosted", (PetscObject)daview));
240: PetscCall(PetscObjectDereference((PetscObject)daview));
241: }
243: /* customize the viewer if present */
244: if (viewer) {
245: DMDAFieldGLVisViewerCtx *ctx;
246: DMDAGhostedGLVisViewerCtx *dactx;
247: char fec[64];
248: Vec xlocal, *Ufield;
249: const char **dafieldname;
250: char **fec_type, **fieldname;
251: PetscInt *nlocal, *bss, *dims;
252: PetscInt dim, M, N, P, dof, s, i, nf;
253: PetscBool bsset;
255: PetscCall(DMDAGetInfo(daview, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));
256: PetscCall(DMGetApplicationContext(daview, &dactx));
257: PetscCall(DMCreateLocalVector(daview, &xlocal));
258: PetscCall(DMDAGetFieldNames(da, (const char *const **)&dafieldname));
259: PetscCall(DMDAGetNumVerticesGhosted(daview, &M, &N, &P));
260: PetscCall(PetscSNPrintf(fec, sizeof(fec), "FiniteElementCollection: H1_%" PetscInt_FMT "D_P1", dim));
261: PetscCall(PetscMalloc6(dof, &fec_type, dof, &nlocal, dof, &bss, dof, &dims, dof, &fieldname, dof, &Ufield));
262: for (i = 0; i < dof; i++) bss[i] = 1;
263: nf = dof;
265: PetscOptionsBegin(PetscObjectComm(oda), oda->prefix, "GLVis PetscViewer DMDA Field options", "PetscViewer");
266: PetscCall(PetscOptionsIntArray("-viewer_glvis_dm_da_bs", "Block sizes for subfields; enable vector representation", NULL, bss, &nf, &bsset));
267: PetscOptionsEnd();
268: if (bsset) {
269: PetscInt t;
270: for (i = 0, t = 0; i < nf; i++) t += bss[i];
271: PetscCheck(t == dof, PetscObjectComm(oda), PETSC_ERR_USER, "Sum of block sizes %" PetscInt_FMT " should equal %" PetscInt_FMT, t, dof);
272: } else nf = dof;
274: for (i = 0, s = 0; i < nf; i++) {
275: PetscCall(PetscStrallocpy(fec, &fec_type[i]));
276: if (bss[i] == 1) {
277: PetscCall(PetscStrallocpy(dafieldname[s], &fieldname[i]));
278: } else {
279: const char prefix[] = "Vector-";
280: const size_t prefix_len = PETSC_STATIC_ARRAY_LENGTH(prefix);
281: const PetscInt bss_i = bss[i];
282: char *fname_i = NULL;
283: size_t tlen = prefix_len, cur_len;
285: for (PetscInt b = 0; b < bss_i; ++b) {
286: size_t len;
288: PetscCall(PetscStrlen(dafieldname[s + b], &len));
289: tlen += len + 1; /* field + "-" */
290: }
291: PetscCall(PetscMalloc1(tlen, &fname_i));
292: PetscCall(PetscArraycpy(fname_i, prefix, prefix_len - 1));
293: cur_len = prefix_len;
294: for (PetscInt b = 0; b < bss_i; ++b) {
295: const char *dafname = dafieldname[s + b];
296: size_t len;
298: PetscCall(PetscStrlen(dafname, &len));
299: PetscCall(PetscArraycpy(fname_i + cur_len, dafname, len + 1));
300: cur_len += len + 1;
301: // don't add for final iteration of the loop
302: if ((b + 1) < bss_i) fname_i[cur_len++] = '-';
303: }
304: fname_i[cur_len] = '\0';
305: fieldname[i] = fname_i;
306: }
307: dims[i] = dim;
308: nlocal[i] = M * N * P * bss[i];
309: s += bss[i];
310: }
312: /* the viewer context takes ownership of xlocal and destroys it in DMDAFieldDestroyGLVisViewerCtx_Private */
313: PetscCall(PetscNew(&ctx));
314: ctx->xlocal = xlocal;
316: /* create work vectors */
317: for (i = 0; i < nf; i++) {
318: PetscCall(VecCreateMPI(PetscObjectComm((PetscObject)da), nlocal[i], PETSC_DECIDE, &Ufield[i]));
319: PetscCall(PetscObjectSetName((PetscObject)Ufield[i], fieldname[i]));
320: PetscCall(VecSetBlockSize(Ufield[i], PetscMax(bss[i], 1)));
321: PetscCall(VecSetDM(Ufield[i], da));
322: }
324: PetscCall(PetscViewerGLVisSetFields(viewer, nf, (const char **)fec_type, dims, DMDASampleGLVisFields_Private, (PetscObject *)Ufield, ctx, DMDAFieldDestroyGLVisViewerCtx_Private));
325: for (i = 0; i < nf; i++) {
326: PetscCall(PetscFree(fec_type[i]));
327: PetscCall(PetscFree(fieldname[i]));
328: PetscCall(VecDestroy(&Ufield[i]));
329: }
330: PetscCall(PetscFree6(fec_type, nlocal, bss, dims, fieldname, Ufield));
331: }
332: PetscFunctionReturn(PETSC_SUCCESS);
333: }
335: static PetscErrorCode DMDAView_GLVis_ASCII(DM dm, PetscViewer viewer)
336: {
337: DM da, cda;
338: Vec xcoorl;
339: PetscMPIInt size;
340: const PetscScalar *array;
341: PetscContainer glvis_container;
342: PetscInt dim, sdim, i, vid[8], mid, cid, cdof;
343: PetscInt sx, sy, sz, ie, je, ke, ien, jen, ken, nel;
344: PetscInt gsx, gsy, gsz, gm, gn, gp, kst, jst, ist;
345: PetscBool enabled = PETSC_TRUE, isascii;
346: const char *fmt;
348: PetscFunctionBegin;
351: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
352: PetscCheck(isascii, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Viewer must be of type VIEWERASCII");
353: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)viewer), &size));
354: PetscCheck(size <= 1, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Use single sequential viewers for parallel visualization");
355: PetscCall(DMGetDimension(dm, &dim));
357: /* get container: determines if a process visualizes is portion of the data or not */
358: PetscCall(PetscObjectQuery((PetscObject)viewer, "_glvis_info_container", (PetscObject *)&glvis_container));
359: PetscCheck(glvis_container, PetscObjectComm((PetscObject)dm), PETSC_ERR_PLIB, "Missing GLVis container");
360: {
361: PetscViewerGLVisInfo glvis_info;
362: PetscCall(PetscContainerGetPointer(glvis_container, (void **)&glvis_info));
363: enabled = glvis_info->enabled;
364: fmt = glvis_info->fmt;
365: }
366: /* this can happen if we are calling DMView outside of VecView_GLVis */
367: PetscCall(PetscObjectQuery((PetscObject)dm, "GLVisGraphicsDMDAGhosted", (PetscObject *)&da));
368: if (!da) PetscCall(DMSetUpGLVisViewer_DMDA((PetscObject)dm, NULL));
369: PetscCall(PetscObjectQuery((PetscObject)dm, "GLVisGraphicsDMDAGhosted", (PetscObject *)&da));
370: PetscCheck(da, PetscObjectComm((PetscObject)dm), PETSC_ERR_PLIB, "Missing GLVis ghosted DMDA");
371: PetscCall(DMGetCoordinateDim(da, &sdim));
373: PetscCall(PetscViewerASCIIPrintf(viewer, "MFEM mesh v1.0\n"));
374: PetscCall(PetscViewerASCIIPrintf(viewer, "\ndimension\n"));
375: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT "\n", dim));
377: if (!enabled) {
378: PetscCall(PetscViewerASCIIPrintf(viewer, "\nelements\n"));
379: PetscCall(PetscViewerASCIIPrintf(viewer, "0\n"));
380: PetscCall(PetscViewerASCIIPrintf(viewer, "\nboundary\n"));
381: PetscCall(PetscViewerASCIIPrintf(viewer, "0\n"));
382: PetscCall(PetscViewerASCIIPrintf(viewer, "\nvertices\n"));
383: PetscCall(PetscViewerASCIIPrintf(viewer, "0\n"));
384: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT "\n", sdim));
385: PetscFunctionReturn(PETSC_SUCCESS);
386: }
388: PetscCall(DMDAGetNumElementsGhosted(da, &ien, &jen, &ken));
389: nel = ien;
390: if (dim > 1) nel *= jen;
391: if (dim > 2) nel *= ken;
392: PetscCall(PetscViewerASCIIPrintf(viewer, "\nelements\n"));
393: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT "\n", nel));
394: switch (dim) {
395: case 1:
396: for (ie = 0; ie < ien; ie++) {
397: vid[0] = ie;
398: vid[1] = ie + 1;
399: mid = 1; /* material id */
400: cid = 1; /* segment */
401: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", mid, cid, vid[0], vid[1]));
402: }
403: break;
404: case 2:
405: for (je = 0; je < jen; je++) {
406: for (ie = 0; ie < ien; ie++) {
407: vid[0] = je * (ien + 1) + ie;
408: vid[1] = je * (ien + 1) + ie + 1;
409: vid[2] = (je + 1) * (ien + 1) + ie + 1;
410: vid[3] = (je + 1) * (ien + 1) + ie;
411: mid = 1; /* material id */
412: cid = 3; /* quad */
413: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", mid, cid, vid[0], vid[1], vid[2], vid[3]));
414: }
415: }
416: break;
417: case 3:
418: for (ke = 0; ke < ken; ke++) {
419: for (je = 0; je < jen; je++) {
420: for (ie = 0; ie < ien; ie++) {
421: vid[0] = ke * (jen + 1) * (ien + 1) + je * (ien + 1) + ie;
422: vid[1] = ke * (jen + 1) * (ien + 1) + je * (ien + 1) + ie + 1;
423: vid[2] = ke * (jen + 1) * (ien + 1) + (je + 1) * (ien + 1) + ie + 1;
424: vid[3] = ke * (jen + 1) * (ien + 1) + (je + 1) * (ien + 1) + ie;
425: vid[4] = (ke + 1) * (jen + 1) * (ien + 1) + je * (ien + 1) + ie;
426: vid[5] = (ke + 1) * (jen + 1) * (ien + 1) + je * (ien + 1) + ie + 1;
427: vid[6] = (ke + 1) * (jen + 1) * (ien + 1) + (je + 1) * (ien + 1) + ie + 1;
428: vid[7] = (ke + 1) * (jen + 1) * (ien + 1) + (je + 1) * (ien + 1) + ie;
429: mid = 1; /* material id */
430: cid = 5; /* hex */
431: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", mid, cid, vid[0], vid[1], vid[2], vid[3], vid[4], vid[5], vid[6], vid[7]));
432: }
433: }
434: }
435: break;
436: default:
437: SETERRQ(PetscObjectComm((PetscObject)da), PETSC_ERR_SUP, "Unsupported dimension %" PetscInt_FMT, dim);
438: }
439: PetscCall(PetscViewerASCIIPrintf(viewer, "\nboundary\n"));
440: PetscCall(PetscViewerASCIIPrintf(viewer, "0\n"));
442: /* vertex coordinates */
443: PetscCall(PetscObjectQuery((PetscObject)da, "GLVisGraphicsCoordsGhosted", (PetscObject *)&xcoorl));
444: PetscCheck(xcoorl, PetscObjectComm((PetscObject)dm), PETSC_ERR_PLIB, "Missing GLVis ghosted coords");
445: PetscCall(DMDAGetNumVerticesGhosted(da, &ien, &jen, &ken));
446: PetscCall(PetscViewerASCIIPrintf(viewer, "\nvertices\n"));
447: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT "\n", ien * jen * ken));
448: if (nel) {
449: PetscCall(VecGetDM(xcoorl, &cda));
450: PetscCall(VecGetArrayRead(xcoorl, &array));
451: if (!cda) { /* HO viz */
452: const char *fecname;
453: PetscInt nc, nl;
455: PetscCall(PetscObjectGetName((PetscObject)xcoorl, &fecname));
456: PetscCall(PetscViewerASCIIPrintf(viewer, "nodes\n"));
457: PetscCall(PetscViewerASCIIPrintf(viewer, "FiniteElementSpace\n"));
458: PetscCall(PetscViewerASCIIPrintf(viewer, "%s\n", fecname));
459: PetscCall(PetscViewerASCIIPrintf(viewer, "VDim: %" PetscInt_FMT "\n", sdim));
460: PetscCall(PetscViewerASCIIPrintf(viewer, "Ordering: 1\n\n")); /*Ordering::byVDIM*/
461: /* L2 coordinates */
462: PetscCall(DMDAGetNumElementsGhosted(da, &ien, &jen, &ken));
463: PetscCall(VecGetLocalSize(xcoorl, &nl));
464: nc = ien * (jen > 0 ? jen : 1) * (ken > 0 ? ken : 1);
465: cdof = nc ? nl / nc : 0;
466: if (!ien) ien++;
467: if (!jen) jen++;
468: if (!ken) ken++;
469: ist = jst = kst = 0;
470: gm = ien;
471: gn = jen;
472: gp = ken;
473: } else {
474: DMDAGhostedGLVisViewerCtx *dactx;
476: PetscCall(DMGetApplicationContext(da, &dactx));
477: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT "\n", sdim));
478: cdof = sdim;
479: PetscCall(DMDAGetCorners(da, &sx, &sy, &sz, NULL, NULL, NULL));
480: PetscCall(DMDAGetGhostCorners(da, &gsx, &gsy, &gsz, &gm, &gn, &gp));
481: if (dactx->ll) {
482: kst = jst = ist = 0;
483: } else {
484: kst = gsz != sz ? 1 : 0;
485: jst = gsy != sy ? 1 : 0;
486: ist = gsx != sx ? 1 : 0;
487: }
488: }
489: for (ke = kst; ke < kst + ken; ke++) {
490: for (je = jst; je < jst + jen; je++) {
491: for (ie = ist; ie < ist + ien; ie++) {
492: PetscInt c;
494: i = ke * gm * gn + je * gm + ie;
495: for (c = 0; c < cdof / sdim; c++) {
496: PetscInt d;
497: for (d = 0; d < sdim; d++) PetscCall(PetscViewerASCIIPrintf(viewer, fmt, PetscRealPart(array[cdof * i + c * sdim + d])));
498: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
499: }
500: }
501: }
502: }
503: PetscCall(VecRestoreArrayRead(xcoorl, &array));
504: }
505: PetscFunctionReturn(PETSC_SUCCESS);
506: }
508: PetscErrorCode DMView_DA_GLVis(DM dm, PetscViewer viewer)
509: {
510: PetscFunctionBegin;
511: PetscCall(DMView_GLVis(dm, viewer, DMDAView_GLVis_ASCII));
512: PetscFunctionReturn(PETSC_SUCCESS);
513: }