1 /*  testSemi.c  */
2 
3 #include "../misc.h"
4 #include "../../SymbFac.h"
5 #include "../../timings.h"
6 
7 /*--------------------------------------------------------------------*/
8 int
main(int argc,char * argv[])9 main ( int argc, char *argv[] )
10 /*
11    ----------------------------------------
12    get statistics for a semi-implicit solve
13 
14    created -- 97dec11, cca
15    ----------------------------------------
16 */
17 {
18 char     *inGraphFileName, *inETreeFileName, *inMapFileName ;
19 double   nA21, nL, nL11, nL22, nPhi, nV, t1, t2 ;
20 ETree    *etree ;
21 int      ii, inside, J, K, msglvl, nfront, nJ,
22          nvtx, rc, sizeJ, v, vsize, w ;
23 int      *adjJ, *frontmap, *map, *nodwghts,
24          *vadj, *vtxToFront, *vwghts ;
25 IV       *mapIV ;
26 IVL      *symbfacIVL ;
27 Graph    *graph ;
28 FILE     *msgFile ;
29 Tree     *tree ;
30 
31 if ( argc != 6 ) {
32    fprintf(stdout,
33      "\n\n usage : %s msglvl msgFile GraphFile ETreeFile mapFile "
34      "\n    msglvl    -- message level"
35      "\n    msgFile   -- message file"
36      "\n    GraphFile -- input graph file, must be *.graphf or *.graphb"
37      "\n    ETreeFile -- input ETree file, must be *.etreef or *.etreeb"
38      "\n    mapFile   -- input map IV file, must be *.ivf or *.ivb"
39      "\n",
40 argv[0]) ;
41    return(0) ;
42 }
43 msglvl = atoi(argv[1]) ;
44 if ( strcmp(argv[2], "stdout") == 0 ) {
45    msgFile = stdout ;
46 } else if ( (msgFile = fopen(argv[2], "a")) == NULL ) {
47    fprintf(stderr, "\n fatal error in %s"
48            "\n unable to open file %s\n",
49            argv[0], argv[2]) ;
50    return(-1) ;
51 }
52 inGraphFileName = argv[3] ;
53 inETreeFileName = argv[4] ;
54 inMapFileName   = argv[5] ;
55 fprintf(msgFile,
56         "\n %s "
57         "\n msglvl        -- %d"
58         "\n msgFile       -- %s"
59         "\n GraphFile     -- %s"
60         "\n ETreeFile     -- %s"
61         "\n mapFile       -- %s"
62         "\n",
63         argv[0], msglvl, argv[2],
64         inGraphFileName, inETreeFileName, inMapFileName) ;
65 fflush(msgFile) ;
66 /*
67    ------------------------
68    read in the Graph object
69    ------------------------
70 */
71 graph = Graph_new() ;
72 if ( strcmp(inGraphFileName, "none") == 0 ) {
73    fprintf(msgFile, "\n no file to read from") ;
74    exit(0) ;
75 }
76 MARKTIME(t1) ;
77 rc = Graph_readFromFile(graph, inGraphFileName) ;
78 MARKTIME(t2) ;
79 fprintf(msgFile, "\n CPU %9.5f : read in graph from file %s",
80         t2 - t1, inGraphFileName) ;
81 nvtx   = graph->nvtx ;
82 vwghts = graph->vwghts ;
83 if ( rc != 1 ) {
84    fprintf(msgFile, "\n return value %d from Graph_readFromFile(%p,%s)",
85            rc, graph, inGraphFileName) ;
86    exit(-1) ;
87 }
88 if ( msglvl > 2 ) {
89    fprintf(msgFile, "\n\n after reading Graph object from file %s",
90            inGraphFileName) ;
91    Graph_writeForHumanEye(graph, msgFile) ;
92    fflush(msgFile) ;
93 }
94 /*
95    ------------------------
96    read in the ETree object
97    ------------------------
98 */
99 etree = ETree_new() ;
100 if ( strcmp(inETreeFileName, "none") == 0 ) {
101    fprintf(msgFile, "\n no file to read from") ;
102    exit(0) ;
103 }
104 MARKTIME(t1) ;
105 rc = ETree_readFromFile(etree, inETreeFileName) ;
106 MARKTIME(t2) ;
107 fprintf(msgFile, "\n CPU %9.5f : read in etree from file %s",
108         t2 - t1, inETreeFileName) ;
109 nfront     = ETree_nfront(etree) ;
110 tree       = ETree_tree(etree) ;
111 vtxToFront = ETree_vtxToFront(etree) ;
112 nodwghts   = ETree_nodwghts(etree) ;
113 nL         = ETree_nFactorEntries(etree, 2) ;
114 if ( rc != 1 ) {
115    fprintf(msgFile, "\n return value %d from ETree_readFromFile(%p,%s)",
116            rc, etree, inETreeFileName) ;
117    exit(-1) ;
118 }
119 if ( msglvl > 2 ) {
120    fprintf(msgFile, "\n\n after reading ETree object from file %s",
121            inETreeFileName) ;
122    ETree_writeForHumanEye(etree, msgFile) ;
123    fflush(msgFile) ;
124 }
125 /*
126    -------------------------
127    read in the map IV object
128    -------------------------
129 */
130 mapIV = IV_new() ;
131 if ( strcmp(inMapFileName, "none") == 0 ) {
132    fprintf(msgFile, "\n no file to read from") ;
133    exit(0) ;
134 }
135 MARKTIME(t1) ;
136 rc = IV_readFromFile(mapIV, inMapFileName) ;
137 MARKTIME(t2) ;
138 fprintf(msgFile, "\n CPU %9.5f : read in mapIV from file %s",
139         t2 - t1, inMapFileName) ;
140 map = IV_entries(mapIV) ;
141 if ( rc != 1 ) {
142    fprintf(msgFile, "\n return value %d from IV_readFromFile(%p,%s)",
143            rc, mapIV, inMapFileName) ;
144    exit(-1) ;
145 }
146 if ( msglvl > 2 ) {
147    fprintf(msgFile, "\n\n after reading IV object from file %s",
148            inMapFileName) ;
149    IV_writeForHumanEye(mapIV, msgFile) ;
150    fflush(msgFile) ;
151 }
152 nV = nPhi = 0 ;
153 if ( vwghts == NULL ) {
154    for ( v = 0 ; v < nvtx ; v++ ) {
155       nV++ ;
156       if ( map[v] == 0 ) {
157          nPhi++ ;
158       }
159    }
160 } else {
161    for ( v = 0 ; v < nvtx ; v++ ) {
162       nV += vwghts[v] ;
163       if ( map[v] == 0 ) {
164          nPhi += vwghts[v] ;
165       }
166    }
167 }
168 fprintf(msgFile, "\n nPhi = %.0f, nV = %.0f", nPhi, nV) ;
169 /*
170    -------------------------
171    get the frontmap[] vector
172    -------------------------
173 */
174 frontmap = IVinit(nfront, -1) ;
175 for ( v = 0 ; v < nvtx ; v++ ) {
176    J = vtxToFront[v] ;
177    if ( frontmap[J] == -1 ) {
178       frontmap[J] = map[v] ;
179    } else if ( frontmap[J] != map[v] ) {
180       fprintf(msgFile, "\n\n error, frontmap[%d] = %d, map[%d] = %d",
181               J, frontmap[J], v, map[v]) ;
182    }
183 }
184 /*
185    ----------------------------------
186    compute the symbolic factorization
187    ----------------------------------
188 */
189 symbfacIVL = SymbFac_initFromGraph(etree, graph) ;
190 if ( msglvl > 2 ) {
191    fprintf(msgFile, "\n\n symbolic factorization") ;
192    IVL_writeForHumanEye(symbfacIVL, msgFile) ;
193    fflush(msgFile) ;
194 }
195 /*
196    --------------------------------------------
197    compute the number of entries in L11 and L22
198    --------------------------------------------
199 */
200 nL11 = nL22 = 0 ;
201 for ( J = Tree_postOTfirst(tree) ;
202       J != -1 ;
203       J = Tree_postOTnext(tree, J) ) {
204    nJ = nodwghts[J] ;
205    if ( msglvl > 3 ) {
206       fprintf(msgFile, "\n\n front %d, nJ = %d", J, nJ) ;
207    }
208    IVL_listAndSize(symbfacIVL, J, &sizeJ, &adjJ) ;
209    for ( ii = 0, inside = 0 ; ii < sizeJ ; ii++ ) {
210       w = adjJ[ii] ;
211       K = vtxToFront[w] ;
212       if ( msglvl > 3 ) {
213          fprintf(msgFile, "\n    w = %d, K = %d", w, K) ;
214       }
215       if ( K > J && frontmap[K] == frontmap[J] ) {
216          inside += (vwghts == NULL) ? 1 : vwghts[w] ;
217          if ( msglvl > 3 ) {
218             fprintf(msgFile, ", inside") ;
219          }
220       }
221    }
222    if ( frontmap[J] != 0 ) {
223       if ( msglvl > 3 ) {
224          fprintf(msgFile, "\n    inside = %d, adding %d to L11",
225                  inside, nJ*nJ + 2*nJ*inside) ;
226       }
227       nL11 += nJ*nJ + 2*nJ*inside ;
228    } else {
229       if ( msglvl > 3 ) {
230          fprintf(msgFile, "\n    inside = %d, adding %d to L22",
231                  inside, nJ*nJ + 2*nJ*inside) ;
232       }
233       nL22 += nJ*nJ + 2*nJ*inside ;
234    }
235 }
236 if ( msglvl > 0 ) {
237    fprintf(msgFile, "\n |L| = %.0f, |L11| = %.0f, |L22| = %.0f",
238            nL, nL11, nL22) ;
239 }
240 /*
241    ------------------------------------
242    compute the number of entries in A21
243    ------------------------------------
244 */
245 nA21 = 0 ;
246 if ( vwghts != NULL ) {
247    for ( v = 0 ; v < nvtx ; v++ ) {
248       if ( map[v] == 0 ) {
249          Graph_adjAndSize(graph, v, &vsize, &vadj) ;
250          for ( ii = 0 ; ii < vsize ; ii++ ) {
251             w = vadj[ii] ;
252             if ( map[v] != map[w] ) {
253                if ( msglvl > 3 ) {
254                   fprintf(msgFile, "\n A21 : v = %d, w = %d", v, w) ;
255                }
256                nA21 += vwghts[v] * vwghts[w] ;
257             }
258          }
259       }
260    }
261 } else {
262    for ( v = 0 ; v < nvtx ; v++ ) {
263       if ( map[v] == 0 ) {
264          Graph_adjAndSize(graph, v, &vsize, &vadj) ;
265          for ( ii = 0 ; ii < vsize ; ii++ ) {
266             w = vadj[ii] ;
267             if ( map[v] != map[w] ) {
268                if ( msglvl > 3 ) {
269                   fprintf(msgFile, "\n A21 : v = %d, w = %d", v, w) ;
270                }
271                nA21++ ;
272             }
273          }
274       }
275    }
276 }
277 if ( msglvl > 0 ) {
278    fprintf(msgFile,
279            "\n |L| = %.0f, |L11| = %.0f, |L22| = %.0f, |A21| = %.0f",
280            nL, nL11, nL22, nA21) ;
281    fprintf(msgFile,
282       "\n storage: explicit = %.0f, semi-implicit = %.0f, ratio = %.3f"
283       "\n opcount: explicit = %.0f, semi-implicit = %.0f, ratio = %.3f",
284       nL, nL11 + nA21 + nL22,
285       nL/(nL11 + nA21 + nL22),
286       2*nL, 4*nL11 + 2*nA21 + 2*nL22,
287       2*nL/(4*nL11 + 2*nA21 + 2*nL22)) ;
288    fprintf(msgFile, "\n ratios %8.3f %8.3f %8.3f",
289            nPhi/nV,
290            nL/(nL11 + nA21 + nL22),
291            2*nL/(4*nL11 + 2*nA21 + 2*nL22)) ;
292 }
293 /*
294    ------------------------
295    free the working storage
296    ------------------------
297 */
298 Graph_free(graph) ;
299 ETree_free(etree) ;
300 IV_free(mapIV) ;
301 IVL_free(symbfacIVL) ;
302 IVfree(frontmap) ;
303 
304 fprintf(msgFile, "\n") ;
305 fclose(msgFile) ;
306 
307 return(1) ; }
308 
309 /*--------------------------------------------------------------------*/
310