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 | }