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) PetscObjectQuery((PetscObject)U, "__PETSc_dm", &dm);
 68:       PetscObjectCompose((PetscObject)Ufield[i], "__PETSc_dm", dm);
 69:     }
 70:   }

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

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

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

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

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