Actual source code: dmglvis.c

  1: /* Routines to visualize DMs through GLVis */

  3: #include <petsc/private/dmimpl.h>
  4: #include <petsc/private/glvisviewerimpl.h>

  6: PetscErrorCode DMView_GLVis(DM dm, PetscViewer viewer, PetscErrorCode (*DMView_GLVis_ASCII)(DM, PetscViewer))
  7: {
  8:   PetscBool isglvis, isascii;

 10:   PetscFunctionBegin;
 13:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERGLVIS, &isglvis));
 14:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
 15:   PetscCheck(isglvis || isascii, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Viewer must be of type VIEWERGLVIS or VIEWERASCII");
 16:   if (isglvis) {
 17:     PetscViewerGLVisType type;
 18:     PetscViewer          view;

 20:     PetscCall(PetscViewerGLVisGetType_Private(viewer, &type));
 21:     PetscCall(PetscViewerGLVisGetDMWindow_Private(viewer, &view));
 22:     if (!view) PetscFunctionReturn(PETSC_SUCCESS); /* socket window has been closed */
 23:     if (type == PETSC_VIEWER_GLVIS_SOCKET) {
 24:       PetscMPIInt size, rank;
 25:       PetscInt    sdim;
 26:       const char *name;

 28:       PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)dm), &size));
 29:       PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
 30:       PetscCall(DMGetCoordinateDim(dm, &sdim));
 31:       PetscCall(PetscObjectGetName((PetscObject)dm, &name));

 33:       PetscCall(PetscGLVisCollectiveBegin(PetscObjectComm((PetscObject)dm), &view));
 34:       PetscCall(PetscViewerASCIIPrintf(view, "parallel %d %d\nmesh\n", size, rank));
 35:       PetscCall(DMView_GLVis_ASCII(dm, view));
 36:       PetscCall(PetscViewerGLVisInitWindow_Private(view, PETSC_TRUE, sdim, name));
 37:       PetscCall(PetscGLVisCollectiveEnd(PetscObjectComm((PetscObject)dm), &view));
 38:     } else {
 39:       PetscCall(DMView_GLVis_ASCII(dm, view));
 40:     }
 41:     PetscCall(PetscViewerGLVisRestoreDMWindow_Private(viewer, &view));
 42:   } else {
 43:     PetscCall(DMView_GLVis_ASCII(dm, viewer));
 44:   }
 45:   PetscFunctionReturn(PETSC_SUCCESS);
 46: }