1 //------------------------------------------------------------------------------
2 // GraphBLAS/Demo/Program/mis_demo.c: maximal independent set
3 //------------------------------------------------------------------------------
4 
5 // SuiteSparse:GraphBLAS, Timothy A. Davis, (c) 2017-2021, All Rights Reserved.
6 // SPDX-License-Identifier: Apache-2.0
7 
8 //------------------------------------------------------------------------------
9 
10 // Read a graph from a file and compute a maximal independent set. Usage:
11 //
12 //  mis_demo < infile
13 //  mis_demo 0 nrows ncols ntuples method
14 //  mis_demo 1 nx ny method
15 //
16 // Where infile has one line per edge in the graph; these have the form
17 //
18 //  i j x
19 //
20 // where A(i,j)=x is performed by GrB_Matrix_build, to construct the matrix.
21 // The dimensions of A are assumed to be the largest row and column indices,
22 // plus one (in read_matrix.c).
23 //
24 // For the second usage (mis_demo 0 ...), a random symmetric matrix is created
25 // of size nrows-by-ncols with ntuples edges (some will be duplicates so the
26 // actual number of edges will be slightly less).  The method is 0 for
27 // setElement and 1 for build.  The matrix will not have any self-edges, which
28 // cause the mis method to fail.
29 //
30 // The 3rd usage (mis_demo 1 ...) creates a finite-element matrix on an
31 // nx-by-ny grid.  Method is 0 to 3; refer to wathen.c for details.
32 
33 // macro used by OK(...) to free workspace if an error occurs
34 #define FREE_ALL                        \
35     GrB_Matrix_free (&A) ;              \
36     GrB_Matrix_free (&C) ;              \
37     GrB_Vector_free (&iset1) ;          \
38     GrB_Vector_free (&iset2) ;          \
39     GrB_Vector_free (&e) ;              \
40     GrB_Descriptor_free (&dt) ;         \
41     GrB_Monoid_free (&Lor) ;            \
42     GrB_Semiring_free (&Boolean) ;      \
43     if (I2 != NULL) free (I2) ;         \
44     if (J2 != NULL) free (J2) ;         \
45     if (X2 != NULL) free (X2) ;         \
46     if (I != NULL) free (I) ;           \
47     if (X != NULL) free (X) ;
48 
49 #include "graphblas_demos.h"
50 
51 int64_t isize1, isize2 ;
52 GrB_Index n ;
53 GrB_Matrix A = NULL, C = NULL ;
54 GrB_Vector iset1 = NULL, iset2 = NULL, e = NULL ;
55 GrB_Descriptor dt = NULL ;
56 GrB_Index *I = NULL, *I2 = NULL, *J2 = NULL ;
57 float *X = NULL ;
58 bool *X2 = NULL ;
59 GrB_Monoid Lor = NULL ;
60 GrB_Semiring Boolean = NULL ;
61 GrB_Info info ;
62 
63 //------------------------------------------------------------------------------
64 // mis_check_results: test if iset is a maximal independent set
65 //------------------------------------------------------------------------------
66 
mis_check_results(int64_t * p_isize,GrB_Vector iset,double t)67 GrB_Info mis_check_results
68 (
69     int64_t *p_isize,
70     GrB_Vector iset,
71     double t
72 )
73 {
74 
75     //--------------------------------------------------------------------------
76     // report the results
77     //--------------------------------------------------------------------------
78 
79     int64_t isize ;
80     OK (GrB_Vector_reduce_INT64 (&isize, NULL, GrB_PLUS_MONOID_INT64,
81         iset, NULL)) ;
82     (*p_isize) = isize ;
83 
84     GrB_Index nvals ;
85     OK (GrB_Vector_nvals (&nvals, iset)) ;
86     I = (GrB_Index *) malloc (nvals * sizeof (GrB_Index)) ;
87     X = (float *) malloc (nvals * sizeof (float)) ;
88 
89     if (I == NULL || X == NULL)
90     {
91         printf ("out of memory\n") ;
92         FREE_ALL ;
93         exit (1) ;
94     }
95 
96     OK (GrB_Vector_extractTuples_FP32 (I, X, &nvals, iset)) ;
97 
98     // I [0..isize-1] is the independent set
99     isize = 0 ;
100     for (int64_t k = 0 ; k < nvals ; k++)
101     {
102         if (X [k])
103         {
104             // printf ("I [%lld] = %lld\n", isize, I [k]) ;
105             I [isize++] = I [k] ;
106         }
107     }
108 
109     free (X) ; X = NULL ;
110 
111     printf ("independent set found: %.16g of %.16g nodes\n",
112         (double) isize, (double) n) ;
113 
114     //--------------------------------------------------------------------------
115     // verify the result
116     //--------------------------------------------------------------------------
117 
118     OK (GrB_Matrix_new (&C, GrB_BOOL, isize, isize)) ;
119     OK (GrB_Matrix_extract (C, NULL, NULL, A, I, isize, I, isize, NULL)) ;
120     OK (GrB_Matrix_nvals (&nvals, C)) ;
121 
122     I2 = (GrB_Index *) malloc (nvals * sizeof (GrB_Index)) ;
123     J2 = (GrB_Index *) malloc (nvals * sizeof (GrB_Index)) ;
124     X2 = (bool *) malloc (nvals * sizeof (bool)) ;
125     if (I2 == NULL || J2 == NULL || X2 == NULL)
126     {
127         printf ("out of memory\n") ;
128         FREE_ALL ;
129         exit (1) ;
130     }
131 
132     // could do this with a mask instead of extractTuples.
133     OK (GrB_Matrix_extractTuples_BOOL (I2, J2, X2, &nvals, C)) ;
134     GrB_Matrix_free (&C) ;
135 
136     for (int64_t k = 0 ; k < nvals ; k++)
137     {
138         if (X2 [k] && I2 [k] != J2 [k])
139         {
140             printf ("error!  A(I,I) has an edge!\n") ;
141             FREE_ALL ;
142             exit (1) ;
143         }
144     }
145 
146     free (I2) ; I2 = NULL ;
147     free (J2) ; J2 = NULL ;
148     free (X2) ; X2 = NULL ;
149 
150     // now check if all other nodes are adjacent to the iset
151 
152     // e = iset
153     OK (GrB_Vector_dup (&e, iset)) ;
154     // e = (e || A*iset), using the Boolean semiring
155     OK (GrB_vxm (e, NULL, GrB_LOR, GxB_LOR_LAND_BOOL, iset, A, NULL)) ;
156     OK (GrB_Vector_nvals (&nvals, e)) ;
157     GrB_Vector_free (&e) ;
158     if (nvals != n)
159     {
160         fprintf (stderr, "error! A (I,I is not maximal!\n") ;
161         exit (1) ;
162     }
163 
164     free (I) ; I = NULL ;
165 
166     fprintf (stderr, "maximal independent set OK %.16g of %.16g nodes"
167         " time: %g\n", (double) isize, (double) n, t) ;
168     return (GrB_SUCCESS) ;
169 }
170 
171 //------------------------------------------------------------------------------
172 // mis_demo main
173 //------------------------------------------------------------------------------
174 
main(int argc,char ** argv)175 int main (int argc, char **argv)
176 {
177 
178     double tic [2], t1, t2 ;
179     OK (GrB_init (GrB_NONBLOCKING)) ;
180     int nthreads ;
181     OK (GxB_Global_Option_get (GxB_GLOBAL_NTHREADS, &nthreads)) ;
182     fprintf (stderr, "\nmis_demo: nthreads: %d\n", nthreads) ;
183 
184     //--------------------------------------------------------------------------
185     // get a symmetric matrix with no self edges
186     //--------------------------------------------------------------------------
187 
188     OK (get_matrix (&A, argc, argv, true, true, true)) ;
189     OK (GrB_Matrix_nrows (&n, A)) ;
190 
191     //--------------------------------------------------------------------------
192     // convert A to boolean
193     //--------------------------------------------------------------------------
194 
195     OK (GrB_Descriptor_new (&dt)) ;
196     OK (GxB_Desc_set (dt, GrB_INP0, GrB_TRAN)) ;
197     OK (GrB_Matrix_new (&C, GrB_BOOL, n, n)) ;
198     OK (GrB_transpose (C, NULL, NULL, A, dt)) ;
199     GrB_Descriptor_free (&dt) ;
200     GrB_Matrix_free (&A) ;
201     A = C ;
202     C = NULL ;
203 
204     //--------------------------------------------------------------------------
205     // compute maximal independent set
206     //--------------------------------------------------------------------------
207 
208     for (int64_t seed = 1 ; seed <= 2 ; seed++)
209     {
210 
211         simple_tic (tic) ;
212         OK (mis (&iset1, A, seed)) ;
213         t1 = simple_toc (tic) ;
214         printf ("MIS time in seconds: %14.6f\n", t1) ;
215 
216         // also test the version that checks every error condition
217         simple_tic (tic) ;
218         OK (mis_check (&iset2, A, seed)) ;
219         t2 = simple_toc (tic) ;
220         printf ("MIS time in seconds: %14.6f\n", t2) ;
221 
222         //----------------------------------------------------------------------
223         // compare results
224         //----------------------------------------------------------------------
225 
226         mis_check_results (&isize1, iset1, t1) ;
227         mis_check_results (&isize2, iset2, t2) ;
228 
229         printf ("isize: %.16g %.16g\n", (double) isize1, (double) isize2) ;
230 
231         if (isize1 != isize2)
232         {
233             fprintf (stderr, "=============================================="
234                 "======size differs!\n") ;
235             printf ("size differs!\n") ;
236         }
237         GrB_Vector_free (&iset1) ;
238         GrB_Vector_free (&iset2) ;
239     }
240 
241     FREE_ALL ;
242     GrB_finalize ( ) ;
243 }
244 
245