1 /* Routines to visualize DMs through GLVis */ 2 3 #include <petsc/private/dmimpl.h> 4 #include <petsc/private/glvisviewerimpl.h> 5 6 PetscErrorCode DMView_GLVis(DM dm, PetscViewer viewer, PetscErrorCode (*DMView_GLVis_ASCII)(DM,PetscViewer)) 7 { 8 PetscErrorCode ierr; 9 PetscBool isglvis,isascii; 10 11 PetscFunctionBegin; 12 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 13 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 14 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERGLVIS,&isglvis);CHKERRQ(ierr); 15 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);CHKERRQ(ierr); 16 if (!isglvis && !isascii) SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Viewer must be of type VIEWERGLVIS or VIEWERASCII"); 17 if (isglvis) { 18 PetscViewerGLVisType type; 19 PetscViewer view; 20 21 ierr = PetscViewerGLVisGetType_Private(viewer,&type);CHKERRQ(ierr); 22 ierr = PetscViewerGLVisGetDMWindow_Private(viewer,&view);CHKERRQ(ierr); 23 if (!view) PetscFunctionReturn(0); /* socket window has been closed */ 24 if (type == PETSC_VIEWER_GLVIS_SOCKET) { 25 PetscMPIInt size,rank; 26 PetscInt sdim; 27 const char* name; 28 29 ierr = MPI_Comm_size(PetscObjectComm((PetscObject)dm),&size);CHKERRMPI(ierr); 30 ierr = MPI_Comm_rank(PetscObjectComm((PetscObject)dm),&rank);CHKERRMPI(ierr); 31 ierr = DMGetCoordinateDim(dm,&sdim);CHKERRQ(ierr); 32 ierr = PetscObjectGetName((PetscObject)dm,&name);CHKERRQ(ierr); 33 34 ierr = PetscGLVisCollectiveBegin(PetscObjectComm((PetscObject)dm),&view);CHKERRQ(ierr); 35 ierr = PetscViewerASCIIPrintf(view,"parallel %d %d\nmesh\n",size,rank);CHKERRQ(ierr); 36 ierr = DMView_GLVis_ASCII(dm,view);CHKERRQ(ierr); 37 ierr = PetscViewerGLVisInitWindow_Private(view,PETSC_TRUE,sdim,name);CHKERRQ(ierr); 38 ierr = PetscGLVisCollectiveEnd(PetscObjectComm((PetscObject)dm),&view);CHKERRQ(ierr); 39 } else { 40 ierr = DMView_GLVis_ASCII(dm,view);CHKERRQ(ierr); 41 } 42 ierr = PetscViewerGLVisRestoreDMWindow_Private(viewer,&view);CHKERRQ(ierr); 43 } else { 44 ierr = DMView_GLVis_ASCII(dm,viewer);CHKERRQ(ierr); 45 } 46 PetscFunctionReturn(0); 47 } 48