Actual source code: glvis.c

petsc-3.8.0 2017-09-26
Report Typos and Errors
  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: };
 29: typedef struct _n_PetscViewerGLVis *PetscViewerGLVis;

 31: /*@
 32:      PetscViewerGLVisSetSnapId - Set the snapshot id. Only relevant when the viewer is of type PETSC_VIEWER_GLVIS_DUMP

 34:   Logically Collective on PetscViewer

 36:   Input Parameters:
 37: +  viewer - the PetscViewer
 38: -  id     - the current snapshot id in a time-dependent simulation

 40:   Level: beginner

 42: .seealso: PetscViewerGLVisOpen(), PetscViewerCreate(), PetscViewerSetType()
 43: @*/
 44: PetscErrorCode PetscViewerGLVisSetSnapId(PetscViewer viewer, PetscInt id)
 45: {

 51:   PetscTryMethod(viewer,"PetscViewerGLVisSetSnapId_C",(PetscViewer,PetscInt),(viewer,id));
 52:   return(0);
 53: }

 55: static PetscErrorCode PetscViewerGLVisSetSnapId_GLVis(PetscViewer viewer, PetscInt id)
 56: {
 57:   PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;

 60:   socket->snapid = id;;
 61:   return(0);
 62: }

 64: /*
 65:      PetscViewerGLVisSetFields - Sets the required information to visualize different fields from within a vector.

 67:   Logically Collective on PetscViewer

 69:   Input Parameters:
 70: +  viewer     - the PetscViewer
 71: .  nf         - number of fields to be visualized
 72: .  namefield  - optional name for each field
 73: .  fec_type   - the type of finite element to be used to visualize the data (see FiniteElementCollection::Name() in MFEM)
 74: .  nlocal     - array of local sizes for field vectors
 75: .  bs         - array of block sizes for field vectors
 76: .  dim        - array of topological dimension for field vectors
 77: .  g2lfields  - User routine to compute the local field vectors to be visualized; PetscObject is used in place of Vec on the prototype
 78: .  ctx        - User context to store the relevant data to apply g2lfields
 79: -  destroyctx - Destroy function for userctx

 81:   Notes: g2lfields is called on the Vec V to be visualized, in order to extract the relevant dofs to be put in Vfield[], as
 82: .vb
 83:   g2lfields((PetscObject)V,nfields,(PetscObject*)Vfield[],ctx). Misses Fortran binding
 84: .ve

 86:   Level: intermediate

 88: .seealso: PetscViewerGLVisOpen(), PetscViewerCreate(), PetscViewerSetType()
 89: */
 90: 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*))
 91: {

 98:   if (!fec_type) SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"You need to provide the FiniteElementCollection names for the fields");
103:   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));
104:   return(0);
105: }

107: 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*))
108: {
109:   PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
110:   PetscInt         i;
111:   PetscErrorCode   ierr;

114:   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);
115:   if (!socket->nwindow) {
116:     socket->nwindow = nfields;

118:     PetscCalloc5(nfields,&socket->window,nfields,&socket->windowtitle,nfields,&socket->fec_type,3*nfields,&socket->locandbs,nfields,&socket->Ufield);
119:     for (i=0;i<nfields;i++) {
120:       if (!namefield || !namefield[i]) {
121:         char name[16];

123:         PetscSNPrintf(name,16,"Field%d",i);
124:         PetscStrallocpy(name,&socket->windowtitle[i]);
125:       } else {
126:         PetscStrallocpy(namefield[i],&socket->windowtitle[i]);
127:       }
128:       PetscStrallocpy(fec_type[i],&socket->fec_type[i]);
129:       socket->locandbs[3*i  ] = nlocal[i];
130:       socket->locandbs[3*i+1] = bs[i];
131:       socket->locandbs[3*i+2] = dim[i];
132:     }
133:   }
134:   /* number of fields are not allowed to vary */
135:   if (nfields != socket->nwindow) SETERRQ2(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Cannot visualize %D fields using %D socket windows",nfields,socket->nwindow);
136:   socket->g2lfield = g2l;
137:   if (socket->destroyctx && socket->userctx) { (*socket->destroyctx)(socket->userctx); }
138:   socket->userctx = ctx;
139:   socket->destroyctx = destroyctx;
140:   return(0);
141: }

143: static PetscErrorCode PetscViewerGLVisInfoDestroy_Private(void *ptr)
144: {
145:   PetscViewerGLVisInfo info = (PetscViewerGLVisInfo)ptr;
146:   PetscErrorCode       ierr;

149:   PetscFree(info);
150:   return(0);
151: }

153: /* we can decide to prevent specific processes from using the viewer */
154: static PetscErrorCode PetscViewerGLVisAttachInfo_Private(PetscViewer viewer, PetscViewer window)
155: {
156:   PetscViewerGLVis     socket = (PetscViewerGLVis)viewer->data;
157:   PetscErrorCode       ierr;
158:   PetscContainer       container;
159:   PetscViewerGLVisInfo info;

162:   PetscNew(&info);
163:   info->enabled = PETSC_TRUE;
164:   info->init    = PETSC_FALSE;
165:   info->pause   = socket->pause;
166:   PetscContainerCreate(PetscObjectComm((PetscObject)window),&container);
167:   PetscContainerSetPointer(container,(void*)info);
168:   PetscContainerSetUserDestroy(container,PetscViewerGLVisInfoDestroy_Private);
169:   PetscObjectCompose((PetscObject)window,"_glvis_info_container",(PetscObject)container);
170:   PetscContainerDestroy(&container);
171:   return(0);
172: }

174: static PetscErrorCode PetscViewerGLVisGetNewWindow_Private(PetscViewer viewer,PetscViewer *view)
175: {
176:   PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
177:   PetscViewer      window = NULL;
178:   PetscBool        ldis,dis;
179:   PetscErrorCode   ierr;

182:   PetscViewerASCIISocketOpen(PETSC_COMM_SELF,socket->name,socket->port,&window);
183:   /* if we could not estabilish a connection the first time,
184:      we disable the socket viewer */
185:   ldis = ierr ? PETSC_TRUE : PETSC_FALSE;
186:   MPI_Allreduce(&ldis,&dis,1,MPIU_BOOL,MPI_LOR,PetscObjectComm((PetscObject)viewer));
187:   if (dis) {
188:     socket->status = PETSCVIEWERGLVIS_DISABLED;
189:     PetscViewerDestroy(&window);
190:   }
191:   *view = window;
192:   return(0);
193: }

195: PetscErrorCode PetscViewerGLVisPause_Private(PetscViewer viewer)
196: {
197:   PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
198:   PetscErrorCode   ierr;

201:   if (socket->pause > 0) {
202:     PetscSleep(socket->pause);
203:   }
204:   return(0);
205: }

207: /* DM specific support */
208: PetscErrorCode PetscViewerGLVisSetDM_Private(PetscViewer viewer, PetscObject dm)
209: {
210:   PetscErrorCode   ierr;
211:   PetscViewerGLVis socket  = (PetscViewerGLVis)viewer->data;

214:   if (socket->dm && socket->dm != dm) SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Cannot change DM associated with the GLVis viewer");
215:   if (!socket->dm) {
216:     PetscErrorCode (*setupwithdm)(PetscObject,PetscViewer) = NULL;

218:     PetscObjectQueryFunction(dm,"DMSetUpGLVisViewer_C",&setupwithdm);
219:     if (setupwithdm) {
220:       (*setupwithdm)(dm,viewer);
221:     } else SETERRQ1(PetscObjectComm(dm),PETSC_ERR_SUP,"No support for DM type %s",dm->type_name);
222:     PetscObjectReference(dm);
223:     socket->dm = dm;
224:   }
225:   return(0);
226: }

228: PetscErrorCode PetscViewerGLVisGetDMWindow_Private(PetscViewer viewer,PetscViewer* view)
229: {
230:   PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
231:   PetscErrorCode   ierr;

234:   if (!socket->meshwindow) {
235:     if (socket->type == PETSC_VIEWER_GLVIS_SOCKET) {
236:       PetscViewerGLVisGetNewWindow_Private(viewer,&socket->meshwindow);
237:     } else {
238:       PetscMPIInt rank;
239:       char        filename[PETSC_MAX_PATH_LEN];

241:       MPI_Comm_rank(PetscObjectComm((PetscObject)viewer),&rank);
242:       PetscSNPrintf(filename,PETSC_MAX_PATH_LEN,"%s-mesh.%06d",socket->name,rank);
243:       PetscViewerASCIIOpen(PETSC_COMM_SELF,filename,&socket->meshwindow);
244:     }
245:     if (socket->meshwindow) {
246:       PetscViewerGLVisAttachInfo_Private(viewer,socket->meshwindow);
247:       PetscViewerPushFormat(socket->meshwindow,PETSC_VIEWER_ASCII_GLVIS);
248:     }
249:   }
250:   *view = socket->meshwindow;
251:   return(0);
252: }

254: PetscErrorCode PetscViewerGLVisGetType_Private(PetscViewer viewer,PetscViewerGLVisType *type)
255: {
256:   PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;

260:   *type = socket->type;
261:   return(0);
262: }

264: /* 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 */
265: PetscErrorCode PetscViewerGLVisGetStatus_Private(PetscViewer viewer, PetscViewerGLVisStatus *sockstatus)
266: {
267:   PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;

271:   if (socket->type == PETSC_VIEWER_GLVIS_DUMP) {
272:     socket->status = PETSCVIEWERGLVIS_DISCONNECTED;
273:   } else if (socket->status == PETSCVIEWERGLVIS_DISCONNECTED && socket->nwindow) {
274:     PetscInt       i;
275:     PetscBool      lconn,conn;

278:     for (i=0,lconn=PETSC_TRUE;i<socket->nwindow;i++)
279:       if (!socket->window[i])
280:         lconn = PETSC_FALSE;

282:     MPI_Allreduce(&lconn,&conn,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)viewer));
283:     if (conn) socket->status = PETSCVIEWERGLVIS_CONNECTED;
284:   }
285:   *sockstatus = socket->status;
286:   return(0);
287: }

289: PetscErrorCode PetscViewerGLVisGetDM_Private(PetscViewer viewer, PetscObject* dm)
290: {
291:   PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;

294:   *dm = socket->dm;
295:   return(0);
296: }

298: 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)
299: {
300:   PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;

303:   if (nfield)   *nfield   = socket->nwindow;
304:   if (names)    *names    = (const char**)socket->windowtitle;
305:   if (fec)      *fec      = (const char**)socket->fec_type;
306:   if (locandbs) *locandbs = socket->locandbs;
307:   if (g2lfield) *g2lfield = socket->g2lfield;
308:   if (Ufield)   *Ufield   = socket->Ufield;
309:   if (userctx)  *userctx  = socket->userctx;
310:   return(0);
311: }

313: /* accessor routines for the viewer windows:
314:    PETSC_VIEWER_GLVIS_DUMP   : it returns a new viewer every time
315:    PETSC_VIEWER_GLVIS_SOCKET : it returns the socket, and creates it if not yet done.
316: */
317: PetscErrorCode PetscViewerGLVisGetWindow_Private(PetscViewer viewer,PetscInt wid,PetscViewer* view)
318: {
319:   PetscViewerGLVis       socket = (PetscViewerGLVis)viewer->data;
320:   PetscViewerGLVisStatus status;
321:   PetscErrorCode         ierr;

326:   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);
327:   status = socket->status;
328:   if (socket->type == PETSC_VIEWER_GLVIS_DUMP && socket->window[wid]) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Window %D is already in use",wid);
329:   switch (status) {
330:     case PETSCVIEWERGLVIS_DISCONNECTED:
331:       if (socket->window[wid]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"This should not happen");
332:       if (socket->type == PETSC_VIEWER_GLVIS_DUMP) {
333:         PetscMPIInt rank;
334:         char        filename[PETSC_MAX_PATH_LEN];

336:         MPI_Comm_rank(PetscObjectComm((PetscObject)viewer),&rank);
337:         PetscSNPrintf(filename,PETSC_MAX_PATH_LEN,"%s-%s-%d.%06d",socket->name,socket->windowtitle[wid],socket->snapid,rank);
338:         PetscViewerASCIIOpen(PETSC_COMM_SELF,filename,&socket->window[wid]);
339:       } else {
340:         PetscViewerGLVisGetNewWindow_Private(viewer,&socket->window[wid]);
341:       }
342:       if (socket->window[wid]) {
343:         PetscViewerGLVisAttachInfo_Private(viewer,socket->window[wid]);
344:         PetscViewerPushFormat(socket->window[wid],PETSC_VIEWER_ASCII_GLVIS);
345:       }
346:       *view = socket->window[wid];
347:       break;
348:     case PETSCVIEWERGLVIS_CONNECTED:
349:       *view = socket->window[wid];
350:       break;
351:     case PETSCVIEWERGLVIS_DISABLED:
352:       *view = NULL;
353:       break;
354:     default:
355:       SETERRQ1(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Unhandled socket status %d\n",(int)status);
356:       break;
357:   }
358:   return(0);
359: }

361: /* Restore the window viewer
362:    PETSC_VIEWER_GLVIS_DUMP  : destroys the temporary created ASCII viewer used for dumping
363:    PETSC_VIEWER_GLVIS_SOCKET: - if the returned window viewer is not NULL, just zeros the pointer.
364:                  - it the returned window viewer is NULL, assumes something went wrong
365:                    with the socket (i.e. SIGPIPE when a user closes the popup window)
366:                    and that the caller already handled it (see VecView_GLVis).
367: */
368: PetscErrorCode PetscViewerGLVisRestoreWindow_Private(PetscViewer viewer,PetscInt wid, PetscViewer* view)
369: {
370:   PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
371:   PetscErrorCode   ierr;

376:   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);
377:   if (*view && *view != socket->window[wid]) SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_USER,"Viewer was not obtained from PetscViewerGLVisGetWindow");
378:   if (*view) {
379:     PetscViewerFlush(*view);
380:     PetscBarrier((PetscObject)viewer);
381:   }
382:   if (socket->type == PETSC_VIEWER_GLVIS_DUMP) { /* destroy the viewer, as it is associated with a single time step */
383:     PetscViewerDestroy(&socket->window[wid]);
384:   } else if (!*view) { /* something went wrong (SIGPIPE) so we just zero the private pointer */
385:     socket->window[wid] = NULL;
386:   }
387:   *view = NULL;
388:   return(0);
389: }

391: /* default window appearance in the PETSC_VIEWER_GLVIS_SOCKET case */
392: PetscErrorCode PetscViewerGLVisInitWindow_Private(PetscViewer viewer, PetscBool mesh, PetscInt dim, const char *name)
393: {
394:   PetscErrorCode       ierr;
395:   PetscViewerGLVisInfo info;
396:   PetscContainer       container;

399:   PetscObjectQuery((PetscObject)viewer,"_glvis_info_container",(PetscObject*)&container);
400:   if (!container) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Viewer was not obtained from PetscGLVisViewerGetNewWindow_Private");
401:   PetscContainerGetPointer(container,(void**)&info);
402:   if (info->init) {
403:     if (info->pause < 0) {
404:       PetscViewerASCIIPrintf(viewer,"pause\n"); /* pause */
405:     }
406:     return(0);
407:   }
408:   PetscViewerASCIIPrintf(viewer,"window_size 800 800\n");
409:   if (name) {
410:     PetscViewerASCIIPrintf(viewer,"window_title '%s'\n",name);
411:   }
412:   if (mesh) {
413:     switch (dim) {
414:     case 1:
415:       break;
416:     case 2:
417:       PetscViewerASCIIPrintf(viewer,"keys cmeeppppp\n"); /* show colorbar, mesh and ranks */
418:       break;
419:     case 3: /* TODO: decide default view in 3D */
420:       break;
421:     }
422:   } else {
423:     PetscViewerASCIIPrintf(viewer,"keys cm\n"); /* show colorbar and mesh */
424:     switch (dim) {
425:     case 1:
426:       PetscViewerASCIIPrintf(viewer,"keys RRj\n"); /* set to 1D (side view) and turn off perspective */
427:       break;
428:     case 2:
429:       PetscViewerASCIIPrintf(viewer,"keys Rjl\n"); /* set to 2D (top view), turn off perspective and light */
430:       break;
431:     case 3:
432:       break;
433:     }
434:     PetscViewerASCIIPrintf(viewer,"autoscale value\n"); /* update value-range; keep mesh-extents fixed */
435:     if (info->pause == 0) {
436:       PetscViewerASCIIPrintf(viewer,"pause\n"); /* pause */
437:     }
438:   }
439:   info->init = PETSC_TRUE;
440:   return(0);
441: }

443: static PetscErrorCode PetscViewerDestroy_GLVis(PetscViewer viewer)
444: {
445:   PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
446:   PetscInt         i;
447:   PetscErrorCode   ierr;

450:   for (i=0;i<socket->nwindow;i++) {
451:     PetscViewerDestroy(&socket->window[i]);
452:     PetscFree(socket->windowtitle[i]);
453:     PetscFree(socket->fec_type[i]);
454:     PetscObjectDestroy(&socket->Ufield[i]);
455:   }
456:   PetscFree(socket->name);
457:   PetscFree5(socket->window,socket->windowtitle,socket->fec_type,socket->locandbs,socket->Ufield);
458:   PetscViewerDestroy(&socket->meshwindow);
459:   PetscObjectDestroy(&socket->dm);
460:   if (socket->destroyctx && socket->userctx) { (*socket->destroyctx)(socket->userctx); }

462:   PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetSnapId_C",NULL);
463:   PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetFields_C",NULL);
464:   PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileSetName_C",NULL);
465:   PetscFree(socket);
466:   viewer->data = NULL;
467:   return(0);
468: }

470: static PetscErrorCode PetscViewerSetFromOptions_GLVis(PetscOptionItems *PetscOptionsObject,PetscViewer v)
471: {
472:   PetscErrorCode   ierr;
473:   PetscViewerGLVis socket = (PetscViewerGLVis)v->data;

476:   PetscOptionsHead(PetscOptionsObject,"GLVis PetscViewer Options");
477:   PetscOptionsReal("-viewer_glvis_pause","-1 to pause after each visualization, otherwise sleeps for given seconds",NULL,socket->pause,&socket->pause,NULL);
478:   PetscOptionsTail();
479:   return(0);
480: }

482: static PetscErrorCode PetscViewerSetFileName_GLVis(PetscViewer viewer, const char name[])
483: {
484:   char             *sport;
485:   PetscViewerGLVis socket = (PetscViewerGLVis)viewer->data;
486:   PetscErrorCode   ierr;

489:   socket->type = PETSC_VIEWER_GLVIS_DUMP;
490:   /* we accept localhost^port */
491:   PetscFree(socket->name);
492:   PetscStrallocpy(name,&socket->name);
493:   PetscStrchr(socket->name,'^',&sport);
494:   if (sport) {
495:     PetscInt port = 19916;
496:     size_t   len;

498:     *sport++ = 0;
499:     PetscStrlen(sport,&len);
500:     if (len) PetscOptionsStringToInt(sport,&port);
501:     if (!ierr) {
502:       socket->port = (port != PETSC_DECIDE && port != PETSC_DEFAULT) ? port : 19916;
503:     } else {
504:       socket->port = 19916;
505:     }
506:     socket->type = PETSC_VIEWER_GLVIS_SOCKET;
507:   }
508:   return(0);
509: }

511: /*
512:      PetscViewerGLVisOpen - Opens a GLVis type viewer

514:   Collective on comm

516:   Input Parameters:
517: +  comm      - the MPI communicator
518: .  type      - the viewer type: PETSC_VIEWER_GLVIS_SOCKET for real-time visualization or PETSC_VIEWER_GLVIS_DUMP for dumping to disk
519: .  name      - either the hostname where the GLVis server is running or the base filename for dumping the data for subsequent visualizations
520: -  port      - socket port where the GLVis server is listening. Not referenced when type is PETSC_VIEWER_GLVIS_DUMP

522:   Output Parameters:
523: -  viewer    - the PetscViewer object

525:   Notes: misses Fortran binding

527:   Level: beginner

529: .seealso: PetscViewerCreate(), PetscViewerSetType(), PetscViewerGLVisType
530: */
531: PETSC_EXTERN PetscErrorCode PetscViewerGLVisOpen(MPI_Comm comm, PetscViewerGLVisType type, const char* name, PetscInt port, PetscViewer* viewer)
532: {
533:   PetscViewerGLVis socket;
534:   PetscErrorCode   ierr;

537:   PetscViewerCreate(comm,viewer);
538:   PetscViewerSetType(*viewer,PETSCVIEWERGLVIS);

540:   socket       = (PetscViewerGLVis)((*viewer)->data);
541:   PetscFree(socket->name);
542:   PetscStrallocpy(name,&socket->name);
543:   socket->type = type;
544:   socket->port = port;

546:   PetscViewerSetFromOptions(*viewer);
547:   return(0);
548: }

550: /*
551:   PETSC_VIEWER_GLVIS_ - Creates an GLVIS PetscViewer shared by all processors in a communicator.

553:   Collective on MPI_Comm

555:   Input Parameter:
556: . comm - the MPI communicator to share the GLVIS PetscViewer

558:   Level: intermediate

560:   Notes: misses Fortran bindings

562:   Environmental variables:
563: + PETSC_VIEWER_GLVIS_FILENAME : output filename (if specified dump to disk, and takes precedence on PETSC_VIEWER_GLVIS_HOSTNAME)
564: . PETSC_VIEWER_GLVIS_HOSTNAME : machine where the GLVis server is listening (defaults to localhost)
565: - PETSC_VIEWER_GLVIS_PORT     : port opened by the GLVis server (defaults to 19916)

567:   Notes:
568:   Unlike almost all other PETSc routines, PETSC_VIEWER_GLVIS_ does not return
569:   an error code.  The GLVIS PetscViewer is usually used in the form
570: $       XXXView(XXX object, PETSC_VIEWER_GLVIS_(comm));

572: .seealso: PetscViewerGLVISOpen(), PetscViewerGLVisType, PetscViewerCreate(), PetscViewerDestroy()
573: */
574: PETSC_EXTERN PetscViewer PETSC_VIEWER_GLVIS_(MPI_Comm comm)
575: {
576:   PetscErrorCode       ierr;
577:   PetscBool            flg;
578:   PetscViewer          viewer;
579:   PetscViewerGLVisType type;
580:   char                 fname[PETSC_MAX_PATH_LEN],sport[16];
581:   PetscInt             port = 19916; /* default for GLVis */

584:   PetscOptionsGetenv(comm,"PETSC_VIEWER_GLVIS_FILENAME",fname,PETSC_MAX_PATH_LEN,&flg);
585:   if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
586:   if (!flg) {
587:     type = PETSC_VIEWER_GLVIS_SOCKET;
588:     PetscOptionsGetenv(comm,"PETSC_VIEWER_GLVIS_HOSTNAME",fname,PETSC_MAX_PATH_LEN,&flg);
589:     if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
590:     if (!flg) {
591:       PetscStrcpy(fname,"localhost");
592:       if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
593:     }
594:     PetscOptionsGetenv(comm,"PETSC_VIEWER_GLVIS_PORT",sport,16,&flg);
595:     if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
596:     if (flg) {
597:       PetscOptionsStringToInt(sport,&port);
598:       if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
599:     }
600:   } else {
601:     type = PETSC_VIEWER_GLVIS_DUMP;
602:   }
603:   PetscViewerGLVisOpen(comm,type,fname,port,&viewer);
604:   if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
605:   PetscObjectRegisterDestroy((PetscObject)viewer);
606:   if (ierr) {PetscError(PETSC_COMM_SELF,__LINE__,"PETSC_VIEWER_GLVIS_",__FILE__,PETSC_ERR_PLIB,PETSC_ERROR_INITIAL," ");return(0);}
607:   PetscFunctionReturn(viewer);
608: }

610: PETSC_EXTERN PetscErrorCode PetscViewerCreate_GLVis(PetscViewer viewer)
611: {
612:   PetscViewerGLVis socket;
613:   PetscErrorCode   ierr;

616:   PetscNewLog(viewer,&socket);

618:   /* defaults to socket viewer */
619:   PetscStrallocpy("localhost",&socket->name);
620:   socket->port  = 19916; /* GLVis default listening port */
621:   socket->type  = PETSC_VIEWER_GLVIS_SOCKET;
622:   socket->pause = 0; /* just pause the first time */

624:   viewer->data                = (void*)socket;
625:   viewer->ops->destroy        = PetscViewerDestroy_GLVis;
626:   viewer->ops->setfromoptions = PetscViewerSetFromOptions_GLVis;

628:   PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetSnapId_C",PetscViewerGLVisSetSnapId_GLVis);
629:   PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerGLVisSetFields_C",PetscViewerGLVisSetFields_GLVis);
630:   PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileSetName_C",PetscViewerSetFileName_GLVis);
631:   return(0);
632: }

634: /* this is a private implementation of a SOCKET with ASCII data format
635:    GLVis does not currently handle binary socket streams */
636: #if defined(PETSC_HAVE_UNISTD_H)
637: #include <unistd.h>
638: #endif

640: #if !defined(PETSC_HAVE_WINDOWS_H)
641: static PetscErrorCode (*PetscViewerDestroy_ASCII)(PetscViewer);

643: static PetscErrorCode PetscViewerDestroy_ASCII_Socket(PetscViewer viewer)
644: {
645:   FILE *stream;
646:   PetscErrorCode 0;
648:   PetscViewerASCIIGetPointer(viewer,&stream);
649:   if (stream) {
650:     fclose(stream);
651:     if (ierr) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SYS,"fclose() failed on stream");
652:   }
653:   PetscViewerDestroy_ASCII(viewer);
654:   return(0);
655: }
656: #endif

658: static PetscErrorCode PetscViewerASCIISocketOpen(MPI_Comm comm,const char* hostname,PetscInt port,PetscViewer* viewer)
659: {
660: #if defined(PETSC_HAVE_WINDOWS_H)
662:   SETERRQ(comm,PETSC_ERR_SUP,"Not implemented for Windows");
663: #else
664:   FILE           *stream = NULL;
665:   int            fd;

671: #if defined(PETSC_USE_SOCKET_VIEWER)
672:   PetscOpenSocket(hostname,port,&fd);
673: #else
674:   SETERRQ(comm,PETSC_ERR_SUP,"Missing Socket viewer");
675: #endif
676:   if (ierr) {
677:     PetscInt sierr;
678:     char     err[1024];

680:     PetscSNPrintf(err,1024,"Cannot connect to socket on %s:%D. Socket visualization is disabled\n",hostname,port);
681:     PetscInfo(NULL,err);
682:     *viewer = NULL;
683:     PetscFunctionReturn(sierr);
684:   } else {
685:     char msg[1024];

687:     PetscSNPrintf(msg,1024,"Successfully connect to socket on %s:%D. Socket visualization is enabled\n",hostname,port);
688:     PetscInfo(NULL,msg);
689:   }
690:   stream = fdopen(fd,"w"); /* Not possible on Windows */
691:   if (!stream) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_SYS,"Cannot open stream from socket %s:%d",hostname,port);
692:   PetscViewerASCIIOpenWithFILE(PETSC_COMM_SELF,stream,viewer);
693:   PetscViewerDestroy_ASCII = (*viewer)->ops->destroy;
694:   (*viewer)->ops->destroy = PetscViewerDestroy_ASCII_Socket;
695: #endif
696:   return(0);
697: }