1 #include <petsc/private/glvisviewerimpl.h> 2 #include <petsc/private/glvisvecimpl.h> 3 4 static PetscErrorCode PetscViewerGLVisVecInfoDestroy_Private(void *ptr) 5 { 6 PetscViewerGLVisVecInfo info = (PetscViewerGLVisVecInfo)ptr; 7 8 PetscFunctionBeginUser; 9 PetscCall(PetscFree(info->fec_type)); 10 PetscCall(PetscFree(info)); 11 PetscFunctionReturn(PETSC_SUCCESS); 12 } 13 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; 25 26 PetscFunctionBegin; 27 PetscCall(PetscViewerGLVisGetStatus_Private(viewer, &sockstatus)); 28 if (sockstatus == PETSCVIEWERGLVIS_DISABLED) PetscFunctionReturn(PETSC_SUCCESS); 29 /* if the user did not customize the viewer through the API, we need extra data that can be attached to the Vec */ 30 PetscCall(PetscViewerGLVisGetFields_Private(viewer, &nfields, NULL, NULL, NULL, NULL, NULL)); 31 if (!nfields) { 32 PetscObject dm; 33 34 PetscCall(PetscObjectQuery((PetscObject)U, "__PETSc_dm", &dm)); 35 if (dm) { 36 PetscCall(PetscViewerGLVisSetDM_Private(viewer, dm)); 37 } else SETERRQ(PetscObjectComm((PetscObject)U), PETSC_ERR_SUP, "You need to provide a DM or use PetscViewerGLVisSetFields()"); 38 } 39 PetscCall(PetscViewerGLVisGetFields_Private(viewer, &nfields, &fec_type, &spacedim, &g2lfields, (PetscObject **)&Ufield, &userctx)); 40 if (!nfields) PetscFunctionReturn(PETSC_SUCCESS); 41 42 PetscCall(PetscViewerGLVisGetType_Private(viewer, &socktype)); 43 44 for (i = 0; i < nfields; i++) { 45 PetscObject fdm; 46 PetscContainer container; 47 48 /* attach visualization info to the vector */ 49 PetscCall(PetscObjectQuery((PetscObject)Ufield[i], "_glvis_info_container", (PetscObject *)&container)); 50 if (!container) { 51 PetscViewerGLVisVecInfo info; 52 53 PetscCall(PetscNew(&info)); 54 PetscCall(PetscStrallocpy(fec_type[i], &info->fec_type)); 55 PetscCall(PetscContainerCreate(PetscObjectComm((PetscObject)U), &container)); 56 PetscCall(PetscContainerSetPointer(container, (void *)info)); 57 PetscCall(PetscContainerSetUserDestroy(container, PetscViewerGLVisVecInfoDestroy_Private)); 58 PetscCall(PetscObjectCompose((PetscObject)Ufield[i], "_glvis_info_container", (PetscObject)container)); 59 PetscCall(PetscContainerDestroy(&container)); 60 } 61 /* attach the mesh to the viz vectors */ 62 PetscCall(PetscObjectQuery((PetscObject)Ufield[i], "__PETSc_dm", &fdm)); 63 if (!fdm) { 64 PetscObject dm; 65 66 PetscCall(PetscViewerGLVisGetDM_Private(viewer, &dm)); 67 if (!dm) PetscCall(PetscObjectQuery((PetscObject)U, "__PETSc_dm", &dm)); 68 PetscCheck(dm, PetscObjectComm((PetscObject)U), PETSC_ERR_SUP, "Mesh not present"); 69 PetscCall(PetscObjectCompose((PetscObject)Ufield[i], "__PETSc_dm", dm)); 70 } 71 } 72 73 /* user-provided sampling */ 74 if (g2lfields) { 75 PetscCall((*g2lfields)((PetscObject)U, nfields, (PetscObject *)Ufield, userctx)); 76 } else { 77 PetscCheck(nfields <= 1, PetscObjectComm((PetscObject)U), PETSC_ERR_SUP, "Don't know how to sample %" PetscInt_FMT " fields", nfields); 78 PetscCall(VecCopy(U, Ufield[0])); 79 } 80 81 /* TODO callback to user routine to disable/enable subdomains */ 82 for (i = 0; i < nfields; i++) { 83 PetscObject dm; 84 PetscViewer view; 85 86 PetscCall(PetscObjectQuery((PetscObject)Ufield[i], "__PETSc_dm", &dm)); 87 PetscCall(PetscViewerGLVisGetWindow_Private(viewer, i, &view)); 88 if (!view) continue; /* socket window has been closed */ 89 if (socktype == PETSC_VIEWER_GLVIS_SOCKET) { 90 PetscMPIInt size, rank; 91 const char *name; 92 93 PetscCallMPI(MPI_Comm_size(PetscObjectComm(dm), &size)); 94 PetscCallMPI(MPI_Comm_rank(PetscObjectComm(dm), &rank)); 95 PetscCall(PetscObjectGetName((PetscObject)Ufield[i], &name)); 96 97 PetscCall(PetscGLVisCollectiveBegin(PetscObjectComm(dm), &view)); 98 PetscCall(PetscViewerASCIIPrintf(view, "parallel %d %d\nsolution\n", size, rank)); 99 PetscCall(PetscObjectView(dm, view)); 100 PetscCall(VecView(Ufield[i], view)); 101 PetscCall(PetscViewerGLVisInitWindow_Private(view, PETSC_FALSE, spacedim[i], name)); 102 PetscCall(PetscGLVisCollectiveEnd(PetscObjectComm(dm), &view)); 103 if (view) pause = PETSC_TRUE; /* at least one window is connected */ 104 } else { 105 PetscCall(PetscObjectView(dm, view)); 106 PetscCall(VecView(Ufield[i], view)); 107 } 108 PetscCall(PetscViewerGLVisRestoreWindow_Private(viewer, i, &view)); 109 } 110 if (pause) PetscCall(PetscViewerGLVisPause_Private(viewer)); 111 PetscFunctionReturn(PETSC_SUCCESS); 112 } 113