Actual source code: vecglvis.c

petsc-3.11.4 2019-09-28
Report Typos and Errors
  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;
  7:   PetscErrorCode          ierr;

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

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

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

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

 44:   PetscViewerGLVisGetType_Private(viewer,&socktype);

 46:   for (i=0;i<nfields;i++) {
 47:     PetscObject    fdm;
 48:     PetscContainer container;

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

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

 68:       PetscViewerGLVisGetDM_Private(viewer,&dm);
 69:       if (!dm) {
 70:         PetscObjectQuery((PetscObject)U, "__PETSc_dm",&dm);
 71:       }
 72:       if (!dm) SETERRQ(PetscObjectComm((PetscObject)U),PETSC_ERR_SUP,"Mesh not present");
 73:       PetscObjectCompose((PetscObject)Ufield[i], "__PETSc_dm",dm);
 74:     }
 75:   }

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

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

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

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

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