1 | /***************************************
2 | $Header$
3 |
4 | This is the petsc.c main file. It has all of the PETSc-dependent functions.
5 | ***************************************/
6 |
7 |
8 | #include <stdlib.h> /* For exit() */
9 | #include <petscda.h>
10 | #include "config.h" /* esp. for inline */
11 | #include "illuminator.h" /* Just to make sure the interface is "right" */
12 | #include "structures.h" /* For the isurface structure definition */
13 |
14 |
15 | #undef __FUNCT__
16 | #define __FUNCT__ "DATriangulateRange"
17 |
18 | /*++++++++++++++++++++++++++++++++++++++
19 | Calculate vertices of isoquant triangles in a 3-D distributed array. This
20 | takes a PETSc DA object, does some sanity checks, calculates array sizes,
21 | gets the local vector and array, and then calls
22 | +latex+{\tt DATriangulateLocal()}
23 | +html+ <tt>DATriangulateLocal()</tt>
24 | to do the rest. Note that global array access (i.e. this function) is
25 | necessary for using default isoquant values, since we need to be able to
26 | calculate the maximum and minimum on the global array.
27 |
28 | int DATriangulateRange Returns 0 or an error code.
29 |
30 | ISurface Surf ISurface object into which to draw this DA's isoquants.
31 |
32 | DA theda The PETSc distributed array object.
33 |
34 | Vec globalX PETSc global vector object associated with the DA with data we'd
35 | like to graph.
36 |
37 | int this Index of the field we'd like to draw.
38 |
39 | PetscScalar *minmax Position of block corners: xmin, xmax, ymin, ymax, zmin,
40 | zmax.
41 |
42 | int n_quants Number of isoquant surfaces to draw (isoquant values), or
43 | +latex+{\tt PETSC\_DECIDE}
44 | +html+ <tt>PETSC_DECIDE</tt>
45 | to use red, yellow, green, blue at 0.2, 0.4, 0.6 and 0.8 between the vector's
46 | global minimum and maximum values.
47 |
48 | PetscScalar *isoquants Array of function values at which to draw isoquants,
49 | +latex+or {\tt PETSC\_NULL} if {\tt n\_quants=PETSC\_DECIDE}.
50 | +html+ or <tt>PETSC_NULL</tt> if <tt>n_quants=PETSC\_DECIDE</tt>.
51 |
52 | PetscScalar *colors Array of color R,G,B,A quads for each isoquant, or
53 | +latex+{\tt PETSC\_NULL} if {\tt n\_quants=PETSC\_DECIDE}.
54 | +html+ <tt>PETSC_NULL</tt> if <tt>n_quants=PETSC\_DECIDE</tt>.
55 |
56 | int xmin Smallest grid x-coordinate to render.
57 |
58 | int xmax Largest grid x-coordinate to render, -1 goes to full x maximum, -2
59 | in periodic systems goes to one short of x maximum.
60 |
61 | int ymin Smallest grid y-coordinate to render.
62 |
63 | int ymax Largest grid y-coordinate to render, -1 goes to full y maximum, -2
64 | in periodic systems goes to one short of y maximum.
65 |
66 | int zmin Smallest grid z-coordinate to render.
67 |
68 | int zmax Largest grid z-coordinate to render, -1 goes to full z maximum, -2
69 | in periodic systems goes to one short of z maximum.
70 | ++++++++++++++++++++++++++++++++++++++*/
71 |
72 | int DATriangulateRange
73 | (ISurface Surf, DA theda, Vec globalX, int this, PetscScalar *minmax,
74 | int n_quants, PetscScalar *isoquants, PetscScalar *colors, int xmin,int xmax,
75 | int ymin,int ymax, int zmin,int zmax)
76 | {
77 | Vec localX;
78 | PetscScalar isoquantdefaults[4],
79 | colordefaults[16] = { 1,0,0,.5, 1,1,0,.5, 0,1,0,.5, 0,0,1,.5 };
80 | PetscReal themax, themin;
81 | int ierr;
82 |
83 | /* Default isoquants */
84 | if (n_quants == PETSC_DECIDE) {
85 | ierr = VecStrideMin (globalX, this, PETSC_NULL, &themin); CHKERRQ (ierr);
86 | ierr = VecStrideMax (globalX, this, PETSC_NULL, &themax); CHKERRQ (ierr);
87 | /* Red, yellow, green, blue at 0.2, 0.4, 0.6, 0.8, all with alpha=0.5 */
88 | n_quants = 4;
89 | isoquantdefaults[0] = themin + 0.2*(themax-themin);
90 | isoquantdefaults[1] = themin + 0.4*(themax-themin);
91 | isoquantdefaults[2] = themin + 0.6*(themax-themin);
92 | isoquantdefaults[3] = themin + 0.8*(themax-themin);
93 | isoquants = isoquantdefaults;
94 | colors = colordefaults;
95 | }
96 |
97 | /* Get the local array of points, with ghosts */
98 | ierr = DACreateLocalVector (theda, &localX); CHKERRQ (ierr);
99 | ierr = DAGlobalToLocalBegin (theda, globalX, INSERT_VALUES, localX); CHKERRQ(ierr);
100 | ierr = DAGlobalToLocalEnd (theda, globalX, INSERT_VALUES, localX); CHKERRQ (ierr);
101 |
102 | /* Use DATriangulateLocalRange() to do the work */
103 | ierr = DATriangulateLocalRange (Surf, theda, localX, this, minmax, n_quants,
104 | isoquants, colors, xmin,xmax, ymin,ymax,
105 | zmin,zmax); CHKERRQ (ierr);
106 |
107 | ierr = VecDestroy (localX); CHKERRQ (ierr);
108 |
109 | return 0;
110 | }
111 |
112 |
113 | #undef __FUNCT__
114 | #define __FUNCT__ "DATriangulateLocalRange"
115 |
116 | /*++++++++++++++++++++++++++++++++++++++
117 | Calculate vertices of isoquant triangles in a 3-D distributed array. This
118 | takes a PETSc DA object, does some sanity checks, calculates array sizes, and
119 | then gets array and passes it to Draw3DBlock for triangulation.
120 |
121 | int DATriangulateLocalRange Returns 0 or an error code.
122 |
123 | ISurface Surf ISurface object into which to draw this DA's isoquants.
124 |
125 | DA theda The PETSc distributed array object.
126 |
127 | Vec localX PETSc local vector object associated with the DA with data we'd
128 | like to graph.
129 |
130 | int this Index of the field we'd like to draw.
131 |
132 | PetscScalar *minmax Position of block corners: xmin, xmax, ymin, ymax, zmin,
133 | zmax.
134 |
135 | int n_quants Number of isoquant surfaces to draw (isoquant values). Note
136 | +latex+{\tt PETSC\_DECIDE}
137 | +html+ <tt>PETSC_DECIDE</tt>
138 | is not a valid option here, because it's impossible to know the global
139 | maximum and minimum and have consistent contours without user-supplied
140 | information.
141 |
142 | PetscScalar *isoquants Array of function values at which to draw isoquants.
143 |
144 | PetscScalar *colors Array of color R,G,B,A quads for each isoquant.
145 |
146 | int xmin Smallest grid x-coordinate to render.
147 |
148 | int xmax Largest grid x-coordinate to render, -1 goes to full x maximum, -2
149 | in periodic systems goes to one short of x maximum.
150 |
151 | int ymin Smallest grid y-coordinate to render.
152 |
153 | int ymax Largest grid y-coordinate to render, -1 goes to full y maximum, -2
154 | in periodic systems goes to one short of y maximum.
155 |
156 | int zmin Smallest grid z-coordinate to render.
157 |
158 | int zmax Largest grid z-coordinate to render, -1 goes to full z maximum, -2
159 | in periodic systems goes to one short of z maximum.
160 | ++++++++++++++++++++++++++++++++++++++*/
161 |
162 | int DATriangulateLocalRange
163 | (ISurface Surf, DA theda, Vec localX, int this, PetscScalar *minmax,
164 | int n_quants, PetscScalar *isoquants, PetscScalar *colors, int xmin,int xmax,
165 | int ymin,int ymax, int zmin,int zmax)
166 | {
167 | PetscScalar *x, isoquantdefaults[4], localminmax[6],
168 | colordefaults[16] = { 1,0,0,.5, 1,1,0,.5, 0,1,0,.5, 0,0,1,.5 };
169 | DAStencilType stencil;
170 | int i, ierr, fields, xw,yw,zw, xs,ys,zs, xm,ym,zm, gxs,gys,gzs, gxm,gym,gzm;
171 | struct isurface *thesurf = (struct isurface *) Surf;
172 |
173 | /* Default isoquant error */
174 | if (n_quants == PETSC_DECIDE)
175 | SETERRQ (PETSC_ERR_ARG_WRONGSTATE, "Can't get default isoquants for local array");
176 |
177 | /* Get global and local grid boundaries */
178 | ierr = DAGetInfo (theda, &i, &xw,&yw,&zw, PETSC_NULL,PETSC_NULL,PETSC_NULL,
179 | &fields, PETSC_NULL, PETSC_NULL, &stencil);CHKERRQ(ierr);
180 | if (i!=3)
181 | SETERRQ (PETSC_ERR_ARG_WRONGSTATE, "DA must be 3-dimensional");
182 | if (stencil!=DA_STENCIL_BOX)
183 | SETERRQ (PETSC_ERR_ARG_WRONGSTATE, "DA must have a box stencil");
184 |
185 | ierr = DAGetCorners (theda, &xs,&ys,&zs, &xm,&ym,&zm); CHKERRQ (ierr);
186 | ierr = DAGetGhostCorners (theda, &gxs,&gys,&gzs, &gxm,&gym,&gzm);
187 | CHKERRQ (ierr);
188 |
189 | /* Get the local array of points, with ghosts */
190 | ierr = VecGetArray (localX, &x); CHKERRQ (ierr);
191 |
192 | /* If the array is too big, cut it down to size */
193 | if (gxm <= xs-gxs+xm)
194 | xm = gxm-xs+gxs-1;
195 | if (gym <= ys-gys+ym)
196 | ym = gym-ys+gys-1;
197 | if (gzm <= zs-gzs+zm)
198 | zm = gzm-zs+gzs-1;
199 | /* Eliminate final rows/planes if *cut and periodic. */
200 | if (xmax == -2 && xs+xm>=xw)
201 | xm--;
202 | if (ymax == -2 && ys+ym>=yw)
203 | ym--;
204 | if (zmax == -2 && zs+zm>=zw)
205 | zm--;
206 |
207 | /* Reframe to range */
208 | xs = PetscMax (xs, xmin);
209 | ys = PetscMax (ys, xmin);
210 | zs = PetscMax (zs, xmin);
211 | xm = (xmax > 0) ? PetscMin (xm, xmax-xs) : xm;
212 | ym = (ymax > 0) ? PetscMin (ym, ymax-ys) : ym;
213 | zm = (zmax > 0) ? PetscMin (zm, zmax-zs) : zm;
214 |
215 | /* Calculate local physical size */
216 | localminmax[0] = minmax[0] + (minmax[1]-minmax[0])*xs/xw;
217 | localminmax[1] = minmax[2] + (minmax[1]-minmax[0])*(xs+xm)/xw;
218 | localminmax[2] = minmax[2] + (minmax[3]-minmax[2])*ys/yw;
219 | localminmax[3] = minmax[2] + (minmax[3]-minmax[2])*(ys+ym)/yw;
220 | localminmax[4] = minmax[4] + (minmax[5]-minmax[4])*zs/zw;
221 | localminmax[5] = minmax[4] + (minmax[5]-minmax[4])*(zs+zm)/zw;
222 |
223 | /* Let 'er rip! */
224 | ierr = Draw3DBlock (Surf, gxm,gym,gzm, xs-gxs,ys-gys,zs-gzs, xm,ym,zm,
225 | localminmax, x+this, fields, n_quants,isoquants, colors);
226 | CHKERRQ (ierr);
227 | ierr = VecRestoreArray (localX, &x); CHKERRQ (ierr);
228 |
229 | /* Prints the number of triangles on all CPUs */
230 | #ifdef DEBUG
231 | {
232 | int rank;
233 | ierr = MPI_Comm_rank (PETSC_COMM_WORLD, &rank); CHKERRQ (ierr);
234 | ierr = PetscSynchronizedPrintf
235 | (PETSC_COMM_WORLD, "[%d] zs=%d, zm=%d, zmin=%g, zmax=%g\n",
236 | rank, zs, zm, localminmax[4], localminmax[5]); CHKERRQ (ierr);
237 | ierr = PetscSynchronizedFlush (PETSC_COMM_WORLD); CHKERRQ (ierr);
238 | ierr = PetscSynchronizedPrintf (PETSC_COMM_WORLD, "[%d] Triangles: %d\n",
239 | rank, thesurf->num_triangles);
240 | CHKERRQ (ierr);
241 | ierr = PetscSynchronizedFlush (PETSC_COMM_WORLD); CHKERRQ (ierr);
242 | }
243 | #endif
244 |
245 | return 0;
246 | }
247 |
248 |
249 | #undef __FUNCT__
250 | #define __FUNCT__ "IllErrorHandler"
251 |
252 | /*++++++++++++++++++++++++++++++++++++++
253 | Handle errors, in this case the PETSc way.
254 |
255 | int IllErrorHandler Returns the error code supplied.
256 |
257 | int id Index of the error, defined in petscerror.h.
258 |
259 | char *message Text of the error message.
260 | ++++++++++++++++++++++++++++++++++++++*/
261 |
262 | int IllErrorHandler (int id, char *message)
263 | {
264 | int ierr;
265 | ierr = PetscPrintf (PETSC_COMM_WORLD, "Warning: IllErrorHandler is broken and deprecated, use PETSc SETERRQ instead.\n");
266 | SETERRQ (id, message);
267 | exit (1);
268 | }