Line data Source code
1 : //------------------------------------------------------------------------------
2 : // LAGraph_cc_lacc.c
3 : //------------------------------------------------------------------------------
4 :
5 : // LAGraph, (c) 2019-2022 by The LAGraph Contributors, All Rights Reserved.
6 : // SPDX-License-Identifier: BSD-2-Clause
7 : //
8 : // For additional details (including references to third party source code and
9 : // other files) see the LICENSE file or contact permission@sei.cmu.edu. See
10 : // Contributors.txt for a full list of contributors. Created, in part, with
11 : // funding and support from the U.S. Government (see Acknowledgments.txt file).
12 : // DM22-0790
13 :
14 : // Contributed by Scott McMillan, SEI, Carnegie Mellon University
15 :
16 : //------------------------------------------------------------------------------
17 :
18 : /**
19 : * Code is based on the algorithm described in the following paper Azad, Buluc;
20 : * LACC: a linear-algebraic algorithm for finding connected components in
21 : * distributed memory (IPDPS 2019)
22 : **/
23 :
24 : // TODO: is this ready for src? Good, except Reduce_assign.
25 : // TODO: Reduce_assign is slow. See src/algorithm/LG_CC_FastSV6.
26 :
27 : #include "LG_internal.h"
28 : #include <LAGraph.h>
29 : #include <LAGraphX.h>
30 :
31 : //****************************************************************************
32 : // mask = NULL, accumulator = GrB_MIN_UINT64, descriptor = NULL
33 :
34 : #undef LG_FREE_ALL
35 : #define LG_FREE_ALL LAGraph_Free ((void **) &mem, msg) ;
36 :
37 188 : static GrB_Info Reduce_assign (GrB_Vector w,
38 : GrB_Vector src,
39 : GrB_Index *index,
40 : GrB_Index nLocs, char *msg)
41 : {
42 : GrB_Index nw, ns;
43 188 : GrB_Vector_nvals(&nw, w);
44 188 : GrB_Vector_nvals(&ns, src);
45 188 : GrB_Index *mem = NULL ;
46 188 : LG_TRY (LAGraph_Malloc ((void **) &mem, nw*3, sizeof (GrB_Index), msg)) ;
47 188 : GrB_Index *ind = mem, *sval = mem + nw, *wval = sval + nw;
48 188 : GRB_TRY (GrB_Vector_extractTuples(ind, wval, &nw, w)) ;
49 188 : GRB_TRY (GrB_Vector_extractTuples(ind, sval, &ns, src)) ;
50 10346 : for (GrB_Index i = 0; i < nLocs; i++)
51 10158 : if (sval[i] < wval[index[i]])
52 9486 : wval[index[i]] = sval[i];
53 188 : GRB_TRY (GrB_Vector_clear(w)) ;
54 188 : GRB_TRY (GrB_Vector_build(w, ind, wval, nw, GrB_PLUS_UINT64)) ;
55 188 : LG_FREE_ALL ;
56 188 : return GrB_SUCCESS;
57 : }
58 :
59 : //****************************************************************************
60 :
61 : #undef LG_FREE_ALL
62 : #define LG_FREE_ALL \
63 : { \
64 : LAGraph_Free ((void **) &I, msg); \
65 : LAGraph_Free ((void **) &V, msg); \
66 : GrB_free (&S2) ; \
67 : GrB_free (&stars); \
68 : GrB_free (&mask); \
69 : GrB_free (&parents); \
70 : GrB_free (&gp); \
71 : GrB_free (&mnp); \
72 : GrB_free (&hookMNP); \
73 : GrB_free (&hookP); \
74 : GrB_free (&pNonstars); \
75 : GrB_free (&tmp); \
76 : GrB_free (&nsgp); \
77 : }
78 :
79 36 : int LAGraph_cc_lacc
80 : (
81 : GrB_Vector *result, // output: array of component identifiers
82 : GrB_Matrix A, // input matrix
83 : bool sanitize, // if true, ensure A is symmetric
84 : char *msg
85 : )
86 : {
87 : //--------------------------------------------------------------------------
88 : // check inputs
89 : //--------------------------------------------------------------------------
90 :
91 36 : LG_CLEAR_MSG ;
92 36 : if (result == NULL)
93 : {
94 12 : return (GrB_NULL_POINTER) ;
95 : }
96 :
97 : GrB_Info info;
98 :
99 24 : GrB_Vector stars = NULL, mask = NULL;
100 24 : GrB_Vector parents = NULL, gp = NULL, mnp = NULL;
101 24 : GrB_Vector hookMNP = NULL, hookP = NULL;
102 24 : GrB_Vector tmp = NULL, pNonstars = NULL, nsgp = NULL; // temporary
103 24 : GrB_Index *I = NULL;
104 24 : GrB_Index *V = NULL;
105 24 : GrB_Matrix S = NULL, S2 = NULL ;
106 :
107 : GrB_Index n ;
108 24 : GRB_TRY (GrB_Matrix_nrows (&n, A)) ;
109 :
110 24 : if (sanitize)
111 : {
112 12 : GRB_TRY (GrB_Matrix_new (&S2, GrB_BOOL, n, n)) ;
113 12 : GRB_TRY (GrB_eWiseAdd (S2, NULL, NULL, GrB_LOR, A, A, GrB_DESC_T1)) ;
114 12 : S = S2 ;
115 : }
116 : else
117 : {
118 : // Use the input as-is, and assume it is binary and symmetric
119 12 : S = A ;
120 : }
121 :
122 : // vectors
123 24 : GRB_TRY (GrB_Vector_new (&stars, GrB_BOOL, n));
124 24 : GRB_TRY (GrB_Vector_new (&mask, GrB_BOOL, n));
125 24 : GRB_TRY (GrB_Vector_new (&parents, GrB_UINT64, n));
126 24 : GRB_TRY (GrB_Vector_new (&gp, GrB_UINT64, n));
127 24 : GRB_TRY (GrB_Vector_new (&hookMNP, GrB_UINT64, n));
128 24 : GRB_TRY (GrB_Vector_new (&hookP, GrB_UINT64, n));
129 24 : GRB_TRY (GrB_Vector_new (&pNonstars, GrB_UINT64, n));
130 :
131 : // temporary arrays
132 24 : LG_TRY (LAGraph_Malloc ((void **) &I, n, sizeof (GrB_Index), msg)) ;
133 24 : LG_TRY (LAGraph_Malloc ((void **) &V, n, sizeof (GrB_Index), msg)) ;
134 :
135 : // prepare the vectors
136 16284 : for (GrB_Index i = 0 ; i < n ; i++)
137 16260 : I[i] = V[i] = i;
138 24 : GRB_TRY (GrB_Vector_build (parents, I, V, n, GrB_PLUS_UINT64));
139 24 : GRB_TRY (GrB_Vector_dup (&mnp, parents));
140 24 : GRB_TRY (GrB_assign (stars, 0, 0, true, GrB_ALL, 0, 0)) ;
141 :
142 : // main computation
143 : GrB_Index nHooks, nStars, nNonstars;
144 : while (true) {
145 : // ---------------------------------------------------------
146 : // CondHook(A, parents, stars);
147 : // ---------------------------------------------------------
148 94 : GRB_TRY (GrB_mxv (mnp, 0, 0, GrB_MIN_SECOND_SEMIRING_UINT64,
149 : S, parents, 0));
150 94 : GRB_TRY (GrB_Vector_clear (mask));
151 94 : GRB_TRY (GrB_eWiseMult(mask, stars, 0, GrB_LT_UINT64, mnp, parents, 0));
152 94 : GRB_TRY (GrB_assign (hookMNP, mask, 0, mnp, GrB_ALL, n, 0));
153 94 : GRB_TRY (GrB_eWiseMult (hookP, 0, 0, GrB_SECOND_UINT64, hookMNP, parents, 0));
154 94 : GRB_TRY (GrB_Vector_clear (mnp));
155 94 : GRB_TRY (GrB_Vector_nvals (&nHooks, hookP));
156 94 : GRB_TRY (GrB_Vector_extractTuples (I, V, &nHooks, hookP));
157 94 : GRB_TRY (GrB_Vector_new (&tmp, GrB_UINT64, nHooks));
158 94 : GRB_TRY (GrB_extract (tmp, 0, 0, hookMNP, I, nHooks, 0));
159 94 : LG_TRY (Reduce_assign (parents, tmp, V, nHooks, msg));
160 94 : GRB_TRY (GrB_Vector_clear (tmp));
161 : // modify the stars vector
162 94 : GRB_TRY (GrB_assign (stars, 0, 0, false, V, nHooks, 0));
163 94 : GRB_TRY (GrB_extract (tmp, 0, 0, parents, V, nHooks, 0)); // extract modified parents
164 94 : GRB_TRY (GrB_Vector_extractTuples (I, V, &nHooks, tmp));
165 94 : GRB_TRY (GrB_assign (stars, 0, 0, false, V, nHooks, 0));
166 94 : GRB_TRY (GrB_Vector_extractTuples (I, V, &n, parents));
167 94 : GRB_TRY (GrB_extract (mask, 0, 0, stars, V, n, 0));
168 94 : GRB_TRY (GrB_assign (stars, 0, GrB_LAND, mask, GrB_ALL, 0, 0));
169 : // clean up
170 94 : GRB_TRY (GrB_Vector_clear (hookMNP));
171 94 : GRB_TRY (GrB_Vector_clear (hookP));
172 94 : GRB_TRY (GrB_free (&tmp));
173 : // ---------------------------------------------------------
174 : // UnCondHook(A, parents, stars);
175 : // ---------------------------------------------------------
176 94 : GRB_TRY (GrB_assign (pNonstars, 0, 0, parents, GrB_ALL, 0, 0));
177 94 : GRB_TRY (GrB_assign (pNonstars, stars, 0, n, GrB_ALL, 0, 0));
178 94 : GRB_TRY (GrB_mxv (hookMNP, stars, 0, GrB_MIN_SECOND_SEMIRING_UINT64,
179 : S, pNonstars, 0));
180 : // select the valid elemenets (<n) of hookMNP
181 94 : GRB_TRY (GrB_assign (pNonstars, 0, 0, n, GrB_ALL, 0, 0));
182 94 : GRB_TRY (GrB_eWiseMult (mask, 0, 0, GrB_LT_UINT64, hookMNP, pNonstars, 0));
183 94 : GRB_TRY (GrB_eWiseMult (hookP, mask, 0, GrB_SECOND_UINT64, hookMNP, parents, 0));
184 94 : GRB_TRY (GrB_Vector_nvals (&nHooks, hookP));
185 94 : GRB_TRY (GrB_Vector_extractTuples (I, V, &nHooks, hookP));
186 94 : GRB_TRY (GrB_Vector_new (&tmp, GrB_UINT64, nHooks));
187 94 : GRB_TRY (GrB_extract (tmp, 0, 0, hookMNP, I, nHooks, 0));
188 94 : GRB_TRY (GrB_assign (parents, 0, 0, n, V, nHooks, 0)); // !!
189 94 : LG_TRY (Reduce_assign (parents, tmp, V, nHooks, msg));
190 : // modify the star vector
191 94 : GRB_TRY (GrB_assign (stars, 0, 0, false, V, nHooks, 0));
192 94 : GRB_TRY (GrB_Vector_extractTuples (I, V, &n, parents));
193 94 : GRB_TRY (GrB_extract (mask, 0, 0, stars, V, n, 0));
194 94 : GRB_TRY (GrB_assign (stars, 0, GrB_LAND, mask, GrB_ALL, 0, 0));
195 : // check termination
196 94 : GRB_TRY (GrB_reduce (&nStars, 0, GrB_PLUS_MONOID_UINT64, stars, 0));
197 94 : if (nStars == n) break;
198 : // clean up
199 70 : GRB_TRY (GrB_Vector_clear(hookMNP));
200 70 : GRB_TRY (GrB_Vector_clear(hookP));
201 70 : GRB_TRY (GrB_Vector_clear(pNonstars));
202 70 : GRB_TRY (GrB_free (&tmp));
203 : // ---------------------------------------------------------
204 : // Shortcut(parents);
205 : // ---------------------------------------------------------
206 70 : GRB_TRY (GrB_Vector_extractTuples (I, V, &n, parents));
207 70 : GRB_TRY (GrB_extract (gp, 0, 0, parents, V, n, 0));
208 70 : GRB_TRY (GrB_assign (parents, 0, 0, gp, GrB_ALL, 0, 0));
209 : // ---------------------------------------------------------
210 : // StarCheck(parents, stars);
211 : // ---------------------------------------------------------
212 : // calculate grandparents
213 70 : GRB_TRY (GrB_Vector_extractTuples (I, V, &n, parents));
214 70 : GRB_TRY (GrB_extract (gp, 0, 0, parents, V, n, 0));
215 : // identify vertices whose parent and grandparent are different
216 70 : GRB_TRY (GrB_eWiseMult (mask, 0, 0, GrB_NE_UINT64, gp, parents, 0));
217 70 : GRB_TRY (GrB_Vector_new (&nsgp, GrB_UINT64, n));
218 70 : GRB_TRY (GrB_assign (nsgp, mask, 0, gp, GrB_ALL, 0, 0));
219 : // extract indices and values for assign
220 70 : GRB_TRY (GrB_Vector_nvals (&nNonstars, nsgp));
221 70 : GRB_TRY (GrB_Vector_extractTuples (I, V, &nNonstars, nsgp));
222 70 : GRB_TRY (GrB_free (&nsgp));
223 70 : GRB_TRY (GrB_assign (stars, 0, 0, true, GrB_ALL, 0, 0));
224 70 : GRB_TRY (GrB_assign (stars, 0, 0, false, I, nNonstars, 0));
225 70 : GRB_TRY (GrB_assign (stars, 0, 0, false, V, nNonstars, 0));
226 : // extract indices and values for assign
227 70 : GRB_TRY (GrB_Vector_extractTuples (I, V, &n, parents));
228 70 : GRB_TRY (GrB_extract (mask, 0, 0, stars, V, n, 0));
229 70 : GRB_TRY (GrB_assign (stars, 0, GrB_LAND, mask, GrB_ALL, 0, 0));
230 : }
231 24 : *result = parents;
232 24 : parents = NULL ; // return parents (set to NULL so it isn't freed)
233 :
234 24 : LG_FREE_ALL;
235 24 : return GrB_SUCCESS;
236 : }
|