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