Actual source code: glvis.c
1: #define PETSC_DESIRE_FEATURE_TEST_MACROS /* for fdopen() */
3: #include <petsc/private/viewerimpl.h>
4: #include <petsc/private/petscimpl.h>
5: #include <petsc/private/glvisviewerimpl.h>
7: /* we may eventually make this function public */
8: static PetscErrorCode PetscViewerASCIISocketOpen(MPI_Comm,const char*,PetscInt,PetscViewer*);
10: struct _n_PetscViewerGLVis {
11: PetscViewerGLVisStatus status;
12: PetscViewerGLVisType type; /* either PETSC_VIEWER_GLVIS_DUMP or PETSC_VIEWER_GLVIS_SOCKET */
13: char *name; /* prefix for filename, or hostname, depending on the type */
14: PetscInt port; /* used just for the socket case */
15: PetscReal pause; /* if positive, calls PetscSleep(pause) after each VecView_GLVis call */
16: PetscViewer meshwindow; /* used just by the ASCII dumping */
17: PetscObject dm; /* DM as passed by PetscViewerGLVisSetDM_Private(): should contain discretization info */
18: PetscInt nwindow; /* number of windows/fields to be visualized */
19: PetscViewer *window;
20: char **windowtitle;
21: PetscInt windowsizes[2];
22: char **fec_type; /* type of elements to be used for visualization, see FiniteElementCollection::Name() */
23: PetscErrorCode (*g2lfield)(PetscObject,PetscInt,PetscObject[],void*); /* global to local operation for generating dofs to be visualized */
24: PetscInt *spacedim; /* geometrical space dimension (just used to initialize the scene) */
25: PetscObject *Ufield; /* work vectors for visualization */
26: PetscInt snapid; /* snapshot id, use PetscViewerGLVisSetSnapId to change this value*/
27: void *userctx; /* User context, used by g2lfield */
28: PetscErrorCode (*destroyctx)(void*); /* destroy routine for userctx */
29: char* fmt; /* format string for FP values */
30: };
31: typedef struct _n_PetscViewerGLVis *PetscViewerGLVis;
33: /*@
34: PetscViewerGLVisSetPrecision - Set the number of digits for floating point values
36: Not Collective
38: Input Parameters:
39: + viewer - the PetscViewer
40: - prec - the number of digits required
42: Level: beginner
44: .seealso: PetscViewerGLVisOpen(), PetscViewerGLVisSetFields(), PetscViewerCreate(), PetscViewerSetType()
45: @*/
46: PetscErrorCode PetscViewerGLVisSetPrecision(PetscViewer viewer, PetscInt prec)
47: {
49: PetscTryMethod(viewer,"PetscViewerGLVisSetPrecision_C",(PetscViewer,PetscInt),(viewer,prec));
50: return 0;
51: }
53: static PetscErrorCode PetscViewerGLVisSetPrecision_GLVis(PetscViewer viewer, PetscInt prec)
54: {
55: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
57: PetscFree(socket->fmt);
58: if (prec > 0) {
59: PetscMalloc1(16,&socket->fmt);
60: PetscSNPrintf(socket->fmt,16," %%.%" PetscInt_FMT "e",prec);
61: } else {
62: PetscStrallocpy(" %g",&socket->fmt);
63: }
64: return 0;
65: }
67: /*@
68: PetscViewerGLVisSetSnapId - Set the snapshot id. Only relevant when the viewer is of type PETSC_VIEWER_GLVIS_DUMP
70: Logically Collective on PetscViewer
72: Input Parameters:
73: + viewer - the PetscViewer
74: - id - the current snapshot id in a time-dependent simulation
76: Level: beginner
78: .seealso: PetscViewerGLVisOpen(), PetscViewerGLVisSetFields(), PetscViewerCreate(), PetscViewerSetType()
79: @*/
80: PetscErrorCode PetscViewerGLVisSetSnapId(PetscViewer viewer, PetscInt id)
81: {
84: PetscTryMethod(viewer,"PetscViewerGLVisSetSnapId_C",(PetscViewer,PetscInt),(viewer,id));
85: return 0;
86: }
88: static PetscErrorCode PetscViewerGLVisSetSnapId_GLVis(PetscViewer viewer, PetscInt id)
89: {
90: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
92: socket->snapid = id;
93: return 0;
94: }
96: /*@C
97: PetscViewerGLVisSetFields - Sets the required information to visualize different fields from a vector.
99: Logically Collective on PetscViewer
101: Input Parameters:
102: + viewer - the PetscViewer
103: . nf - number of fields to be visualized
104: . fec_type - the type of finite element to be used to visualize the data (see FiniteElementCollection::Name() in MFEM)
105: . dim - array of space dimension for field vectors (used to initialize the scene)
106: . g2lfields - User routine to compute the local field vectors to be visualized; PetscObject is used in place of Vec on the prototype
107: . Vfield - array of work vectors, one for each field
108: . ctx - User context to store the relevant data to apply g2lfields
109: - destroyctx - Destroy function for userctx
111: Notes:
112: g2lfields is called on the vector V to be visualized in order to extract the relevant dofs to be put in Vfield[], as
113: .vb
114: g2lfields((PetscObject)V,nfields,(PetscObject*)Vfield[],ctx).
115: .ve
116: For vector spaces, the block size of Vfield[i] represents the vector dimension. It misses the Fortran bindings.
117: The names of the Vfield vectors will be displayed in the window title.
119: Level: intermediate
121: .seealso: PetscViewerGLVisOpen(), PetscViewerCreate(), PetscViewerSetType(), PetscObjectSetName()
122: @*/
123: PetscErrorCode PetscViewerGLVisSetFields(PetscViewer viewer, PetscInt nf, const char* fec_type[], PetscInt dim[], PetscErrorCode(*g2l)(PetscObject,PetscInt,PetscObject[],void*), PetscObject Vfield[], void* ctx, PetscErrorCode(*destroyctx)(void*))
124: {
131: PetscTryMethod(viewer,"PetscViewerGLVisSetFields_C",(PetscViewer,PetscInt,const char*[],PetscInt[],PetscErrorCode(*)(PetscObject,PetscInt,PetscObject[],void*),PetscObject[],void*,PetscErrorCode(*)(void*)),(viewer,nf,fec_type,dim,g2l,Vfield,ctx,destroyctx));
132: return 0;
133: }
135: static PetscErrorCode PetscViewerGLVisSetFields_GLVis(PetscViewer viewer, PetscInt nfields, const char* fec_type[], PetscInt dim[], PetscErrorCode(*g2l)(PetscObject,PetscInt,PetscObject[],void*), PetscObject Vfield[], void* ctx, PetscErrorCode(*destroyctx)(void*))
136: {
137: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
138: PetscInt i;
141: if (!socket->nwindow) {
142: socket->nwindow = nfields;
144: PetscCalloc5(nfields,&socket->window,nfields,&socket->windowtitle,nfields,&socket->fec_type,nfields,&socket->spacedim,nfields,&socket->Ufield);
145: for (i=0;i<nfields;i++) {
146: const char *name;
148: PetscObjectGetName(Vfield[i],&name);
149: PetscStrallocpy(name,&socket->windowtitle[i]);
150: PetscStrallocpy(fec_type[i],&socket->fec_type[i]);
151: PetscObjectReference(Vfield[i]);
152: socket->Ufield[i] = Vfield[i];
153: socket->spacedim[i] = dim[i];
154: }
155: }
156: /* number of fields are not allowed to vary */
158: socket->g2lfield = g2l;
159: if (socket->destroyctx && socket->userctx) (*socket->destroyctx)(socket->userctx);
160: socket->userctx = ctx;
161: socket->destroyctx = destroyctx;
162: return 0;
163: }
165: static PetscErrorCode PetscViewerGLVisInfoDestroy_Private(void *ptr)
166: {
167: PetscViewerGLVisInfo info = (PetscViewerGLVisInfo)ptr;
169: PetscFree(info->fmt);
170: PetscFree(info);
171: return 0;
172: }
174: /* we can decide to prevent specific processes from using the viewer */
175: static PetscErrorCode PetscViewerGLVisAttachInfo_Private(PetscViewer viewer, PetscViewer window)
176: {
177: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
178: PetscContainer container;
179: PetscViewerGLVisInfo info;
181: PetscObjectQuery((PetscObject)window,"_glvis_info_container",(PetscObject*)&container);
182: if (!container) {
183: PetscNew(&info);
184: info->enabled = PETSC_TRUE;
185: info->init = PETSC_FALSE;
186: info->size[0] = socket->windowsizes[0];
187: info->size[1] = socket->windowsizes[1];
188: info->pause = socket->pause;
189: PetscContainerCreate(PetscObjectComm((PetscObject)window),&container);
190: PetscContainerSetPointer(container,(void*)info);
191: PetscContainerSetUserDestroy(container,PetscViewerGLVisInfoDestroy_Private);
192: PetscObjectCompose((PetscObject)window,"_glvis_info_container",(PetscObject)container);
193: PetscContainerDestroy(&container);
194: } else {
195: PetscContainerGetPointer(container,(void**)&info);
196: }
197: PetscFree(info->fmt);
198: PetscStrallocpy(socket->fmt,&info->fmt);
199: return 0;
200: }
202: static PetscErrorCode PetscViewerGLVisGetNewWindow_Private(PetscViewer viewer,PetscViewer *view)
203: {
204: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
205: PetscViewer window = NULL;
206: PetscBool ldis,dis;
207: PetscErrorCode ierr;
209: PetscViewerASCIISocketOpen(PETSC_COMM_SELF,socket->name,socket->port,&window);
210: /* if we could not estabilish a connection the first time,
211: we disable the socket viewer */
212: ldis = ierr ? PETSC_TRUE : PETSC_FALSE;
213: MPIU_Allreduce(&ldis,&dis,1,MPIU_BOOL,MPI_LOR,PetscObjectComm((PetscObject)viewer));
214: if (dis) {
215: socket->status = PETSCVIEWERGLVIS_DISABLED;
216: PetscViewerDestroy(&window);
217: }
218: *view = window;
219: return 0;
220: }
222: PetscErrorCode PetscViewerGLVisPause_Private(PetscViewer viewer)
223: {
224: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
226: if (socket->type == PETSC_VIEWER_GLVIS_SOCKET && socket->pause > 0) {
227: PetscSleep(socket->pause);
228: }
229: return 0;
230: }
232: /* DM specific support */
233: PetscErrorCode PetscViewerGLVisSetDM_Private(PetscViewer viewer, PetscObject dm)
234: {
235: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
238: if (!socket->dm) {
239: PetscErrorCode (*setupwithdm)(PetscObject,PetscViewer) = NULL;
241: PetscObjectQueryFunction(dm,"DMSetUpGLVisViewer_C",&setupwithdm);
242: if (setupwithdm) {
243: (*setupwithdm)(dm,viewer);
244: } else SETERRQ(PetscObjectComm(dm),PETSC_ERR_SUP,"No support for DM type %s",dm->type_name);
245: PetscObjectReference(dm);
246: socket->dm = dm;
247: }
248: return 0;
249: }
251: PetscErrorCode PetscViewerGLVisGetDMWindow_Private(PetscViewer viewer,PetscViewer *view)
252: {
253: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
256: if (!socket->meshwindow) {
257: if (socket->type == PETSC_VIEWER_GLVIS_SOCKET) {
258: PetscViewerGLVisGetNewWindow_Private(viewer,&socket->meshwindow);
259: } else {
260: size_t len;
261: PetscBool isstdout;
263: PetscStrlen(socket->name,&len);
264: PetscStrcmp(socket->name,"stdout",&isstdout);
265: if (!socket->name || !len || isstdout) {
266: PetscViewerASCIIOpen(PETSC_COMM_SELF,"stdout",&socket->meshwindow);
267: } else {
268: PetscMPIInt rank;
269: char filename[PETSC_MAX_PATH_LEN];
270: MPI_Comm_rank(PetscObjectComm((PetscObject)viewer),&rank);
271: PetscSNPrintf(filename,PETSC_MAX_PATH_LEN,"%s-mesh.%06d",socket->name,rank);
272: PetscViewerASCIIOpen(PETSC_COMM_SELF,filename,&socket->meshwindow);
273: }
274: }
275: if (socket->meshwindow) {
276: PetscViewerPushFormat(socket->meshwindow,PETSC_VIEWER_ASCII_GLVIS);
277: }
278: }
279: if (socket->meshwindow) {
280: PetscViewerGLVisAttachInfo_Private(viewer,socket->meshwindow);
281: }
282: *view = socket->meshwindow;
283: return 0;
284: }
286: PetscErrorCode PetscViewerGLVisRestoreDMWindow_Private(PetscViewer viewer,PetscViewer *view)
287: {
288: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
292: if (*view) {
293: PetscViewerFlush(*view);
294: PetscBarrier((PetscObject)viewer);
295: }
296: if (socket->type == PETSC_VIEWER_GLVIS_DUMP) { /* destroy the viewer, as it is associated with a single time step */
297: PetscViewerDestroy(&socket->meshwindow);
298: } else if (!*view) { /* something went wrong (SIGPIPE) so we just zero the private pointer */
299: socket->meshwindow = NULL;
300: }
301: *view = NULL;
302: return 0;
303: }
305: PetscErrorCode PetscViewerGLVisGetType_Private(PetscViewer viewer,PetscViewerGLVisType *type)
306: {
307: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
310: *type = socket->type;
311: return 0;
312: }
314: /* This function is only relevant in the SOCKET_GLIVS case. The status is computed the first time it is requested, as GLVis currently has issues when connecting the first time through the socket */
315: PetscErrorCode PetscViewerGLVisGetStatus_Private(PetscViewer viewer, PetscViewerGLVisStatus *sockstatus)
316: {
317: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
320: if (socket->type == PETSC_VIEWER_GLVIS_DUMP) {
321: socket->status = PETSCVIEWERGLVIS_DISCONNECTED;
322: } else if (socket->status == PETSCVIEWERGLVIS_DISCONNECTED && socket->nwindow) {
323: PetscInt i;
324: PetscBool lconn,conn;
326: for (i=0,lconn=PETSC_TRUE;i<socket->nwindow;i++)
327: if (!socket->window[i])
328: lconn = PETSC_FALSE;
330: MPIU_Allreduce(&lconn,&conn,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)viewer));
331: if (conn) socket->status = PETSCVIEWERGLVIS_CONNECTED;
332: }
333: *sockstatus = socket->status;
334: return 0;
335: }
337: PetscErrorCode PetscViewerGLVisGetDM_Private(PetscViewer viewer, PetscObject* dm)
338: {
339: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
341: *dm = socket->dm;
342: return 0;
343: }
345: PetscErrorCode PetscViewerGLVisGetFields_Private(PetscViewer viewer, PetscInt* nfield, const char **fec[], PetscInt *spacedim[], PetscErrorCode(**g2lfield)(PetscObject,PetscInt,PetscObject[],void*), PetscObject *Ufield[], void **userctx)
346: {
347: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
349: if (nfield) *nfield = socket->nwindow;
350: if (fec) *fec = (const char**)socket->fec_type;
351: if (spacedim) *spacedim = socket->spacedim;
352: if (g2lfield) *g2lfield = socket->g2lfield;
353: if (Ufield) *Ufield = socket->Ufield;
354: if (userctx) *userctx = socket->userctx;
355: return 0;
356: }
358: /* accessor routines for the viewer windows:
359: PETSC_VIEWER_GLVIS_DUMP : it returns a new viewer every time
360: PETSC_VIEWER_GLVIS_SOCKET : it returns the socket, and creates it if not yet done.
361: */
362: PetscErrorCode PetscViewerGLVisGetWindow_Private(PetscViewer viewer,PetscInt wid,PetscViewer* view)
363: {
364: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
365: PetscViewerGLVisStatus status;
370: status = socket->status;
372: switch (status) {
373: case PETSCVIEWERGLVIS_DISCONNECTED:
375: else if (socket->type == PETSC_VIEWER_GLVIS_DUMP) {
376: size_t len;
377: PetscBool isstdout;
379: PetscStrlen(socket->name,&len);
380: PetscStrcmp(socket->name,"stdout",&isstdout);
381: if (!socket->name || !len || isstdout) {
382: PetscViewerASCIIOpen(PETSC_COMM_SELF,"stdout",&socket->window[wid]);
383: } else {
384: PetscMPIInt rank;
385: char filename[PETSC_MAX_PATH_LEN];
387: MPI_Comm_rank(PetscObjectComm((PetscObject)viewer),&rank);
388: PetscSNPrintf(filename,PETSC_MAX_PATH_LEN,"%s-%s-%" PetscInt_FMT ".%06d",socket->name,socket->windowtitle[wid],socket->snapid,rank);
389: PetscViewerASCIIOpen(PETSC_COMM_SELF,filename,&socket->window[wid]);
390: }
391: } else {
392: PetscViewerGLVisGetNewWindow_Private(viewer,&socket->window[wid]);
393: }
394: if (socket->window[wid]) {
395: PetscViewerPushFormat(socket->window[wid],PETSC_VIEWER_ASCII_GLVIS);
396: }
397: *view = socket->window[wid];
398: break;
399: case PETSCVIEWERGLVIS_CONNECTED:
400: *view = socket->window[wid];
401: break;
402: case PETSCVIEWERGLVIS_DISABLED:
403: *view = NULL;
404: break;
405: default:
406: SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Unhandled socket status %d",(int)status);
407: }
408: if (*view) {
409: PetscViewerGLVisAttachInfo_Private(viewer,*view);
410: }
411: return 0;
412: }
414: /* Restore the window viewer
415: PETSC_VIEWER_GLVIS_DUMP : destroys the temporary created ASCII viewer used for dumping
416: PETSC_VIEWER_GLVIS_SOCKET: - if the returned window viewer is not NULL, just zeros the pointer.
417: - it the returned window viewer is NULL, assumes something went wrong
418: with the socket (i.e. SIGPIPE when a user closes the popup window)
419: and that the caller already handled it (see VecView_GLVis).
420: */
421: PetscErrorCode PetscViewerGLVisRestoreWindow_Private(PetscViewer viewer,PetscInt wid, PetscViewer* view)
422: {
423: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
430: if (*view) {
431: PetscViewerFlush(*view);
432: PetscBarrier((PetscObject)viewer);
433: }
434: if (socket->type == PETSC_VIEWER_GLVIS_DUMP) { /* destroy the viewer, as it is associated with a single time step */
435: PetscViewerDestroy(&socket->window[wid]);
436: } else if (!*view) { /* something went wrong (SIGPIPE) so we just zero the private pointer */
437: socket->window[wid] = NULL;
438: }
439: *view = NULL;
440: return 0;
441: }
443: /* default window appearance in the PETSC_VIEWER_GLVIS_SOCKET case */
444: PetscErrorCode PetscViewerGLVisInitWindow_Private(PetscViewer viewer, PetscBool mesh, PetscInt dim, const char *name)
445: {
446: PetscViewerGLVisInfo info;
447: PetscContainer container;
449: PetscObjectQuery((PetscObject)viewer,"_glvis_info_container",(PetscObject*)&container);
451: PetscContainerGetPointer(container,(void**)&info);
452: if (info->init) return 0;
454: /* Configure window */
455: if (info->size[0] > 0) {
456: PetscViewerASCIIPrintf(viewer,"window_size %" PetscInt_FMT " %" PetscInt_FMT "\n",info->size[0],info->size[1]);
457: }
458: if (name) {
459: PetscViewerASCIIPrintf(viewer,"window_title '%s'\n",name);
460: }
462: /* Configure default view */
463: if (mesh) {
464: switch (dim) {
465: case 1:
466: PetscViewerASCIIPrintf(viewer,"keys m\n"); /* show mesh */
467: break;
468: case 2:
469: PetscViewerASCIIPrintf(viewer,"keys m\n"); /* show mesh */
470: break;
471: case 3: /* TODO: decide default view in 3D */
472: break;
473: }
474: } else {
475: PetscViewerASCIIPrintf(viewer,"keys cm\n"); /* show colorbar and mesh */
476: switch (dim) {
477: case 1:
478: PetscViewerASCIIPrintf(viewer,"keys RRjl\n"); /* set to 1D (side view), turn off perspective and light */
479: break;
480: case 2:
481: PetscViewerASCIIPrintf(viewer,"keys Rjl\n"); /* set to 2D (top view), turn off perspective and light */
482: break;
483: case 3:
484: break;
485: }
486: PetscViewerASCIIPrintf(viewer,"autoscale value\n"); /* update value-range; keep mesh-extents fixed */
487: }
489: { /* Additional keys and commands */
490: char keys[256] = "", cmds[2*PETSC_MAX_PATH_LEN] = "";
491: PetscOptions opt = ((PetscObject)viewer)->options;
492: const char *pre = ((PetscObject)viewer)->prefix;
494: PetscOptionsGetString(opt,pre,"-glvis_keys",keys,sizeof(keys),NULL);
495: PetscOptionsGetString(opt,pre,"-glvis_exec",cmds,sizeof(cmds),NULL);
496: if (keys[0]) PetscViewerASCIIPrintf(viewer,"keys %s\n",keys);
497: if (cmds[0]) PetscViewerASCIIPrintf(viewer,"%s\n",cmds);
498: }
500: /* Pause visualization */
501: if (!mesh && info->pause == -1) {
502: PetscViewerASCIIPrintf(viewer,"autopause 1\n");
503: }
504: if (!mesh && info->pause == 0) {
505: PetscViewerASCIIPrintf(viewer,"pause\n");
506: }
508: info->init = PETSC_TRUE;
509: return 0;
510: }
512: static PetscErrorCode PetscViewerDestroy_GLVis(PetscViewer viewer)
513: {
514: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
515: PetscInt i;
517: for (i=0;i<socket->nwindow;i++) {
518: PetscViewerDestroy(&socket->window[i]);
519: PetscFree(socket->windowtitle[i]);
520: PetscFree(socket->fec_type[i]);
521: PetscObjectDestroy(&socket->Ufield[i]);
522: }
523: PetscFree(socket->name);
524: PetscFree5(socket->window,socket->windowtitle,socket->fec_type,socket->spacedim,socket->Ufield);
525: PetscFree(socket->fmt);
526: PetscViewerDestroy(&socket->meshwindow);
527: PetscObjectDestroy(&socket->dm);
528: if (socket->destroyctx && socket->userctx) (*socket->destroyctx)(socket->userctx);
530: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetPrecision_C",NULL);
531: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetSnapId_C",NULL);
532: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetFields_C",NULL);
533: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileSetName_C",NULL);
534: PetscFree(socket);
535: viewer->data = NULL;
536: return 0;
537: }
539: static PetscErrorCode PetscViewerSetFromOptions_GLVis(PetscOptionItems *PetscOptionsObject,PetscViewer v)
540: {
541: PetscViewerGLVis socket = (PetscViewerGLVis)v->data;
542: PetscInt nsizes = 2, prec = PETSC_DECIDE;
543: PetscBool set;
545: PetscOptionsHead(PetscOptionsObject,"GLVis PetscViewer Options");
546: PetscOptionsInt("-glvis_precision","Number of digits for floating point values","PetscViewerGLVisSetPrecision",prec,&prec,&set);
547: if (set) PetscViewerGLVisSetPrecision(v,prec);
548: PetscOptionsIntArray("-glvis_size","Window sizes",NULL,socket->windowsizes,&nsizes,&set);
549: if (set && (nsizes == 1 || socket->windowsizes[1] < 0)) socket->windowsizes[1] = socket->windowsizes[0];
550: PetscOptionsReal("-glvis_pause","-1 to pause after each visualization, otherwise sleeps for given seconds",NULL,socket->pause,&socket->pause,NULL);
551: PetscOptionsName("-glvis_keys","Additional keys to configure visualization",NULL,NULL);
552: PetscOptionsName("-glvis_exec","Additional commands to configure visualization",NULL,NULL);
553: PetscOptionsTail();
554: return 0;
555: }
557: static PetscErrorCode PetscViewerFileSetName_GLVis(PetscViewer viewer, const char name[])
558: {
559: char *sport;
560: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
562: socket->type = PETSC_VIEWER_GLVIS_DUMP;
563: /* we accept localhost^port */
564: PetscFree(socket->name);
565: PetscStrallocpy(name,&socket->name);
566: PetscStrchr(socket->name,'^',&sport);
567: if (sport) {
568: PetscInt port = 19916;
569: size_t len;
572: *sport++ = 0;
573: PetscStrlen(sport,&len);
574: PetscOptionsStringToInt(sport,&port);
575: if (PetscUnlikely(ierr)) {
576: socket->port = 19916;
577: } else {
578: socket->port = (port != PETSC_DECIDE && port != PETSC_DEFAULT) ? port : 19916;
579: }
580: socket->type = PETSC_VIEWER_GLVIS_SOCKET;
581: }
582: return 0;
583: }
585: /*@C
586: PetscViewerGLVisOpen - Opens a GLVis type viewer
588: Collective
590: Input Parameters:
591: + comm - the MPI communicator
592: . type - the viewer type: PETSC_VIEWER_GLVIS_SOCKET for real-time visualization or PETSC_VIEWER_GLVIS_DUMP for dumping to disk
593: . name - either the hostname where the GLVis server is running or the base filename for dumping the data for subsequent visualizations
594: - port - socket port where the GLVis server is listening. Not referenced when type is PETSC_VIEWER_GLVIS_DUMP
596: Output Parameters:
597: - viewer - the PetscViewer object
599: Options Database Keys:
600: + -glvis_precision <precision> - Sets number of digits for floating point values
601: . -glvis_size <width,height> - Sets the window size (in pixels)
602: . -glvis_pause <pause> - Sets time (in seconds) that the program pauses after each visualization
603: (0 is default, -1 implies every visualization)
604: . -glvis_keys - Additional keys to configure visualization
605: - -glvis_exec - Additional commands to configure visualization
607: Notes:
608: misses Fortran binding
610: Level: beginner
612: .seealso: PetscViewerCreate(), PetscViewerSetType(), PetscViewerGLVisType
613: @*/
614: PetscErrorCode PetscViewerGLVisOpen(MPI_Comm comm, PetscViewerGLVisType type, const char name[], PetscInt port, PetscViewer *viewer)
615: {
616: PetscViewerGLVis socket;
618: PetscViewerCreate(comm,viewer);
619: PetscViewerSetType(*viewer,PETSCVIEWERGLVIS);
621: socket = (PetscViewerGLVis)((*viewer)->data);
622: socket->type = type;
623: if (type == PETSC_VIEWER_GLVIS_DUMP || name) {
624: PetscFree(socket->name);
625: PetscStrallocpy(name,&socket->name);
626: }
627: socket->port = (!port || port == PETSC_DETERMINE || port == PETSC_DECIDE) ? 19916 : port;
629: PetscViewerSetFromOptions(*viewer);
630: return 0;
631: }
633: /*
634: PETSC_VIEWER_GLVIS_ - Creates an GLVIS PetscViewer shared by all processors in a communicator.
636: Collective
638: Input Parameter:
639: . comm - the MPI communicator to share the GLVIS PetscViewer
641: Level: intermediate
643: Notes:
644: misses Fortran bindings
646: Environmental variables:
647: + PETSC_VIEWER_GLVIS_FILENAME : output filename (if specified dump to disk, and takes precedence on PETSC_VIEWER_GLVIS_HOSTNAME)
648: . PETSC_VIEWER_GLVIS_HOSTNAME : machine where the GLVis server is listening (defaults to localhost)
649: - PETSC_VIEWER_GLVIS_PORT : port opened by the GLVis server (defaults to 19916)
651: Notes:
652: Unlike almost all other PETSc routines, PETSC_VIEWER_GLVIS_ does not return
653: an error code. The GLVIS PetscViewer is usually used in the form
654: $ XXXView(XXX object, PETSC_VIEWER_GLVIS_(comm));
656: .seealso: PetscViewerGLVISOpen(), PetscViewerGLVisType, PetscViewerCreate(), PetscViewerDestroy()
657: */
658: PetscViewer PETSC_VIEWER_GLVIS_(MPI_Comm comm)
659: {
660: PetscErrorCode ierr;
661: PetscBool flg;
662: PetscViewer viewer;
663: PetscViewerGLVisType type;
664: char fname[PETSC_MAX_PATH_LEN],sport[16];
665: PetscInt port = 19916; /* default for GLVis */
667: PetscOptionsGetenv(comm,"PETSC_VIEWER_GLVIS_FILENAME",fname,PETSC_MAX_PATH_LEN,&flg);
668: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return NULL;}
669: if (!flg) {
670: type = PETSC_VIEWER_GLVIS_SOCKET;
671: PetscOptionsGetenv(comm,"PETSC_VIEWER_GLVIS_HOSTNAME",fname,PETSC_MAX_PATH_LEN,&flg);
672: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return NULL;}
673: if (!flg) {
674: PetscStrcpy(fname,"localhost");
675: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return NULL;}
676: }
677: PetscOptionsGetenv(comm,"PETSC_VIEWER_GLVIS_PORT",sport,16,&flg);
678: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return NULL;}
679: if (flg) {
680: PetscOptionsStringToInt(sport,&port);
681: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return NULL;}
682: }
683: } else {
684: type = PETSC_VIEWER_GLVIS_DUMP;
685: }
686: PetscViewerGLVisOpen(comm,type,fname,port,&viewer);
687: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return NULL;}
688: PetscObjectRegisterDestroy((PetscObject)viewer);
689: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return NULL;}
690: return viewer;
691: }
693: PETSC_EXTERN PetscErrorCode PetscViewerCreate_GLVis(PetscViewer viewer)
694: {
695: PetscViewerGLVis socket;
697: PetscNewLog(viewer,&socket);
699: /* defaults to socket viewer */
700: PetscStrallocpy("localhost",&socket->name);
701: socket->port = 19916; /* GLVis default listening port */
702: socket->type = PETSC_VIEWER_GLVIS_SOCKET;
703: socket->pause = 0; /* just pause the first time */
705: socket->windowsizes[0] = 600;
706: socket->windowsizes[1] = 600;
708: /* defaults to full precision */
709: PetscStrallocpy(" %g",&socket->fmt);
711: viewer->data = (void*)socket;
712: viewer->ops->destroy = PetscViewerDestroy_GLVis;
713: viewer->ops->setfromoptions = PetscViewerSetFromOptions_GLVis;
715: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetPrecision_C",PetscViewerGLVisSetPrecision_GLVis);
716: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetSnapId_C",PetscViewerGLVisSetSnapId_GLVis);
717: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetFields_C",PetscViewerGLVisSetFields_GLVis);
718: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileSetName_C",PetscViewerFileSetName_GLVis);
719: return 0;
720: }
722: /* this is a private implementation of a SOCKET with ASCII data format
723: GLVis does not currently handle binary socket streams */
724: #if defined(PETSC_HAVE_UNISTD_H)
725: #include <unistd.h>
726: #endif
728: #if !defined(PETSC_HAVE_WINDOWS_H)
729: static PetscErrorCode (*PetscViewerDestroy_ASCII)(PetscViewer);
731: static PetscErrorCode PetscViewerDestroy_ASCII_Socket(PetscViewer viewer)
732: {
733: FILE *stream;
734: PetscErrorCode 0;
735: PetscViewerASCIIGetPointer(viewer,&stream);
736: if (stream) {
737: fclose(stream);
739: }
740: PetscViewerDestroy_ASCII(viewer);
741: return 0;
742: }
743: #endif
745: static PetscErrorCode PetscViewerASCIISocketOpen(MPI_Comm comm,const char* hostname,PetscInt port,PetscViewer* viewer)
746: {
747: #if defined(PETSC_HAVE_WINDOWS_H)
748: SETERRQ(comm,PETSC_ERR_SUP,"Not implemented for Windows");
749: #else
750: FILE *stream = NULL;
751: int fd=0;
756: #if defined(PETSC_USE_SOCKET_VIEWER)
757: PetscOpenSocket(hostname,port,&fd);
758: #else
759: SETERRQ(comm,PETSC_ERR_SUP,"Missing Socket viewer");
760: #endif
761: if (PetscUnlikely(ierr)) {
762: PetscInt sierr;
763: char err[1024];
765: PetscSNPrintf(err,1024,"Cannot connect to socket on %s:%" PetscInt_FMT ". Socket visualization is disabled\n",hostname,port);
766: PetscInfo(NULL,"%s",err);
767: *viewer = NULL;
768: return sierr;
769: } else {
770: char msg[1024];
772: PetscSNPrintf(msg,1024,"Successfully connect to socket on %s:%" PetscInt_FMT ". Socket visualization is enabled\n",hostname,port);
773: PetscInfo(NULL,"%s",msg);
774: }
775: stream = fdopen(fd,"w"); /* Not possible on Windows */
777: PetscViewerASCIIOpenWithFILE(PETSC_COMM_SELF,stream,viewer);
778: PetscViewerDestroy_ASCII = (*viewer)->ops->destroy;
779: (*viewer)->ops->destroy = PetscViewerDestroy_ASCII_Socket;
780: #endif
781: return 0;
782: }
784: #if !defined(PETSC_MISSING_SIGPIPE)
786: #include <signal.h>
788: #if defined(PETSC_HAVE_WINDOWS_H)
789: #define PETSC_DEVNULL "NUL"
790: #else
791: #define PETSC_DEVNULL "/dev/null"
792: #endif
794: static volatile PetscBool PetscGLVisBrokenPipe = PETSC_FALSE;
796: static void (*PetscGLVisSigHandler_save)(int) = NULL;
798: static void PetscGLVisSigHandler_SIGPIPE(PETSC_UNUSED int sig)
799: {
800: PetscGLVisBrokenPipe = PETSC_TRUE;
801: #if !defined(PETSC_MISSING_SIG_IGN)
802: signal(SIGPIPE,SIG_IGN);
803: #endif
804: }
806: PetscErrorCode PetscGLVisCollectiveBegin(PETSC_UNUSED MPI_Comm comm,PETSC_UNUSED PetscViewer *win)
807: {
809: PetscGLVisBrokenPipe = PETSC_FALSE;
810: PetscGLVisSigHandler_save = signal(SIGPIPE,PetscGLVisSigHandler_SIGPIPE);
811: return 0;
812: }
814: PetscErrorCode PetscGLVisCollectiveEnd(MPI_Comm comm,PetscViewer *win)
815: {
816: PetscBool flag,brokenpipe;
818: flag = PetscGLVisBrokenPipe;
819: MPIU_Allreduce(&flag,&brokenpipe,1,MPIU_BOOL,MPI_LOR,comm);
820: if (brokenpipe) {
821: FILE *sock, *null = fopen(PETSC_DEVNULL,"w");
822: PetscViewerASCIIGetPointer(*win,&sock);
823: PetscViewerASCIISetFILE(*win,null);
824: PetscViewerDestroy(win);
825: if (sock) (void)fclose(sock);
826: }
827: (void)signal(SIGPIPE,PetscGLVisSigHandler_save);
828: PetscGLVisSigHandler_save = NULL;
829: PetscGLVisBrokenPipe = PETSC_FALSE;
830: return 0;
831: }
833: #else
835: PetscErrorCode PetscGLVisCollectiveBegin(PETSC_UNUSED MPI_Comm comm,PETSC_UNUSED PetscViewer *win)
836: {
837: return 0;
838: }
840: PetscErrorCode PetscGLVisCollectiveEnd(PETSC_UNUSED MPI_Comm comm,PETSC_UNUSED PetscViewer *win)
841: {
842: return 0;
843: }
845: #endif