Actual source code: vecglvis.c

  1: #include <petsc/private/glvisviewerimpl.h>
  2: #include <petsc/private/glvisvecimpl.h>

  4: static PetscErrorCode PetscViewerGLVisVecInfoDestroy_Private(void *ptr)
  5: {
  6:   PetscViewerGLVisVecInfo info = (PetscViewerGLVisVecInfo)ptr;

  9:   PetscFree(info->fec_type);
 10:   PetscFree(info);
 11:   return 0;
 12: }

 14: /* the main function to visualize vectors using GLVis */
 15: PetscErrorCode VecView_GLVis(Vec U,PetscViewer viewer)
 16: {
 17:   PetscErrorCode         (*g2lfields)(PetscObject,PetscInt,PetscObject[],void*);
 18:   Vec                    *Ufield;
 19:   const char             **fec_type;
 20:   PetscViewerGLVisStatus sockstatus;
 21:   PetscViewerGLVisType   socktype;
 22:   void                   *userctx;
 23:   PetscInt               i,nfields,*spacedim;
 24:   PetscBool              pause = PETSC_FALSE;

 26:   PetscViewerGLVisGetStatus_Private(viewer,&sockstatus);
 27:   if (sockstatus == PETSCVIEWERGLVIS_DISABLED) return 0;
 28:   /* if the user did not customize the viewer through the API, we need extra data that can be attached to the Vec */
 29:   PetscViewerGLVisGetFields_Private(viewer,&nfields,NULL,NULL,NULL,NULL,NULL);
 30:   if (!nfields) {
 31:     PetscObject dm;

 33:     PetscObjectQuery((PetscObject)U, "__PETSc_dm",&dm);
 34:     if (dm) {
 35:       PetscViewerGLVisSetDM_Private(viewer,dm);
 36:     } else SETERRQ(PetscObjectComm((PetscObject)U),PETSC_ERR_SUP,"You need to provide a DM or use PetscViewerGLVisSetFields()");
 37:   }
 38:   PetscViewerGLVisGetFields_Private(viewer,&nfields,&fec_type,&spacedim,&g2lfields,(PetscObject**)&Ufield,&userctx);
 39:   if (!nfields) return 0;

 41:   PetscViewerGLVisGetType_Private(viewer,&socktype);

 43:   for (i=0;i<nfields;i++) {
 44:     PetscObject    fdm;
 45:     PetscContainer container;

 47:     /* attach visualization info to the vector */
 48:     PetscObjectQuery((PetscObject)Ufield[i],"_glvis_info_container",(PetscObject*)&container);
 49:     if (!container) {
 50:       PetscViewerGLVisVecInfo info;

 52:       PetscNew(&info);
 53:       PetscStrallocpy(fec_type[i],&info->fec_type);
 54:       PetscContainerCreate(PetscObjectComm((PetscObject)U),&container);
 55:       PetscContainerSetPointer(container,(void*)info);
 56:       PetscContainerSetUserDestroy(container,PetscViewerGLVisVecInfoDestroy_Private);
 57:       PetscObjectCompose((PetscObject)Ufield[i],"_glvis_info_container",(PetscObject)container);
 58:       PetscContainerDestroy(&container);
 59:     }
 60:     /* attach the mesh to the viz vectors */
 61:     PetscObjectQuery((PetscObject)Ufield[i], "__PETSc_dm",&fdm);
 62:     if (!fdm) {
 63:       PetscObject dm;

 65:       PetscViewerGLVisGetDM_Private(viewer,&dm);
 66:       if (!dm) {
 67:         PetscObjectQuery((PetscObject)U, "__PETSc_dm",&dm);
 68:       }
 70:       PetscObjectCompose((PetscObject)Ufield[i], "__PETSc_dm",dm);
 71:     }
 72:   }

 74:   /* user-provided sampling */
 75:   if (g2lfields) {
 76:     (*g2lfields)((PetscObject)U,nfields,(PetscObject*)Ufield,userctx);
 77:   } else {
 79:     VecCopy(U,Ufield[0]);
 80:   }

 82:   /* TODO callback to user routine to disable/enable subdomains */
 83:   for (i=0; i<nfields; i++) {
 84:     PetscObject dm;
 85:     PetscViewer view;

 87:     PetscObjectQuery((PetscObject)Ufield[i], "__PETSc_dm",&dm);
 88:     PetscViewerGLVisGetWindow_Private(viewer,i,&view);
 89:     if (!view) continue; /* socket window has been closed */
 90:     if (socktype == PETSC_VIEWER_GLVIS_SOCKET) {
 91:       PetscMPIInt size,rank;
 92:       const char *name;

 94:       MPI_Comm_size(PetscObjectComm(dm),&size);
 95:       MPI_Comm_rank(PetscObjectComm(dm),&rank);
 96:       PetscObjectGetName((PetscObject)Ufield[i],&name);

 98:       PetscGLVisCollectiveBegin(PetscObjectComm(dm),&view);
 99:       PetscViewerASCIIPrintf(view,"parallel %d %d\nsolution\n",size,rank);
100:       PetscObjectView(dm,view);
101:       VecView(Ufield[i],view);
102:       PetscViewerGLVisInitWindow_Private(view,PETSC_FALSE,spacedim[i],name);
103:       PetscGLVisCollectiveEnd(PetscObjectComm(dm),&view);
104:       if (view) pause = PETSC_TRUE; /* at least one window is connected */
105:     } else {
106:       PetscObjectView(dm,view);
107:       VecView(Ufield[i],view);
108:     }
109:     PetscViewerGLVisRestoreWindow_Private(viewer,i,&view);
110:   }
111:   if (pause) PetscViewerGLVisPause_Private(viewer);
112:   return 0;
113: }