Actual source code: vecglvis.c
petsc-3.9.2 2018-05-20
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: #if defined(PETSC_HAVE_SETJMP_H) && !defined(PETSC_MISSING_SIGPIPE)
16: #include <setjmp.h>
17: #include <signal.h>
19: #if defined(PETSC_HAVE_WINDOWS_H)
20: #define DEV_NULL "NUL"
21: #else
22: #define DEV_NULL "/dev/null"
23: #endif
25: static jmp_buf PetscGLVisSigPipeJmpBuf;
27: static void PetscGLVisSigPipeHandler(PETSC_UNUSED int sig)
28: {
29: longjmp(PetscGLVisSigPipeJmpBuf,1);
30: }
31: #endif
33: /* the main function to visualize vectors using GLVis */
34: PetscErrorCode VecView_GLVis(Vec U,PetscViewer viewer)
35: {
36: PetscErrorCode (*g2lfields)(PetscObject,PetscInt,PetscObject[],void*);
37: Vec *Ufield;
38: const char **fec_type;
39: PetscViewerGLVisStatus sockstatus;
40: PetscViewerGLVisType socktype;
41: void *userctx;
42: PetscInt i,nfields,*spacedim;
43: PetscErrorCode ierr;
46: PetscViewerGLVisGetStatus_Private(viewer,&sockstatus);
47: if (sockstatus == PETSCVIEWERGLVIS_DISABLED) return(0);
48: /* if the user did not customize the viewer through the API, we need extra data that can be attached to the Vec */
49: PetscViewerGLVisGetFields_Private(viewer,&nfields,NULL,NULL,NULL,NULL,NULL);
50: if (!nfields) {
51: PetscObject dm;
53: PetscObjectQuery((PetscObject)U, "__PETSc_dm",&dm);
54: if (dm) {
55: PetscViewerGLVisSetDM_Private(viewer,dm);
56: } else SETERRQ(PetscObjectComm((PetscObject)U),PETSC_ERR_SUP,"You need to provide a DM or use PetscViewerGLVisSetFields()");
57: }
58: PetscViewerGLVisGetFields_Private(viewer,&nfields,&fec_type,&spacedim,&g2lfields,(PetscObject**)&Ufield,&userctx);
59: if (!nfields) return(0);
61: PetscViewerGLVisGetType_Private(viewer,&socktype);
63: for (i=0;i<nfields;i++) {
64: PetscObject fdm;
65: PetscContainer container;
67: /* attach visualization info to the vector */
68: PetscObjectQuery((PetscObject)Ufield[i],"_glvis_info_container",(PetscObject*)&container);
69: if (!container) {
70: PetscViewerGLVisVecInfo info;
72: PetscNew(&info);
73: PetscStrallocpy(fec_type[i],&info->fec_type);
74: PetscContainerCreate(PetscObjectComm((PetscObject)U),&container);
75: PetscContainerSetPointer(container,(void*)info);
76: PetscContainerSetUserDestroy(container,PetscViewerGLVisVecInfoDestroy_Private);
77: PetscObjectCompose((PetscObject)Ufield[i],"_glvis_info_container",(PetscObject)container);
78: PetscContainerDestroy(&container);
79: }
80: /* attach the mesh to the viz vectors */
81: PetscObjectQuery((PetscObject)Ufield[i], "__PETSc_dm",&fdm);
82: if (!fdm) {
83: PetscObject dm;
85: PetscViewerGLVisGetDM_Private(viewer,&dm);
86: if (!dm) {
87: PetscObjectQuery((PetscObject)U, "__PETSc_dm",&dm);
88: }
89: if (!dm) SETERRQ(PetscObjectComm((PetscObject)U),PETSC_ERR_SUP,"Mesh not present");
90: PetscObjectCompose((PetscObject)Ufield[i], "__PETSc_dm",dm);
91: }
92: }
94: /* user-provided sampling */
95: if (g2lfields) {
96: (*g2lfields)((PetscObject)U,nfields,(PetscObject*)Ufield,userctx);
97: } else {
98: VecCopy(U,Ufield[0]);
99: }
101: /* TODO callback to user routine to disable/enable subdomains */
102: for (i=0;i<nfields;i++) {
103: PetscObject dm;
104: PetscViewer view;
106: PetscObjectQuery((PetscObject)Ufield[i], "__PETSc_dm",&dm);
107: PetscViewerGLVisGetWindow_Private(viewer,i,&view);
108: if (!view) continue; /* socket window has been closed */
109: if (socktype == PETSC_VIEWER_GLVIS_DUMP) {
110: PetscObjectView(dm,view);
111: VecView(Ufield[i],view);
112: } else {
113: /* It may happen that the user has closed the GLVis window */
114: #if defined(PETSC_HAVE_SETJMP_H) && !defined(PETSC_MISSING_SIGPIPE)
115: void (*sighdl)(int) = signal(SIGPIPE,PetscGLVisSigPipeHandler);
116: if (!setjmp(PetscGLVisSigPipeJmpBuf)) {
117: #endif
118: PetscMPIInt size,rank;
119: const char *name;
121: MPI_Comm_size(PetscObjectComm(dm),&size);
122: MPI_Comm_rank(PetscObjectComm(dm),&rank);
123: PetscViewerASCIIPrintf(view,"parallel %D %D\nsolution\n",size,rank);
124: PetscObjectView(dm,view);
125: VecView(Ufield[i],view);
126: PetscObjectGetName((PetscObject)Ufield[i],&name);
127: PetscViewerGLVisInitWindow_Private(view,PETSC_FALSE,spacedim[i],name);
128: #if defined(PETSC_HAVE_SETJMP_H) && !defined(PETSC_MISSING_SIGPIPE)
129: } else {
130: FILE *sock,*null = fopen(DEV_NULL,"w");
131: PetscInt readonly;
133: VecLockGet(Ufield[i],&readonly);
134: if (readonly) {
135: VecLockPop(Ufield[i]);
136: }
137: PetscViewerASCIIGetPointer(view,&sock);
138: PetscViewerASCIISetFILE(view,null);
139: PetscViewerDestroy(&view);
140: (void)fclose(sock);
141: }
142: (void)signal(SIGPIPE,sighdl);
143: #endif
144: }
145: PetscViewerGLVisRestoreWindow_Private(viewer,i,&view);
146: }
147: PetscViewerGLVisPause_Private(viewer);
148: return(0);
149: }