1    | /***************************************
2    |   $Header: /cvsroot/petscgraphics/petsc.c,v 1.8 2003/05/20 13:58:11 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   | 
53   | int DATriangulate
54   | (DA theda, Vec globalX, int this, PetscScalar *minmax, int n_quants,
55   |  PetscScalar *isoquants, PetscScalar *colors)
56   | {
57   |   Vec localX;
58   |   PetscScalar isoquantdefaults[4],
59   |     colordefaults[16] = { 1,0,0,.5, 1,1,0,.5, 0,1,0,.5, 0,0,1,.5 };
60   |   PetscReal themax, themin;
61   |   int ierr;
62   | 
63   |   /* Default isoquants */
64   |   if (n_quants == PETSC_DECIDE) {
65   |     ierr = VecStrideMin (globalX, this, PETSC_NULL, &themin); CHKERRQ (ierr);
66   |     ierr = VecStrideMax (globalX, this, PETSC_NULL, &themax); CHKERRQ (ierr);
67   |     /* Red, yellow, green, blue at 0.2, 0.4, 0.6, 0.8, all with alpha=0.5 */
68   |     n_quants = 4;
69   |     isoquantdefaults[0] = themin + 0.2*(themax-themin);
70   |     isoquantdefaults[1] = themin + 0.4*(themax-themin);
71   |     isoquantdefaults[2] = themin + 0.6*(themax-themin);
72   |     isoquantdefaults[3] = themin + 0.8*(themax-themin);
73   |     isoquants = isoquantdefaults;
74   |     colors = colordefaults;
75   |   }
76   | 
77   |   /* Get the local array of points, with ghosts */
78   |   ierr = DACreateLocalVector (theda, &localX); CHKERRQ (ierr);
79   |   ierr = DAGlobalToLocalBegin (theda, globalX, INSERT_VALUES, localX); CHKERRQ(ierr);
80   |   ierr = DAGlobalToLocalEnd (theda, globalX, INSERT_VALUES, localX); CHKERRQ (ierr);
81   | 
82   |   /* Use DATriangulateLocal() to do the work */
83   |   ierr = DATriangulateLocal (theda, localX, this, minmax, n_quants, isoquants,
84   | 			     colors); CHKERRQ (ierr);
85   | 
86   |   ierr = VecDestroy (localX); CHKERRQ (ierr);
87   | 
88   |   return 0;
89   | }
90   | 
91   | 
92   | #undef __FUNCT__
93   | #define __FUNCT__ "DATriangulateLocal"
94   | 
95   | /*++++++++++++++++++++++++++++++++++++++
96   |   Calculate vertices of isoquant triangles in a 3-D distributed array.  This
97   |   takes a PETSc DA object, does some sanity checks, calculates array sizes, and
98   |   then gets array and passes it to Draw3DBlock for triangulation.
99   | 
100  |   int DATriangulateLocal Returns 0 or an error code.
101  | 
102  |   DA theda The PETSc distributed array object.
103  | 
104  |   Vec localX PETSc local vector object associated with the DA with data we'd
105  |   like to graph.
106  | 
107  |   int this Index of the field we'd like to draw.
108  | 
109  |   PetscScalar *minmax Position of block corners: xmin, xmax, ymin, ymax, zmin,
110  |   zmax.
111  | 
112  |   int n_quants Number of isoquant surfaces to draw (isoquant values).  Note
113  |   +latex+{\tt PETSC\_DECIDE}
114  |   +html+ <tt>PETSC\_DECIDE</tt>
115  |   is not a valid option here, because it's impossible to know the global
116  |   maximum and minimum and have consistent contours without user-supplied
117  |   information.
118  | 
119  |   PetscScalar *isoquants Array of function values at which to draw isoquants.
120  | 
121  |   PetscScalar *colors Array of color R,G,B,A quads for each isoquant.
122  |   ++++++++++++++++++++++++++++++++++++++*/
123  | 
124  | int DATriangulateLocal
125  | (DA theda, Vec localX, int this, PetscScalar *minmax, int n_quants,
126  |  PetscScalar *isoquants, PetscScalar *colors)
127  | {
128  |   PetscScalar *x, isoquantdefaults[4], localminmax[6],
129  |     colordefaults[16] = { 1,0,0,.5, 1,1,0,.5, 0,1,0,.5, 0,0,1,.5 };
130  |   DAStencilType stencil;
131  |   int i, ierr, fields, xw,yw,zw, xs,ys,zs, xm,ym,zm, gxs,gys,gzs, gxm,gym,gzm;
132  | 
133  |   /* Default isoquant error */
134  |   if (n_quants == PETSC_DECIDE)
135  |     {
136  |       SETERRQ (PETSC_ERR_ARG_WRONGSTATE, "Can't get default isoquants for local array");
137  |     }
138  | 
139  |   /* Get global and local grid boundaries */
140  |   ierr = DAGetInfo (theda, &i, &xw,&yw,&zw, PETSC_NULL,PETSC_NULL,PETSC_NULL,
141  | 		    &fields, PETSC_NULL, PETSC_NULL, &stencil);CHKERRQ(ierr);
142  |   if (i!=3)
143  |     {
144  |       SETERRQ (PETSC_ERR_ARG_WRONGSTATE, "DA must be 3-dimensional");
145  |     }
146  |   if (stencil!=DA_STENCIL_BOX)
147  |     {
148  |       SETERRQ (PETSC_ERR_ARG_WRONGSTATE, "DA must have a box stencil");
149  |     }
150  |   ierr = DAGetCorners (theda, &xs,&ys,&zs, &xm,&ym,&zm); CHKERRQ (ierr);
151  |   ierr = DAGetGhostCorners (theda, &gxs,&gys,&gzs, &gxm,&gym,&gzm);
152  |   CHKERRQ (ierr);
153  | 
154  |   /* Get the local array of points, with ghosts */
155  |   ierr = VecGetArray (localX, &x); CHKERRQ (ierr);
156  | 
157  |   /* Calculate local physical size */
158  |   localminmax[0] = minmax[0] + (minmax[1]-minmax[0])*xs/xw;
159  |   localminmax[1] = minmax[0] + (minmax[1]-minmax[0])*(xs+xm)/xw;
160  |   localminmax[2] = minmax[2] + (minmax[3]-minmax[2])*ys/yw;
161  |   localminmax[3] = minmax[2] + (minmax[3]-minmax[2])*(ys+ym)/yw;
162  |   localminmax[4] = minmax[4] + (minmax[5]-minmax[4])*zs/zw;
163  |   localminmax[5] = minmax[4] + (minmax[5]-minmax[4])*(zs+zm)/zw;
164  | 
165  |   /* If the array is too big, cut it down to size */
166  |   if (gxm <= xs-gxs+xm)
167  |     xm = gxm-xs+gxs-1;
168  |   if (gym <= ys-gys+ym)
169  |     ym = gym-ys+gys-1;
170  |   if (gzm <= zs-gzs+zm)
171  |     zm = gzm-zs+gzs-1;
172  | 
173  |   /* Let 'er rip! */
174  |   ierr = Draw3DBlock (gxm,gym,gzm, xs-gxs,ys-gys,zs-gzs, xm,ym,zm, localminmax,
175  | 		      x+this, fields, n_quants, isoquants, colors);
176  |   CHKERRQ (ierr);
177  |   ierr = VecRestoreArray (localX, &x); CHKERRQ (ierr);
178  | 
179  |   /* This debugging information may be useful to someone at some point;
180  |      note that it requires "extern int num_triangles;" somewhere above. */
181  |   /* {
182  |     int rank;
183  |     ierr = MPI_Comm_rank (PETSC_COMM_WORLD, &rank); CHKERRQ (ierr);
184  |     ierr = PetscSynchronizedPrintf
185  |       (PETSC_COMM_WORLD, "[%d] zs=%d, zm=%d, zmin=%g, zmax=%g\n",
186  |        rank, zs, zm, localminmax[4], localminmax[5]); CHKERRQ (ierr);
187  |     ierr = PetscSynchronizedFlush (PETSC_COMM_WORLD); CHKERRQ (ierr);
188  |     ierr = PetscSynchronizedPrintf (PETSC_COMM_WORLD, "[%d] Triangles: %d\n",
189  | 				    rank, num_triangles); CHKERRQ (ierr);
190  |     ierr = PetscSynchronizedFlush (PETSC_COMM_WORLD); CHKERRQ (ierr);
191  |     } */
192  | 
193  |   return 0;
194  | }
195  | 
196  | 
197  | #undef __FUNCT__
198  | #define __FUNCT__ "IllErrorHandler"
199  | 
200  | /*++++++++++++++++++++++++++++++++++++++
201  |   Handle errors, in this case the PETSc way.
202  | 
203  |   int IllErrorHandler Returns the error code supplied.
204  | 
205  |   int id Index of the error, defined in petscerror.h.
206  | 
207  |   char *message Text of the error message.
208  |   ++++++++++++++++++++++++++++++++++++++*/
209  | 
210  | int IllErrorHandler (int id, char *message)
211  | {
212  |   SETERRQ (id, message);
213  |   exit (1);
214  | }