Line data Source code
1 : //------------------------------------------------------------------------------
2 : // LAGraph_MaximalMatching: maximal matching using an adaptation of Luby's MIS algorithm
3 : // on a line graph. Derived from LAGraph_MaximalIndependentSet
4 : //------------------------------------------------------------------------------
5 :
6 : // LAGraph, (c) 2022 by The LAGraph Contributors, All Rights Reserved.
7 : // SPDX-License-Identifier: BSD-2-Clause
8 :
9 : // Contributed by Vidith Madhu, Texas A&M University
10 :
11 : //------------------------------------------------------------------------------
12 :
13 : // TODO: ready for src? need a vanilla non-GxB, and incidence graphs.
14 :
15 : /*
16 : Uses a modified version of Luby's MIS algorithm
17 : Major algorithm steps:
18 : - Compute score for each edge
19 : - Find max score neighbor of each edge (*)
20 : - Retain edges with score == max score neighbor (*)
21 : - Add retained edges to result
22 : - Remove retained edges and their neighbors from the graph (*)
23 :
24 : (*): these steps involve what can be thought as a "2-hop" process that involves two
25 : GrB_mxv's: the first to go from edges to vertices, and the second from vertices back to edges.
26 : Tying both steps together yields a single BFS-like step in the line graph. A important side effect
27 : of this is that the source edge gets included in the result of this 2-hop step, which cannot be avoided
28 : since we do not compute E'E explicitly.
29 :
30 : The input to this method is an incidence matrix E, of size n-by-e where the
31 : undirected graph G has n nodes and e edges. If the kth edge of G is the edge
32 : (i,j), then the column E(:,k) contains two entries: E(i,k) and E(j,k), which
33 : have the same value. If the graph G is weighted, then both E(i,k) and E(j,k)
34 : are equal to the weight of the (i,j) edge. If G is unweighted, then both are
35 : equal to 1 (and the matrix E is thus iso-valued).
36 :
37 : The output is vector 'matching' of size e, where matching(k) is present (and
38 : equal to true) if the kth edge appears in the maximal matching. If (i,j) is
39 : a matched edge, then no other edges of G that are incident on nodes i and j
40 : appear in the matching.
41 :
42 : This method requires O(e) space for an undirected graph with e edges
43 : */
44 :
45 : #include "LG_internal.h"
46 : #include "LAGraphX.h"
47 :
48 : // #define dbg
49 :
50 : #undef LG_FREE_ALL
51 : #undef LG_FREE_WORK
52 :
53 : #define OPTIMIZE_PUSH_PULL
54 :
55 : #define LG_FREE_WORK \
56 : { \
57 : GrB_free(&score) ; \
58 : GrB_free(&candidates) ; \
59 : GrB_free(&Seed) ; \
60 : GrB_free(&node_degree) ; \
61 : GrB_free(°ree) ; \
62 : GrB_free(&max_node_neighbor) ; \
63 : GrB_free(&max_neighbor) ; \
64 : GrB_free(&new_members) ; \
65 : GrB_free(&new_neighbors) ; \
66 : GrB_free(&new_members_nodes) ; \
67 : GrB_free(&new_members_node_degree) ; \
68 : GrB_free(&empty) ; \
69 : GrB_free(&weight) ; \
70 : } \
71 :
72 : #define LG_FREE_ALL \
73 : { \
74 : LG_FREE_WORK ; \
75 : GrB_free(&result) ; \
76 : } \
77 :
78 : #define MAX_FAILURES 50
79 :
80 213 : int LAGraph_MaximalMatching
81 : (
82 : // outputs:
83 : GrB_Vector *matching, // pointer to output vector
84 : // inputs:
85 : GrB_Matrix E, // incidence
86 : GrB_Matrix E_t, // incidence transposed (if incorrect, results are undefined)
87 : LAGraph_Matching_kind matching_type, // 0 (random), 1 (heavy weight matching), 2 (light weight matching)
88 : uint64_t seed, // random number seed
89 : char *msg
90 : )
91 : {
92 : #if LAGRAPH_SUITESPARSE
93 213 : LG_CLEAR_MSG ;
94 :
95 213 : if ((matching == NULL) || (E == NULL) || (E_t == NULL)) {
96 3 : return GrB_NULL_POINTER ;
97 : }
98 :
99 210 : GrB_Vector score = NULL ; // score for each edge. Computed according to matching_type
100 210 : GrB_Vector weight = NULL ; // weight of each edge
101 210 : GrB_Vector candidates = NULL ; // set of candidate edges
102 210 : GrB_Vector Seed = NULL ; // random number seed vector
103 210 : GrB_Vector node_degree = NULL ; // intermediate result for computing edge degree; degree of the node. Only computed once
104 210 : GrB_Vector degree = NULL ; // edge degree; number of incident edges. Only computed once
105 210 : GrB_Vector max_node_neighbor = NULL ; // intermediate result for computing max edge neighbor; max edge touching a node
106 210 : GrB_Vector max_neighbor = NULL ; // max neighbor of an edge (including itself)
107 210 : GrB_Vector new_members = NULL ; // new edges to include in the matching
108 210 : GrB_Vector new_neighbors = NULL ; // union of new members and their neighbor edges
109 210 : GrB_Vector new_members_nodes = NULL ; // nodes touching an edge in new_members
110 210 : GrB_Vector new_members_node_degree = NULL ; // node degree considering only new members
111 210 : GrB_Vector result = NULL ; // resulting matching
112 210 : GrB_Vector empty = NULL ; // empty vector
113 :
114 : GrB_Index num_edges ;
115 : GrB_Index num_nodes ;
116 :
117 : char typename[LAGRAPH_MAX_NAME_LEN] ;
118 : GrB_Type type ;
119 210 : LG_TRY (LAGraph_Matrix_TypeName (typename, E, msg)) ;
120 210 : LG_TRY (LAGraph_TypeFromName (&type, typename, msg)) ;
121 :
122 :
123 210 : GRB_TRY (GrB_Matrix_nrows (&num_nodes, E)) ;
124 210 : GRB_TRY (GrB_Matrix_ncols (&num_edges, E)) ;
125 :
126 210 : GRB_TRY (GrB_Vector_new (&candidates, GrB_BOOL, num_edges)) ;
127 210 : GRB_TRY (GrB_Vector_new (&Seed, GrB_UINT64, num_edges)) ;
128 210 : GRB_TRY (GrB_Vector_new (&score, GrB_FP64, num_edges)) ;
129 210 : GRB_TRY (GrB_Vector_new (&weight, type, num_edges)) ;
130 210 : GRB_TRY (GrB_Vector_new (&node_degree, GrB_UINT64, num_nodes)) ;
131 210 : GRB_TRY (GrB_Vector_new (°ree, GrB_UINT64, num_edges)) ;
132 210 : GRB_TRY (GrB_Vector_new (&max_node_neighbor, GrB_FP64, num_nodes)) ;
133 210 : GRB_TRY (GrB_Vector_new (&max_neighbor, GrB_FP64, num_edges)) ;
134 210 : GRB_TRY (GrB_Vector_new (&new_members, GrB_BOOL, num_edges)) ;
135 210 : GRB_TRY (GrB_Vector_new (&new_neighbors, GrB_BOOL, num_edges)) ;
136 210 : GRB_TRY (GrB_Vector_new (&new_members_nodes, GrB_BOOL, num_nodes)) ;
137 210 : GRB_TRY (GrB_Vector_new (&new_members_node_degree, GrB_UINT64, num_nodes)) ;
138 210 : GRB_TRY (GrB_Vector_new (&result, GrB_BOOL, num_edges)) ;
139 210 : GRB_TRY (GrB_Vector_new (&empty, GrB_BOOL, num_edges)) ;
140 :
141 210 : GRB_TRY (GrB_assign (Seed, NULL, NULL, 0, GrB_ALL, num_edges, NULL)) ;
142 :
143 210 : LG_TRY (LAGraph_Random_Seed (Seed, seed, msg)) ;
144 :
145 : // initially all edges are considered
146 210 : GRB_TRY (GrB_assign (candidates, NULL, NULL, true, GrB_ALL,
147 : num_edges, NULL)) ;
148 :
149 : GrB_Index ncandidates ;
150 210 : GrB_Index nfailures = 0 ; // counts how many iterations have failed due to invalid matchings
151 :
152 210 : GRB_TRY (GrB_Vector_nvals(&ncandidates, candidates)) ;
153 :
154 : // for each node, counts incident edges
155 210 : GRB_TRY (GrB_mxv (node_degree, NULL, NULL, LAGraph_plus_one_uint64, E, candidates, NULL)) ;
156 :
157 : // for each edge, sums incident edges for each node. Each edge has an excess of 2 degree, but it doesn't matter since
158 : // we care about relative degree
159 210 : GRB_TRY (GrB_mxv (degree, NULL, NULL, LAGraph_plus_second_uint64, E_t, node_degree, NULL)) ;
160 :
161 210 : GRB_TRY (GrB_reduce (weight, NULL, NULL, GrB_MAX_MONOID_FP64, E_t, NULL)) ; // use ANY ?
162 :
163 210 : double sparsity_thresh =
164 : #ifdef OPTIMIZE_PUSH_PULL
165 : 0.04 ;
166 : #else
167 : 1.0;
168 : #endif
169 :
170 : #if defined ( COVERAGE )
171 210 : int kount = 0 ;
172 : #endif
173 :
174 1161 : while (ncandidates > 0) {
175 : // first just generate the scores again
176 951 : GRB_TRY (GrB_eWiseMult (score, candidates, NULL, GrB_DIV_FP64, Seed, degree, GrB_DESC_RS)) ;
177 :
178 : // for light matching, can multiply scores by 1 / (edge weight)
179 951 : if (matching_type == LAGraph_Matching_heavy) {
180 : // heavy
181 225 : GRB_TRY (GrB_eWiseMult (score, NULL, NULL, GrB_TIMES_FP64, score, weight, NULL)) ;
182 726 : } else if (matching_type == LAGraph_Matching_light) {
183 : // light
184 217 : GRB_TRY (GrB_eWiseMult (score, NULL, NULL, GrB_DIV_FP64, score, weight, NULL)) ;
185 : }
186 :
187 : // the actual edge selection is common regardless of matching type
188 :
189 : // intermediate result. Max score edge touching each node
190 : // don't need to clear this out first because we populate the result for all nodes
191 951 : if (ncandidates > sparsity_thresh * num_edges) {
192 586 : GRB_TRY (LG_SET_FORMAT_HINT (score, LG_BITMAP)) ;
193 586 : GRB_TRY (GrB_mxv (max_node_neighbor, NULL, NULL, GrB_MAX_SECOND_SEMIRING_FP64, E, score, NULL)) ;
194 : } else {
195 365 : GRB_TRY (LG_SET_FORMAT_HINT (score, LG_SPARSE)) ;
196 365 : GRB_TRY (GrB_vxm (max_node_neighbor, NULL, NULL, GrB_MAX_FIRST_SEMIRING_FP64, score, E_t, NULL)) ;
197 : }
198 :
199 : GrB_Index node_nvals ;
200 951 : GRB_TRY (GrB_Vector_nvals (&node_nvals, max_node_neighbor)) ;
201 :
202 : // Max edge touching each candidate edge, including itself
203 951 : if (node_nvals > sparsity_thresh * num_nodes) {
204 760 : GRB_TRY (LG_SET_FORMAT_HINT (max_node_neighbor, LG_BITMAP)) ;
205 760 : GRB_TRY (GrB_mxv (max_neighbor, candidates, NULL, GrB_MAX_SECOND_SEMIRING_FP64, E_t, max_node_neighbor, GrB_DESC_RS)) ;
206 : } else {
207 191 : GRB_TRY (LG_SET_FORMAT_HINT (max_node_neighbor, LG_SPARSE)) ;
208 191 : GRB_TRY (GrB_vxm (max_neighbor, candidates, NULL, GrB_MAX_FIRST_SEMIRING_FP64, max_node_neighbor, E, GrB_DESC_RS)) ;
209 : }
210 : // Note that we are using the GE operator and not G, since max_neighbor includes the self score
211 : // correctness: both score and max_neighbor only have entries for candidates, so no non-candidate members are produced
212 : // GRB_TRY (GrB_assign (new_members, NULL, NULL, empty, GrB_ALL, num_edges, NULL)) ; // just experimenting
213 951 : GRB_TRY (GrB_eWiseAdd (new_members, NULL, NULL, GrB_GE_FP64, score, max_neighbor, NULL)) ;
214 :
215 : // makes new_members structural
216 951 : GRB_TRY (GrB_select (new_members, NULL, NULL, GrB_VALUEEQ_BOOL, new_members, true, NULL)) ;
217 : #ifdef dbg
218 : printf("new members for ncandidates = %lld:\n", ncandidates);
219 : LAGRAPH_TRY (LAGraph_Vector_Print (new_members, LAGraph_SHORT, stdout, msg)) ;
220 : #endif
221 :
222 : GrB_Index new_members_nvals ;
223 951 : GRB_TRY (GrB_Vector_nvals (&new_members_nvals, new_members)) ;
224 :
225 : // check if any node has > 1 edge touching it.
226 951 : if (new_members_nvals > sparsity_thresh * num_edges) {
227 271 : GRB_TRY (LG_SET_FORMAT_HINT (new_members, LG_BITMAP)) ;
228 271 : GRB_TRY (GrB_mxv (new_members_node_degree, NULL, NULL, LAGraph_plus_one_uint64, E, new_members, NULL)) ;
229 : } else {
230 680 : GRB_TRY (LG_SET_FORMAT_HINT (new_members, LG_SPARSE)) ;
231 680 : GRB_TRY (GrB_vxm (new_members_node_degree, NULL, NULL, LAGraph_plus_one_uint64, new_members, E_t, NULL)) ;
232 : }
233 :
234 : GrB_Index max_degree ;
235 951 : GRB_TRY (GrB_reduce (&max_degree, NULL, GrB_MAX_MONOID_UINT64, new_members_node_degree, NULL)) ;
236 :
237 : #if defined ( COVERAGE )
238 951 : if (num_nodes == 20 && kount++ == 1) max_degree = 2 ;
239 951 : if (num_nodes == 30 && kount++ == 0) max_degree = 2 ;
240 : #endif
241 :
242 951 : if (max_degree > 1) {
243 10 : nfailures++ ;
244 10 : LG_ASSERT_MSG (nfailures <= MAX_FAILURES, LAGRAPH_CONVERGENCE_FAILURE,
245 : "method has stalled") ;
246 : // regen seed and seed vector
247 10 : LG_TRY (LAGraph_Random_Seed (Seed, seed + nfailures, msg)) ;
248 10 : continue ;
249 : }
250 : // add new members to result and remove from candidates
251 : // also want to remove all adjacent edges in new_members from candidates
252 941 : GRB_TRY (GrB_assign (result, new_members, NULL, true, GrB_ALL, num_edges, GrB_DESC_S)) ;
253 : // to include neighbor edges, need to compute new_neighbors
254 : // to do this, we need to compute the intermediate result new_members_nodes
255 941 : if (new_members_nvals > sparsity_thresh * num_edges) {
256 261 : GRB_TRY (LG_SET_FORMAT_HINT (new_members, LG_BITMAP)) ;
257 261 : GRB_TRY (GrB_mxv (new_members_nodes, NULL, NULL, LAGraph_any_one_bool, E, new_members, NULL)) ;
258 : } else {
259 680 : GRB_TRY (LG_SET_FORMAT_HINT (new_members, LG_SPARSE)) ;
260 680 : GRB_TRY (GrB_vxm (new_members_nodes, NULL, NULL, LAGraph_any_one_bool, new_members, E_t, NULL)) ;
261 : }
262 :
263 941 : GRB_TRY (GrB_Vector_nvals (&node_nvals, new_members_nodes)) ;
264 :
265 941 : if (node_nvals > sparsity_thresh * num_nodes) {
266 682 : GRB_TRY (LG_SET_FORMAT_HINT (new_members_nodes, LG_BITMAP)) ;
267 682 : GRB_TRY (GrB_mxv (new_neighbors, NULL, NULL, LAGraph_any_one_bool, E_t, new_members_nodes, NULL)) ;
268 : } else {
269 259 : GRB_TRY (LG_SET_FORMAT_HINT (new_members_nodes, LG_SPARSE)) ;
270 259 : GRB_TRY (GrB_vxm (new_neighbors, NULL, NULL, LAGraph_any_one_bool, new_members_nodes, E, NULL)) ;
271 : }
272 :
273 : #ifdef dbg
274 : LAGRAPH_TRY (LAGraph_Vector_Print (new_neighbors, LAGraph_SHORT, stdout, msg)) ;
275 : #endif
276 : // removes the union of new_members and their neighbors
277 941 : GRB_TRY (GrB_assign (candidates, new_neighbors, NULL, empty, GrB_ALL, num_edges, GrB_DESC_S)) ;
278 :
279 : #ifdef dbg
280 : printf("candidates:\n");
281 : LAGRAPH_TRY (LAGraph_Vector_Print (candidates, LAGraph_SHORT, stdout, msg)) ;
282 : #endif
283 941 : GrB_Index last_ncandidates = ncandidates ;
284 :
285 941 : GrB_Vector_nvals(&ncandidates, candidates) ;
286 :
287 : // advance seed vector
288 941 : LG_TRY (LAGraph_Random_Next (Seed, msg)) ;
289 :
290 : #if defined ( COVERAGE )
291 941 : if (num_nodes == 50 && kount++ == 0)
292 : {
293 : // hack the Seed vector
294 20 : GRB_TRY (GrB_assign (Seed, NULL, NULL, 42, GrB_ALL, num_edges, NULL)) ;
295 : }
296 : #endif
297 : }
298 :
299 210 : (*matching) = result ;
300 :
301 210 : LG_FREE_WORK ;
302 210 : return (GrB_SUCCESS) ;
303 : #else
304 : return (GrB_NOT_IMPLEMENTED) ;
305 : #endif
306 : }
|
