Actual source code: glvis.c
petsc-3.8.4 2018-03-24
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: char **fec_type; /* type of elements to be used for visualization, see FiniteElementCollection::Name() */
22: PetscErrorCode (*g2lfield)(PetscObject,PetscInt,PetscObject[],void*); /* global to local operation for generating dofs to be visualized */
23: PetscInt *locandbs; /* local and block sizes for work vectors */
24: PetscObject *Ufield; /* work vectors for visualization */
25: PetscInt snapid; /* snapshot id, use PetscViewerGLVisSetSnapId to change this value*/
26: void *userctx; /* User context, used by g2lfield */
27: PetscErrorCode (*destroyctx)(void*); /* destroy routine for userctx */
28: char* fmt; /* format string for FP values */
29: };
30: typedef struct _n_PetscViewerGLVis *PetscViewerGLVis;
32: /*@
33: PetscViewerGLVisSetPrecision - Set the number of digits for floating point values
35: Not Collective
37: Input Parameters:
38: + viewer - the PetscViewer
39: - prec - the number of digits required
41: Level: beginner
43: .seealso: PetscViewerGLVisOpen(), PetscViewerGLVisSetFields(), PetscViewerCreate(), PetscViewerSetType()
44: @*/
45: PetscErrorCode PetscViewerGLVisSetPrecision(PetscViewer viewer, PetscInt prec)
46: {
51: PetscTryMethod(viewer,"PetscViewerGLVisSetPrecision_C",(PetscViewer,PetscInt),(viewer,prec));
52: return(0);
53: }
55: static PetscErrorCode PetscViewerGLVisSetPrecision_GLVis(PetscViewer viewer, PetscInt prec)
56: {
57: PetscErrorCode ierr;
58: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
61: PetscFree(socket->fmt);
62: if (prec > 0) {
63: PetscMalloc1(16,&socket->fmt);
64: PetscSNPrintf(socket->fmt,16," %%.%De",prec);
65: } else {
66: PetscStrallocpy(" %g",&socket->fmt);
67: }
68: return(0);
69: }
71: /*@
72: PetscViewerGLVisSetSnapId - Set the snapshot id. Only relevant when the viewer is of type PETSC_VIEWER_GLVIS_DUMP
74: Logically Collective on PetscViewer
76: Input Parameters:
77: + viewer - the PetscViewer
78: - id - the current snapshot id in a time-dependent simulation
80: Level: beginner
82: .seealso: PetscViewerGLVisOpen(), PetscViewerGLVisSetFields(), PetscViewerCreate(), PetscViewerSetType()
83: @*/
84: PetscErrorCode PetscViewerGLVisSetSnapId(PetscViewer viewer, PetscInt id)
85: {
91: PetscTryMethod(viewer,"PetscViewerGLVisSetSnapId_C",(PetscViewer,PetscInt),(viewer,id));
92: return(0);
93: }
95: static PetscErrorCode PetscViewerGLVisSetSnapId_GLVis(PetscViewer viewer, PetscInt id)
96: {
97: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
100: socket->snapid = id;
101: return(0);
102: }
104: /*
105: PetscViewerGLVisSetFields - Sets the required information to visualize different fields from within a vector.
107: Logically Collective on PetscViewer
109: Input Parameters:
110: + viewer - the PetscViewer
111: . nf - number of fields to be visualized
112: . namefield - optional name for each field
113: . fec_type - the type of finite element to be used to visualize the data (see FiniteElementCollection::Name() in MFEM)
114: . nlocal - array of local sizes for field vectors
115: . bs - array of block sizes for field vectors
116: . dim - array of topological dimension for field vectors
117: . g2lfields - User routine to compute the local field vectors to be visualized; PetscObject is used in place of Vec on the prototype
118: . ctx - User context to store the relevant data to apply g2lfields
119: - destroyctx - Destroy function for userctx
121: Notes: g2lfields is called on the Vec V to be visualized, in order to extract the relevant dofs to be put in Vfield[], as
122: .vb
123: g2lfields((PetscObject)V,nfields,(PetscObject*)Vfield[],ctx). Misses Fortran binding
124: .ve
126: Level: intermediate
128: .seealso: PetscViewerGLVisOpen(), PetscViewerCreate(), PetscViewerSetType()
129: */
130: PetscErrorCode PetscViewerGLVisSetFields(PetscViewer viewer, PetscInt nf, const char* namefield[], const char* fec_type[], PetscInt nlocal[], PetscInt bs[], PetscInt dim[], PetscErrorCode(*g2l)(PetscObject,PetscInt,PetscObject[],void*), void* ctx, PetscErrorCode(*destroyctx)(void*))
131: {
138: if (!fec_type) SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"You need to provide the FiniteElementCollection names for the fields");
143: PetscTryMethod(viewer,"PetscViewerGLVisSetFields_C",(PetscViewer,PetscInt,const char*[],const char*[],PetscInt[],PetscInt[],PetscInt[],PetscErrorCode(*)(PetscObject,PetscInt,PetscObject[],void*),void*,PetscErrorCode(*)(void*)),(viewer,nf,namefield,fec_type,nlocal,bs,dim,g2l,ctx,destroyctx));
144: return(0);
145: }
147: static PetscErrorCode PetscViewerGLVisSetFields_GLVis(PetscViewer viewer, PetscInt nfields, const char* namefield[], const char* fec_type[], PetscInt nlocal[], PetscInt bs[], PetscInt dim[], PetscErrorCode(*g2l)(PetscObject,PetscInt,PetscObject[],void*),void* ctx, PetscErrorCode(*destroyctx)(void*))
148: {
149: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
150: PetscInt i;
151: PetscErrorCode ierr;
154: if (socket->nwindow && socket->nwindow != nfields) SETERRQ2(PetscObjectComm((PetscObject)viewer),PETSC_ERR_USER,"Cannot set number of fields %D with number of windows %D",nfields,socket->nwindow);
155: if (!socket->nwindow) {
156: socket->nwindow = nfields;
158: PetscCalloc5(nfields,&socket->window,nfields,&socket->windowtitle,nfields,&socket->fec_type,3*nfields,&socket->locandbs,nfields,&socket->Ufield);
159: for (i=0;i<nfields;i++) {
160: if (!namefield || !namefield[i]) {
161: char name[16];
163: PetscSNPrintf(name,16,"Field%d",i);
164: PetscStrallocpy(name,&socket->windowtitle[i]);
165: } else {
166: PetscStrallocpy(namefield[i],&socket->windowtitle[i]);
167: }
168: PetscStrallocpy(fec_type[i],&socket->fec_type[i]);
169: socket->locandbs[3*i ] = nlocal[i];
170: socket->locandbs[3*i+1] = bs[i];
171: socket->locandbs[3*i+2] = dim[i];
172: }
173: }
174: /* number of fields are not allowed to vary */
175: if (nfields != socket->nwindow) SETERRQ2(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Cannot visualize %D fields using %D socket windows",nfields,socket->nwindow);
176: socket->g2lfield = g2l;
177: if (socket->destroyctx && socket->userctx) { (*socket->destroyctx)(socket->userctx); }
178: socket->userctx = ctx;
179: socket->destroyctx = destroyctx;
180: return(0);
181: }
183: static PetscErrorCode PetscViewerGLVisInfoDestroy_Private(void *ptr)
184: {
185: PetscViewerGLVisInfo info = (PetscViewerGLVisInfo)ptr;
186: PetscErrorCode ierr;
189: PetscFree(info->fmt);
190: PetscFree(info);
191: return(0);
192: }
194: /* we can decide to prevent specific processes from using the viewer */
195: static PetscErrorCode PetscViewerGLVisAttachInfo_Private(PetscViewer viewer, PetscViewer window)
196: {
197: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
198: PetscErrorCode ierr;
199: PetscContainer container;
200: PetscViewerGLVisInfo info;
203: PetscObjectQuery((PetscObject)window,"_glvis_info_container",(PetscObject*)&container);
204: if (!container) {
205: PetscNew(&info);
206: info->enabled = PETSC_TRUE;
207: info->init = PETSC_FALSE;
208: info->pause = socket->pause;
209: PetscContainerCreate(PetscObjectComm((PetscObject)window),&container);
210: PetscContainerSetPointer(container,(void*)info);
211: PetscContainerSetUserDestroy(container,PetscViewerGLVisInfoDestroy_Private);
212: PetscObjectCompose((PetscObject)window,"_glvis_info_container",(PetscObject)container);
213: PetscContainerDestroy(&container);
214: } else {
215: PetscContainerGetPointer(container,(void**)&info);
216: }
217: PetscFree(info->fmt);
218: PetscStrallocpy(socket->fmt,&info->fmt);
219: return(0);
220: }
222: static PetscErrorCode PetscViewerGLVisGetNewWindow_Private(PetscViewer viewer,PetscViewer *view)
223: {
224: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
225: PetscViewer window = NULL;
226: PetscBool ldis,dis;
227: PetscErrorCode ierr;
230: PetscViewerASCIISocketOpen(PETSC_COMM_SELF,socket->name,socket->port,&window);
231: /* if we could not estabilish a connection the first time,
232: we disable the socket viewer */
233: ldis = ierr ? PETSC_TRUE : PETSC_FALSE;
234: MPI_Allreduce(&ldis,&dis,1,MPIU_BOOL,MPI_LOR,PetscObjectComm((PetscObject)viewer));
235: if (dis) {
236: socket->status = PETSCVIEWERGLVIS_DISABLED;
237: PetscViewerDestroy(&window);
238: }
239: *view = window;
240: return(0);
241: }
243: PetscErrorCode PetscViewerGLVisPause_Private(PetscViewer viewer)
244: {
245: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
246: PetscErrorCode ierr;
249: if (socket->pause > 0) {
250: PetscSleep(socket->pause);
251: }
252: return(0);
253: }
255: /* DM specific support */
256: PetscErrorCode PetscViewerGLVisSetDM_Private(PetscViewer viewer, PetscObject dm)
257: {
258: PetscErrorCode ierr;
259: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
262: if (socket->dm && socket->dm != dm) SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Cannot change DM associated with the GLVis viewer");
263: if (!socket->dm) {
264: PetscErrorCode (*setupwithdm)(PetscObject,PetscViewer) = NULL;
266: PetscObjectQueryFunction(dm,"DMSetUpGLVisViewer_C",&setupwithdm);
267: if (setupwithdm) {
268: (*setupwithdm)(dm,viewer);
269: } else SETERRQ1(PetscObjectComm(dm),PETSC_ERR_SUP,"No support for DM type %s",dm->type_name);
270: PetscObjectReference(dm);
271: socket->dm = dm;
272: }
273: return(0);
274: }
276: PetscErrorCode PetscViewerGLVisGetDMWindow_Private(PetscViewer viewer,PetscViewer* view)
277: {
278: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
279: PetscErrorCode ierr;
282: if (!socket->meshwindow) {
283: if (socket->type == PETSC_VIEWER_GLVIS_SOCKET) {
284: PetscViewerGLVisGetNewWindow_Private(viewer,&socket->meshwindow);
285: } else {
286: PetscMPIInt rank;
287: char filename[PETSC_MAX_PATH_LEN];
289: MPI_Comm_rank(PetscObjectComm((PetscObject)viewer),&rank);
290: PetscSNPrintf(filename,PETSC_MAX_PATH_LEN,"%s-mesh.%06d",socket->name,rank);
291: PetscViewerASCIIOpen(PETSC_COMM_SELF,filename,&socket->meshwindow);
292: }
293: if (socket->meshwindow) {
294: PetscViewerPushFormat(socket->meshwindow,PETSC_VIEWER_ASCII_GLVIS);
295: }
296: }
297: if (socket->meshwindow) {
298: PetscViewerGLVisAttachInfo_Private(viewer,socket->meshwindow);
299: }
300: *view = socket->meshwindow;
301: return(0);
302: }
304: PetscErrorCode PetscViewerGLVisGetType_Private(PetscViewer viewer,PetscViewerGLVisType *type)
305: {
306: 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;
321: if (socket->type == PETSC_VIEWER_GLVIS_DUMP) {
322: socket->status = PETSCVIEWERGLVIS_DISCONNECTED;
323: } else if (socket->status == PETSCVIEWERGLVIS_DISCONNECTED && socket->nwindow) {
324: PetscInt i;
325: PetscBool lconn,conn;
328: for (i=0,lconn=PETSC_TRUE;i<socket->nwindow;i++)
329: if (!socket->window[i])
330: lconn = PETSC_FALSE;
332: MPI_Allreduce(&lconn,&conn,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)viewer));
333: if (conn) socket->status = PETSCVIEWERGLVIS_CONNECTED;
334: }
335: *sockstatus = socket->status;
336: return(0);
337: }
339: PetscErrorCode PetscViewerGLVisGetDM_Private(PetscViewer viewer, PetscObject* dm)
340: {
341: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
344: *dm = socket->dm;
345: return(0);
346: }
348: PetscErrorCode PetscViewerGLVisGetFields_Private(PetscViewer viewer, PetscInt* nfield, const char** names[], const char **fec[], PetscInt *locandbs[], PetscErrorCode(**g2lfield)(PetscObject,PetscInt,PetscObject[],void*), PetscObject *Ufield[], void **userctx)
349: {
350: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
353: if (nfield) *nfield = socket->nwindow;
354: if (names) *names = (const char**)socket->windowtitle;
355: if (fec) *fec = (const char**)socket->fec_type;
356: if (locandbs) *locandbs = socket->locandbs;
357: if (g2lfield) *g2lfield = socket->g2lfield;
358: if (Ufield) *Ufield = socket->Ufield;
359: if (userctx) *userctx = socket->userctx;
360: return(0);
361: }
363: /* accessor routines for the viewer windows:
364: PETSC_VIEWER_GLVIS_DUMP : it returns a new viewer every time
365: PETSC_VIEWER_GLVIS_SOCKET : it returns the socket, and creates it if not yet done.
366: */
367: PetscErrorCode PetscViewerGLVisGetWindow_Private(PetscViewer viewer,PetscInt wid,PetscViewer* view)
368: {
369: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
370: PetscViewerGLVisStatus status;
371: PetscErrorCode ierr;
376: if (wid < 0 || wid > socket->nwindow-1) SETERRQ2(PetscObjectComm((PetscObject)viewer),PETSC_ERR_USER,"Cannot get window id %D: allowed range [0,%D)",wid,socket->nwindow-1);
377: status = socket->status;
378: if (socket->type == PETSC_VIEWER_GLVIS_DUMP && socket->window[wid]) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Window %D is already in use",wid);
379: switch (status) {
380: case PETSCVIEWERGLVIS_DISCONNECTED:
381: if (socket->window[wid]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"This should not happen");
382: if (socket->type == PETSC_VIEWER_GLVIS_DUMP) {
383: PetscMPIInt rank;
384: char filename[PETSC_MAX_PATH_LEN];
386: MPI_Comm_rank(PetscObjectComm((PetscObject)viewer),&rank);
387: PetscSNPrintf(filename,PETSC_MAX_PATH_LEN,"%s-%s-%d.%06d",socket->name,socket->windowtitle[wid],socket->snapid,rank);
388: PetscViewerASCIIOpen(PETSC_COMM_SELF,filename,&socket->window[wid]);
389: } else {
390: PetscViewerGLVisGetNewWindow_Private(viewer,&socket->window[wid]);
391: }
392: if (socket->window[wid]) {
393: PetscViewerPushFormat(socket->window[wid],PETSC_VIEWER_ASCII_GLVIS);
394: }
395: *view = socket->window[wid];
396: break;
397: case PETSCVIEWERGLVIS_CONNECTED:
398: *view = socket->window[wid];
399: break;
400: case PETSCVIEWERGLVIS_DISABLED:
401: *view = NULL;
402: break;
403: default:
404: SETERRQ1(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Unhandled socket status %d\n",(int)status);
405: break;
406: }
407: if (*view) {
408: PetscViewerGLVisAttachInfo_Private(viewer,*view);
409: }
410: return(0);
411: }
413: /* Restore the window viewer
414: PETSC_VIEWER_GLVIS_DUMP : destroys the temporary created ASCII viewer used for dumping
415: PETSC_VIEWER_GLVIS_SOCKET: - if the returned window viewer is not NULL, just zeros the pointer.
416: - it the returned window viewer is NULL, assumes something went wrong
417: with the socket (i.e. SIGPIPE when a user closes the popup window)
418: and that the caller already handled it (see VecView_GLVis).
419: */
420: PetscErrorCode PetscViewerGLVisRestoreWindow_Private(PetscViewer viewer,PetscInt wid, PetscViewer* view)
421: {
422: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
423: PetscErrorCode ierr;
428: if (wid < 0 || wid > socket->nwindow-1) SETERRQ2(PetscObjectComm((PetscObject)viewer),PETSC_ERR_USER,"Cannot restore window id %D: allowed range [0,%D)",wid,socket->nwindow);
429: if (*view && *view != socket->window[wid]) SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_USER,"Viewer was not obtained from PetscViewerGLVisGetWindow");
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: PetscErrorCode ierr;
447: PetscViewerGLVisInfo info;
448: PetscContainer container;
451: PetscObjectQuery((PetscObject)viewer,"_glvis_info_container",(PetscObject*)&container);
452: if (!container) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Viewer was not obtained from PetscGLVisViewerGetNewWindow_Private");
453: PetscContainerGetPointer(container,(void**)&info);
454: if (info->init) {
455: if (info->pause < 0) {
456: PetscViewerASCIIPrintf(viewer,"pause\n"); /* pause */
457: }
458: return(0);
459: }
460: PetscViewerASCIIPrintf(viewer,"window_size 800 800\n");
461: if (name) {
462: PetscViewerASCIIPrintf(viewer,"window_title '%s'\n",name);
463: }
464: if (mesh) {
465: switch (dim) {
466: case 1:
467: break;
468: case 2:
469: PetscViewerASCIIPrintf(viewer,"keys cmeeppppp\n"); /* show colorbar, mesh and ranks */
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 RRj\n"); /* set to 1D (side view) and turn off perspective */
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: if (info->pause == 0) {
488: PetscViewerASCIIPrintf(viewer,"pause\n"); /* pause */
489: }
490: }
491: info->init = PETSC_TRUE;
492: return(0);
493: }
495: static PetscErrorCode PetscViewerDestroy_GLVis(PetscViewer viewer)
496: {
497: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
498: PetscInt i;
499: PetscErrorCode ierr;
502: for (i=0;i<socket->nwindow;i++) {
503: PetscViewerDestroy(&socket->window[i]);
504: PetscFree(socket->windowtitle[i]);
505: PetscFree(socket->fec_type[i]);
506: PetscObjectDestroy(&socket->Ufield[i]);
507: }
508: PetscFree(socket->name);
509: PetscFree(socket->fmt);
510: PetscFree5(socket->window,socket->windowtitle,socket->fec_type,socket->locandbs,socket->Ufield);
511: PetscViewerDestroy(&socket->meshwindow);
512: PetscObjectDestroy(&socket->dm);
513: if (socket->destroyctx && socket->userctx) { (*socket->destroyctx)(socket->userctx); }
515: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetPrecision_C",NULL);
516: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetSnapId_C",NULL);
517: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetFields_C",NULL);
518: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileSetName_C",NULL);
519: PetscFree(socket);
520: viewer->data = NULL;
521: return(0);
522: }
524: static PetscErrorCode PetscViewerSetFromOptions_GLVis(PetscOptionItems *PetscOptionsObject,PetscViewer v)
525: {
526: PetscErrorCode ierr;
527: PetscViewerGLVis socket = (PetscViewerGLVis)v->data;
530: PetscOptionsHead(PetscOptionsObject,"GLVis PetscViewer Options");
531: PetscOptionsReal("-viewer_glvis_pause","-1 to pause after each visualization, otherwise sleeps for given seconds",NULL,socket->pause,&socket->pause,NULL);
532: PetscOptionsTail();
533: return(0);
534: }
536: static PetscErrorCode PetscViewerSetFileName_GLVis(PetscViewer viewer, const char name[])
537: {
538: char *sport;
539: PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
540: PetscErrorCode ierr;
543: socket->type = PETSC_VIEWER_GLVIS_DUMP;
544: /* we accept localhost^port */
545: PetscFree(socket->name);
546: PetscStrallocpy(name,&socket->name);
547: PetscStrchr(socket->name,'^',&sport);
548: if (sport) {
549: PetscInt port = 19916;
550: size_t len;
552: *sport++ = 0;
553: PetscStrlen(sport,&len);
554: if (len) PetscOptionsStringToInt(sport,&port);
555: if (!ierr) {
556: socket->port = (port != PETSC_DECIDE && port != PETSC_DEFAULT) ? port : 19916;
557: } else {
558: socket->port = 19916;
559: }
560: socket->type = PETSC_VIEWER_GLVIS_SOCKET;
561: }
562: return(0);
563: }
565: /*
566: PetscViewerGLVisOpen - Opens a GLVis type viewer
568: Collective on comm
570: Input Parameters:
571: + comm - the MPI communicator
572: . type - the viewer type: PETSC_VIEWER_GLVIS_SOCKET for real-time visualization or PETSC_VIEWER_GLVIS_DUMP for dumping to disk
573: . name - either the hostname where the GLVis server is running or the base filename for dumping the data for subsequent visualizations
574: - port - socket port where the GLVis server is listening. Not referenced when type is PETSC_VIEWER_GLVIS_DUMP
576: Output Parameters:
577: - viewer - the PetscViewer object
579: Notes: misses Fortran binding
581: Level: beginner
583: .seealso: PetscViewerCreate(), PetscViewerSetType(), PetscViewerGLVisType
584: */
585: PETSC_EXTERN PetscErrorCode PetscViewerGLVisOpen(MPI_Comm comm, PetscViewerGLVisType type, const char* name, PetscInt port, PetscViewer* viewer)
586: {
587: PetscViewerGLVis socket;
588: PetscErrorCode ierr;
591: PetscViewerCreate(comm,viewer);
592: PetscViewerSetType(*viewer,PETSCVIEWERGLVIS);
594: socket = (PetscViewerGLVis)((*viewer)->data);
595: PetscFree(socket->name);
596: PetscStrallocpy(name,&socket->name);
597: socket->type = type;
598: socket->port = port;
600: PetscViewerSetFromOptions(*viewer);
601: return(0);
602: }
604: /*
605: PETSC_VIEWER_GLVIS_ - Creates an GLVIS PetscViewer shared by all processors in a communicator.
607: Collective on MPI_Comm
609: Input Parameter:
610: . comm - the MPI communicator to share the GLVIS PetscViewer
612: Level: intermediate
614: Notes: misses Fortran bindings
616: Environmental variables:
617: + PETSC_VIEWER_GLVIS_FILENAME : output filename (if specified dump to disk, and takes precedence on PETSC_VIEWER_GLVIS_HOSTNAME)
618: . PETSC_VIEWER_GLVIS_HOSTNAME : machine where the GLVis server is listening (defaults to localhost)
619: - PETSC_VIEWER_GLVIS_PORT : port opened by the GLVis server (defaults to 19916)
621: Notes:
622: Unlike almost all other PETSc routines, PETSC_VIEWER_GLVIS_ does not return
623: an error code. The GLVIS PetscViewer is usually used in the form
624: $ XXXView(XXX object, PETSC_VIEWER_GLVIS_(comm));
626: .seealso: PetscViewerGLVISOpen(), PetscViewerGLVisType, PetscViewerCreate(), PetscViewerDestroy()
627: */
628: PETSC_EXTERN PetscViewer PETSC_VIEWER_GLVIS_(MPI_Comm comm)
629: {
630: PetscErrorCode ierr;
631: PetscBool flg;
632: PetscViewer viewer;
633: PetscViewerGLVisType type;
634: char fname[PETSC_MAX_PATH_LEN],sport[16];
635: PetscInt port = 19916; /* default for GLVis */
638: PetscOptionsGetenv(comm,"PETSC_VIEWER_GLVIS_FILENAME",fname,PETSC_MAX_PATH_LEN,&flg);
639: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
640: if (!flg) {
641: type = PETSC_VIEWER_GLVIS_SOCKET;
642: PetscOptionsGetenv(comm,"PETSC_VIEWER_GLVIS_HOSTNAME",fname,PETSC_MAX_PATH_LEN,&flg);
643: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
644: if (!flg) {
645: PetscStrcpy(fname,"localhost");
646: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
647: }
648: PetscOptionsGetenv(comm,"PETSC_VIEWER_GLVIS_PORT",sport,16,&flg);
649: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
650: if (flg) {
651: PetscOptionsStringToInt(sport,&port);
652: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
653: }
654: } else {
655: type = PETSC_VIEWER_GLVIS_DUMP;
656: }
657: PetscViewerGLVisOpen(comm,type,fname,port,&viewer);
658: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
659: PetscObjectRegisterDestroy((PetscObject)viewer);
660: if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
661: PetscFunctionReturn(viewer);
662: }
664: PETSC_EXTERN PetscErrorCode PetscViewerCreate_GLVis(PetscViewer viewer)
665: {
666: PetscViewerGLVis socket;
667: PetscErrorCode ierr;
670: PetscNewLog(viewer,&socket);
672: /* defaults to socket viewer */
673: PetscStrallocpy("localhost",&socket->name);
674: socket->port = 19916; /* GLVis default listening port */
675: socket->type = PETSC_VIEWER_GLVIS_SOCKET;
676: socket->pause = 0; /* just pause the first time */
678: /* defaults to full precision */
679: PetscStrallocpy(" %g",&socket->fmt);
681: viewer->data = (void*)socket;
682: viewer->ops->destroy = PetscViewerDestroy_GLVis;
683: viewer->ops->setfromoptions = PetscViewerSetFromOptions_GLVis;
685: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetPrecision_C",PetscViewerGLVisSetPrecision_GLVis);
686: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetSnapId_C",PetscViewerGLVisSetSnapId_GLVis);
687: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetFields_C",PetscViewerGLVisSetFields_GLVis);
688: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileSetName_C",PetscViewerSetFileName_GLVis);
689: return(0);
690: }
692: /* this is a private implementation of a SOCKET with ASCII data format
693: GLVis does not currently handle binary socket streams */
694: #if defined(PETSC_HAVE_UNISTD_H)
695: #include <unistd.h>
696: #endif
698: #if !defined(PETSC_HAVE_WINDOWS_H)
699: static PetscErrorCode (*PetscViewerDestroy_ASCII)(PetscViewer);
701: static PetscErrorCode PetscViewerDestroy_ASCII_Socket(PetscViewer viewer)
702: {
703: FILE *stream;
704: PetscErrorCode 0;
706: PetscViewerASCIIGetPointer(viewer,&stream);
707: if (stream) {
708: fclose(stream);
709: if (ierr) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SYS,"fclose() failed on stream");
710: }
711: PetscViewerDestroy_ASCII(viewer);
712: return(0);
713: }
714: #endif
716: static PetscErrorCode PetscViewerASCIISocketOpen(MPI_Comm comm,const char* hostname,PetscInt port,PetscViewer* viewer)
717: {
718: #if defined(PETSC_HAVE_WINDOWS_H)
720: SETERRQ(comm,PETSC_ERR_SUP,"Not implemented for Windows");
721: #else
722: FILE *stream = NULL;
723: int fd;
729: #if defined(PETSC_USE_SOCKET_VIEWER)
730: PetscOpenSocket(hostname,port,&fd);
731: #else
732: SETERRQ(comm,PETSC_ERR_SUP,"Missing Socket viewer");
733: #endif
734: if (ierr) {
735: PetscInt sierr;
736: char err[1024];
738: PetscSNPrintf(err,1024,"Cannot connect to socket on %s:%D. Socket visualization is disabled\n",hostname,port);
739: PetscInfo(NULL,err);
740: *viewer = NULL;
741: PetscFunctionReturn(sierr);
742: } else {
743: char msg[1024];
745: PetscSNPrintf(msg,1024,"Successfully connect to socket on %s:%D. Socket visualization is enabled\n",hostname,port);
746: PetscInfo(NULL,msg);
747: }
748: stream = fdopen(fd,"w"); /* Not possible on Windows */
749: if (!stream) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_SYS,"Cannot open stream from socket %s:%d",hostname,port);
750: PetscViewerASCIIOpenWithFILE(PETSC_COMM_SELF,stream,viewer);
751: PetscViewerDestroy_ASCII = (*viewer)->ops->destroy;
752: (*viewer)->ops->destroy = PetscViewerDestroy_ASCII_Socket;
753: #endif
754: return(0);
755: }