1    | /***************************************
2    |   $Header: /cvsroot/petscgraphics/geomview.c,v 1.14 2006/02/05 21:15:52 hazelsct Exp $
3    | 
4    |   This file has the Geomview interface, including the PETSc vector gather
5    |   operations to bring everything to CPU 0.
6    | ***************************************/
7    | 
8    | 
9    | #include <stdio.h>
10   | #include <petscvec.h>
11   | #include "config.h"      /* esp. for inline */
12   | #include "illuminator.h" /* Just to make sure the interface is "right" */
13   | #include "surface.h"     /* For isurface structure definition */
14   | #include <stdlib.h>      /* For malloc() */
15   | #include <unistd.h>      /* For execlp() */
16   | 
17   | /* Build with -DDEBUG for debugging output */
18   | #undef DPRINTF
19   | #ifdef DEBUG
20   | #define DPRINTF(fmt, args...) PetscPrintf (PETSC_COMM_WORLD, "%s: " fmt, __FUNCT__, args)
21   | #else
22   | #define DPRINTF(fmt, args...)
23   | #endif
24   | 
25   | 
26   | #undef __FUNCT__
27   | #define __FUNCT__ "GeomviewBegin"
28   | 
29   | /*++++++++++++++++++++++++++++++++++++++
30   |   Spawn a new geomview process.  Most of this was shamelessly ripped from Ken
31   |   Brakke's Surface Evolver.
32   | 
33   |   int GeomviewBegin Returns 0 or an error code.
34   | 
35   |   MPI_Comm comm MPI communicator for rank information, if NULL it uses
36   |   PETSC_COMM_WORLD.
37   | 
38   |   IDisplay *newdisp Pointer in which to put a new IDisplay object.
39   |   ++++++++++++++++++++++++++++++++++++++*/
40   | 
41   | int GeomviewBegin(MPI_Comm comm, IDisplay *newdisp)
42   | {
43   |   int rank, ierr, gv_pid,gv_pipe[2], to_gv_pipe[2];
44   |   char gv_version[100];
45   |   struct idisplay *thedisp;
46   | 
47   |   if (!(thedisp = (struct idisplay *) malloc (sizeof (struct idisplay))))
48   |     SETERRQ (PETSC_ERR_MEM, "Unable to allocate memory for isurface object");
49   | 
50   |   /* Initialize display image to zero */
51   |   thedisp->rgb = NULL;
52   |   thedisp->rgb_width = thedisp->rgb_height = thedisp->rgb_rowskip =
53   |     thedisp->rgb_bpp = 0;
54   |   thedisp->zbuf = NULL;
55   |   thedisp->zbuf_width = thedisp->zbuf_height = thedisp->zbuf_depth =
56   |     thedisp->zbuf_rowskip = 0;
57   | 
58   |   if (!comm)
59   |     comm = PETSC_COMM_WORLD;
60   |   ierr = MPI_Comm_rank (comm,&rank); CHKERRQ(ierr);
61   |   if (!rank)
62   |     {
63   |       pipe (gv_pipe); /* from geomview stdout */
64   |       pipe (to_gv_pipe); /* to geomview stdin */
65   |       gv_pid = fork ();
66   | 
67   |       if(gv_pid==0) /* child */
68   | 	{
69   | 	  close (0);
70   | 	  dup (to_gv_pipe[0]);
71   | 	  close (to_gv_pipe[0]);
72   | 	  close (to_gv_pipe[1]);
73   | 	  close (1);
74   | 	  dup (gv_pipe[1]);
75   | 	  close (gv_pipe[0]);
76   | 	  close (gv_pipe[1]);
77   | 	  /* signal(SIGINT,SIG_IGN); */
78   | 	  execlp (GEOMVIEW,GEOMVIEW,"-c","(interest (pick world))","-",NULL);
79   | 	  perror (GEOMVIEW); /* only error gets here */
80   | 	  SETERRQ (PETSC_ERR_ARG_WRONGSTATE,"Geomview invocation failed.\n");
81   | 	}
82   | 
83   |       /* PETSc program execution resumes here */
84   |       close (gv_pipe[1]);
85   |       close (to_gv_pipe[0]);
86   |       thedisp->to_geomview = fdopen (to_gv_pipe[1], "w"); /* geomview stdin */
87   |       fprintf (thedisp->to_geomview, "(echo (geomview-version) \"\n\")\n");
88   |       fflush (thedisp->to_geomview);
89   |       read (gv_pipe[0], gv_version, sizeof (gv_version));
90   |     }
91   | 
92   |   *newdisp = (IDisplay) thedisp;
93   |   return 0;
94   | }
95   | 
96   | 
97   | #undef __FUNCT__
98   | #define __FUNCT__ "GeomviewEnd"
99   | 
100  | /*++++++++++++++++++++++++++++++++++++++
101  |   Exit the current running Geomview process and close its pipe.  Based in part
102  |   on Ken Brakke's Surface Evolver.
103  | 
104  |   int GeomviewEnd Returns 0 or an error code.
105  | 
106  |   MPI_Comm comm MPI communicator for rank information, if NULL it uses
107  |   PETSC_COMM_WORLD.
108  | 
109  |   IDisplay Disp IDisplay object whose geomview FILE we're closing.
110  |   ++++++++++++++++++++++++++++++++++++++*/
111  | 
112  | int GeomviewEnd (MPI_Comm comm, IDisplay Disp)
113  | {
114  |   int rank, ierr;
115  |   struct idisplay *thedisp = (struct idisplay *) Disp;
116  | 
117  |   if (!comm)
118  |     comm = PETSC_COMM_WORLD;
119  |   ierr = MPI_Comm_rank (comm, &rank); CHKERRQ (ierr);
120  |   if (!rank)
121  |     {
122  |       if (thedisp->to_geomview)
123  | 	{
124  | 	  fprintf (thedisp->to_geomview, "(exit)");
125  | 	  fclose (thedisp->to_geomview);
126  | 	  thedisp->to_geomview = NULL;
127  | 	}
128  |     }
129  |   return 0;
130  | }
131  | 
132  | 
133  | #undef __FUNCT__
134  | #define __FUNCT__ "GeomviewDisplayTriangulation"
135  | 
136  | /*++++++++++++++++++++++++++++++++++++++
137  |   Pipe the current triangulation to Geomview for display.  Much of this is
138  |   based on Ken Brakke's Surface Evolver.
139  | 
140  |   int GeomviewDisplayTriangulation Returns 0 or an error code.
141  | 
142  |   MPI_Comm comm MPI communicator for rank information, if NULL it uses
143  |   PETSC_COMM_WORLD.
144  | 
145  |   ISurface Surf ISurface object containing triangles to render using Geomview.
146  | 
147  |   IDisplay Disp IDisplay object whose geomview FILE we're closing.
148  | 
149  |   PetscScalar *minmax Position of block corners: xmin, xmax, ymin, ymax, zmin, zmax.
150  | 
151  |   char *name Name to give the Geomview OOGL object which we create.
152  | 
153  |   PetscTruth transparent Geomview transparency flag.
154  |   ++++++++++++++++++++++++++++++++++++++*/
155  | 
156  | int GeomviewDisplayTriangulation
157  | (MPI_Comm comm, ISurface Surf, IDisplay Disp, PetscScalar *minmax, char *name,
158  |  PetscTruth transparent)
159  | {
160  |   int ierr, i, total_entries, start, end, rank;
161  |   PetscScalar *localvals;
162  |   Vec globalverts, localverts;
163  |   VecScatter vecscat;
164  |   IS islocal;
165  |   struct isurface *thesurf = (struct isurface *) Surf;
166  |   struct idisplay *thedisp = (struct idisplay *) Disp;
167  | 
168  |   /* Simple "assertion" check */
169  |   if (!comm)
170  |     comm = PETSC_COMM_WORLD;
171  |   ierr = MPI_Comm_rank (comm, &rank); CHKERRQ (ierr);
172  |   if (!rank && !thedisp->to_geomview)
173  |     SETERRQ (PETSC_ERR_ARG_WRONGSTATE,"Geomview display is not open!");
174  | 
175  |   /*+ First, this creates global and local vectors for all of the triangle
176  |     vertices. +*/
177  |   ierr = VecCreateMPIWithArray (comm, 13*thesurf->num_triangles,
178  | 				PETSC_DETERMINE,
179  | 				thesurf->vertices, &globalverts); CHKERRQ (ierr);
180  |   ierr = VecGetSize (globalverts, &total_entries); CHKERRQ (ierr);
181  |   ierr = VecCreateMPI (comm, (rank == 0) ? total_entries:0, total_entries,
182  | 		       &localverts); CHKERRQ (ierr);
183  |   DPRINTF ("Total triangles: %d\n", total_entries/13);
184  | 
185  |   /* Diagnostics which may be useful to somebody sometime */
186  |   /* ierr = VecGetOwnershipRange (globalverts, &start, &end); CHKERRQ (ierr);
187  |   ierr = PetscSynchronizedPrintf
188  |     (comm, "[%d] Global vector local size: %d, ownership from %d to %d\n",
189  |      rank, end-start, start, end); CHKERRQ (ierr);
190  |   ierr = PetscSynchronizedFlush (comm); CHKERRQ (ierr);
191  | 
192  |   ierr = VecGetOwnershipRange (localverts, &start, &end); CHKERRQ (ierr);
193  |   ierr = PetscSynchronizedPrintf
194  |     (comm, "[%d] Local vector local size: %d, ownership from %d to %d\n", rank,
195  |      end-start, start, end); CHKERRQ (ierr);
196  |   ierr = PetscSynchronizedFlush (comm); CHKERRQ (ierr);
197  |   ierr = PetscPrintf (comm, "Global vector size: %d\n", total_entries);
198  |   CHKERRQ (ierr); */
199  | 
200  |   /*+ It then gathers (``scatters'') all vertex data to processor zero, +*/
201  |   ierr = ISCreateStride (comm, total_entries, 0, 1, &islocal); CHKERRQ (ierr);
202  |   ierr = VecScatterCreate (globalverts, PETSC_NULL, localverts, islocal,
203  | 			   &vecscat); CHKERRQ (ierr);
204  |   DPRINTF ("Starting triangle scatter to head node\n",0);
205  |   ierr = VecScatterBegin (globalverts, localverts, INSERT_VALUES,
206  | 			  SCATTER_FORWARD, vecscat); CHKERRQ (ierr);
207  |   DPRINTF ("Triangle scatter to head node in progress...\n",0);
208  |   ierr = VecScatterEnd (globalverts, localverts, INSERT_VALUES,
209  | 			SCATTER_FORWARD, vecscat); CHKERRQ (ierr);
210  |   DPRINTF ("Triangle scatter to head node complete\n",0);
211  |   ierr = VecScatterDestroy (vecscat); CHKERRQ (ierr);
212  |   ierr = ISDestroy (islocal); CHKERRQ (ierr);
213  | 
214  |   /*+ and puts them in an array. +*/
215  |   ierr = VecGetArray (localverts, &localvals); CHKERRQ (ierr);
216  | 
217  |   /*+ Finally, it sends everything to Geomview, +*/
218  |   if (!rank) {
219  |     fprintf (thedisp->to_geomview, "(geometry \"%s\" { : bor })", name);
220  |     fprintf (thedisp->to_geomview, "(read geometry { define bor \n");
221  |     fprintf (thedisp->to_geomview, "appearance {%s}\nOFF\n",
222  | 	     transparent ? "+transparent" : "");
223  |     fprintf (thedisp->to_geomview, "%d %d 0\n", 3*(total_entries/13) + 2,
224  | 	     total_entries/13);
225  |     for (i=0; i<total_entries; i+=13)
226  |       fprintf (thedisp->to_geomview, "%g %g %g\n%g %g %g\n%g %g %g\n",
227  | 	       localvals[i], localvals[i+1], localvals[i+2], localvals[i+3],
228  | 	       localvals[i+4], localvals[i+5], localvals[i+6], localvals[i+7],
229  | 	       localvals[i+8]);
230  |     fprintf (thedisp->to_geomview, "%g %g %g\n%g %g %g\n", minmax[0],minmax[2],
231  | 	     minmax[4], minmax[1], minmax[3], minmax[5]);
232  |     for (i=0; i<total_entries/13; i++)
233  |       fprintf (thedisp->to_geomview,"3 %d %d %d %g %g %g %g\n",3*i,3*i+1,3*i+2,
234  | 	       localvals[13*i+9], localvals[13*i+10], localvals[13*i+11],
235  | 	       localvals[13*i+12]);
236  |     fprintf (thedisp->to_geomview, "})\n");
237  |     fflush (thedisp->to_geomview);
238  |   }
239  | 
240  |   /*+ and cleans up the mess. +*/
241  |   ierr = VecRestoreArray (localverts, &localvals); CHKERRQ (ierr);
242  |   ierr = VecDestroy (localverts); CHKERRQ (ierr);
243  |   ierr = VecDestroy (globalverts); CHKERRQ (ierr);
244  | 
245  |   return 0;
246  | }