1 | /***************************************
2 | $Header: /cvsroot/petscgraphics/geomview.c,v 1.6 2003/04/30 19:34:55 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 |
14 | /* Build with -DDEBUG for debugging output */
15 | #undef DPRINTF
16 | #ifdef DEBUG
17 | #define DPRINTF(fmt, args...) PetscPrintf (PETSC_COMM_WORLD, "%s: " fmt, __FUNCT__, args)
18 | #else
19 | #define DPRINTF(fmt, args...)
20 | #endif
21 |
22 | /*+ Stream holding the pipe to the Geomview process. It is initialized to
23 | NULL, then set by GeomviewBegin, and cleared by GeomviewEnd. This can be
24 | used as a flag to determine whether the Geomview display is open. +*/
25 | static FILE *pfd = NULL;
26 |
27 | /*+ Declared in illuminator.c, this gives the current number of triangles on
28 | this node. +*/
29 | extern int num_triangles;
30 |
31 | /*+ Declared in illuminator.c, this array stores the corner coordinates and
32 | color information for each triangle. +*/
33 | extern PetscScalar vertices[];
34 |
35 |
36 | #undef __FUNCT__
37 | #define __FUNCT__ "GeomviewBegin"
38 |
39 | /*++++++++++++++++++++++++++++++++++++++
40 | Spawn a new geomview process. Most of this was shamelessly ripped from Ken
41 | Brakke's Surface Evolver.
42 |
43 | int GeomviewBegin Returns 0 or an error code.
44 |
45 | MPI_Comm comm MPI communicator for rank information, if NULL it uses
46 | PETSC_COMM_WORLD.
47 | ++++++++++++++++++++++++++++++++++++++*/
48 |
49 | int GeomviewBegin(MPI_Comm comm)
50 | {
51 | int rank, ierr, gv_pid,gv_pipe[2], to_gv_pipe[2];
52 | char gv_version[100];
53 |
54 | if (!comm)
55 | comm = PETSC_COMM_WORLD;
56 | ierr = MPI_Comm_rank (comm,&rank); CHKERRQ(ierr);
57 | if (!rank)
58 | {
59 | if (pfd)
60 | SETERRQ (PETSC_ERR_ARG_WRONGSTATE,
61 | "Multiple geomview displays not allowed");
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 | pfd = fdopen (to_gv_pipe[1], "w"); /* hooked to stdin of geomview */
87 | fprintf (pfd, "(echo (geomview-version) \"\n\")\n");
88 | fflush (pfd);
89 | read (gv_pipe[0], gv_version, sizeof (gv_version));
90 | }
91 | return 0;
92 | }
93 |
94 |
95 | #undef __FUNCT__
96 | #define __FUNCT__ "GeomviewEnd"
97 |
98 | /*++++++++++++++++++++++++++++++++++++++
99 | Exit the current running Geomview process and close its pipe. Based in part
100 | on Ken Brakke's Surface Evolver.
101 |
102 | int GeomviewEnd Returns 0 or an error code.
103 |
104 | MPI_Comm comm MPI communicator for rank information, if NULL it uses
105 | PETSC_COMM_WORLD.
106 | ++++++++++++++++++++++++++++++++++++++*/
107 |
108 | int GeomviewEnd (MPI_Comm comm)
109 | {
110 | int rank, ierr;
111 |
112 | if (!comm)
113 | comm = PETSC_COMM_WORLD;
114 | ierr = MPI_Comm_rank (comm, &rank); CHKERRQ (ierr);
115 | if (!rank)
116 | {
117 | if (pfd)
118 | {
119 | fprintf (pfd, "(exit)");
120 | fclose (pfd);
121 | pfd = NULL;
122 | }
123 | }
124 | return 0;
125 | }
126 |
127 |
128 | #undef __FUNCT__
129 | #define __FUNCT__ "GeomviewDisplayTriangulation"
130 |
131 | /*++++++++++++++++++++++++++++++++++++++
132 | Pipe the current triangulation to Geomview for display. Much of this is
133 | based on Ken Brakke's Surface Evolver.
134 |
135 | int GeomviewDisplayTriangulation Returns 0 or an error code.
136 |
137 | MPI_Comm comm MPI communicator for rank information, if NULL it uses
138 | PETSC_COMM_WORLD.
139 |
140 | PetscScalar *minmax Position of block corners: xmin, xmax, ymin, ymax, zmin, zmax.
141 |
142 | char *name Name to give the Geomview OOGL object which we create.
143 |
144 | PetscTruth transparent Geomview transparency flag.
145 | ++++++++++++++++++++++++++++++++++++++*/
146 |
147 | int GeomviewDisplayTriangulation
148 | (MPI_Comm comm, PetscScalar *minmax, char *name, PetscTruth transparent)
149 | {
150 | int ierr, i, total_entries, start, end, rank;
151 | PetscScalar *localvals;
152 | Vec globalverts, localverts;
153 | VecScatter vecscat;
154 | IS islocal;
155 |
156 | /* Simple "assertion" check */
157 | if (!comm)
158 | comm = PETSC_COMM_WORLD;
159 | ierr = MPI_Comm_rank (comm, &rank); CHKERRQ (ierr);
160 | if (!rank && !pfd)
161 | SETERRQ (PETSC_ERR_ARG_WRONGSTATE,"Geomview display is not open!");
162 |
163 | /*+ First, this creates global and local vectors for all of the triangle
164 | vertices. +*/
165 | ierr = VecCreateMPIWithArray (comm, 13*num_triangles, PETSC_DETERMINE,
166 | vertices, &globalverts); CHKERRQ (ierr);
167 | ierr = VecGetSize (globalverts, &total_entries); CHKERRQ (ierr);
168 | ierr = VecCreateMPI (comm, (rank == 0) ? total_entries:0, total_entries,
169 | &localverts); CHKERRQ (ierr);
170 | DPRINTF ("Total triangles: %d\n", total_entries/13);
171 |
172 | /* Diagnostics which may be useful to somebody sometime */
173 | /* ierr = VecGetOwnershipRange (globalverts, &start, &end); CHKERRQ (ierr);
174 | ierr = PetscSynchronizedPrintf
175 | (comm, "[%d] Global vector local size: %d, ownership from %d to %d\n",
176 | rank, end-start, start, end); CHKERRQ (ierr);
177 | ierr = PetscSynchronizedFlush (comm); CHKERRQ (ierr);
178 |
179 | ierr = VecGetOwnershipRange (localverts, &start, &end); CHKERRQ (ierr);
180 | ierr = PetscSynchronizedPrintf
181 | (comm, "[%d] Local vector local size: %d, ownership from %d to %d\n", rank,
182 | end-start, start, end); CHKERRQ (ierr);
183 | ierr = PetscSynchronizedFlush (comm); CHKERRQ (ierr);
184 | ierr = PetscPrintf (comm, "Global vector size: %d\n", total_entries);
185 | CHKERRQ (ierr); */
186 |
187 | /*+ It then gathers (``scatters'') all vertex data to processor zero, +*/
188 | ierr = ISCreateStride (comm, total_entries, 0, 1, &islocal); CHKERRQ (ierr);
189 | ierr = VecScatterCreate (globalverts, PETSC_NULL, localverts, islocal,
190 | &vecscat); CHKERRQ (ierr);
191 | DPRINTF ("Starting triangle scatter to head node\n",0);
192 | ierr = VecScatterBegin (globalverts, localverts, INSERT_VALUES,
193 | SCATTER_FORWARD, vecscat); CHKERRQ (ierr);
194 | DPRINTF ("Triangle scatter to head node in progress...\n",0);
195 | ierr = VecScatterEnd (globalverts, localverts, INSERT_VALUES,
196 | SCATTER_FORWARD, vecscat); CHKERRQ (ierr);
197 | DPRINTF ("Triangle scatter to head node complete\n",0);
198 | ierr = VecScatterDestroy (vecscat); CHKERRQ (ierr);
199 | ierr = ISDestroy (islocal); CHKERRQ (ierr);
200 |
201 | /*+ and puts them in an array. +*/
202 | ierr = VecGetArray (localverts, &localvals); CHKERRQ (ierr);
203 |
204 | /*+ Finally, it sends everything to Geomview, +*/
205 | if (!rank) {
206 | fprintf (pfd, "(geometry \"%s\" { : bor })", name);
207 | fprintf (pfd, "(read geometry { define bor \n");
208 | fprintf (pfd, "appearance {%s}\nOFF\n", transparent ? "+transparent" : "");
209 | fprintf (pfd, "%d %d 0\n", 3*(total_entries/13) + 2, total_entries/13);
210 | for (i=0; i<total_entries; i+=13)
211 | fprintf (pfd, "%g %g %g\n%g %g %g\n%g %g %g\n", localvals[i],
212 | localvals[i+1], localvals[i+2], localvals[i+3], localvals[i+4],
213 | localvals[i+5], localvals[i+6], localvals[i+7], localvals[i+8]);
214 | fprintf (pfd, "%g %g %g\n%g %g %g\n", minmax[0], minmax[2], minmax[4],
215 | minmax[1], minmax[3], minmax[5]);
216 | for (i=0; i<total_entries/13; i++)
217 | fprintf (pfd, "3 %d %d %d %g %g %g %g\n", 3*i,3*i+1,3*i+2,
218 | localvals[13*i+9], localvals[13*i+10], localvals[13*i+11],
219 | localvals[13*i+12]);
220 | fprintf (pfd, "})\n");
221 | fflush (pfd);
222 | }
223 |
224 | /*+ and cleans up the mess, resetting the number of triangles to zero. +*/
225 | ierr = VecRestoreArray (localverts, &localvals); CHKERRQ (ierr);
226 | ierr = VecDestroy (localverts); CHKERRQ (ierr);
227 | ierr = VecDestroy (globalverts); CHKERRQ (ierr);
228 | num_triangles = 0;
229 |
230 | return 0;
231 | }