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