1 /* asnhall.c (find bipartite matching of maximum cardinality) */
2
3 /***********************************************************************
4 * This code is part of GLPK (GNU Linear Programming Kit).
5 * Copyright (C) 2009-2016 Free Software Foundation, Inc.
6 * Written by Andrew Makhorin <mao@gnu.org>.
7 *
8 * GLPK is free software: you can redistribute it and/or modify it
9 * under the terms of the GNU General Public License as published by
10 * the Free Software Foundation, either version 3 of the License, or
11 * (at your option) any later version.
12 *
13 * GLPK is distributed in the hope that it will be useful, but WITHOUT
14 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15 * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
16 * License for more details.
17 *
18 * You should have received a copy of the GNU General Public License
19 * along with GLPK. If not, see <http://www.gnu.org/licenses/>.
20 ***********************************************************************/
21
22 #include "env.h"
23 #include "glpk.h"
24 #include "mc21a.h"
25
26 /***********************************************************************
27 * NAME
28 *
29 * glp_asnprob_hall - find bipartite matching of maximum cardinality
30 *
31 * SYNOPSIS
32 *
33 * int glp_asnprob_hall(glp_graph *G, int v_set, int a_x);
34 *
35 * DESCRIPTION
36 *
37 * The routine glp_asnprob_hall finds a matching of maximal cardinality
38 * in the specified bipartite graph G. It uses a version of the Fortran
39 * routine MC21A developed by I.S.Duff [1], which implements Hall's
40 * algorithm [2].
41 *
42 * RETURNS
43 *
44 * The routine glp_asnprob_hall returns the cardinality of the matching
45 * found. However, if the specified graph is incorrect (as detected by
46 * the routine glp_check_asnprob), the routine returns negative value.
47 *
48 * REFERENCES
49 *
50 * 1. I.S.Duff, Algorithm 575: Permutations for zero-free diagonal, ACM
51 * Trans. on Math. Softw. 7 (1981), 387-390.
52 *
53 * 2. M.Hall, "An Algorithm for distinct representatives," Amer. Math.
54 * Monthly 63 (1956), 716-717. */
55
glp_asnprob_hall(glp_graph * G,int v_set,int a_x)56 int glp_asnprob_hall(glp_graph *G, int v_set, int a_x)
57 { glp_vertex *v;
58 glp_arc *a;
59 int card, i, k, loc, n, n1, n2, xij;
60 int *num, *icn, *ip, *lenr, *iperm, *pr, *arp, *cv, *out;
61 if (v_set >= 0 && v_set > G->v_size - (int)sizeof(int))
62 xerror("glp_asnprob_hall: v_set = %d; invalid offset\n",
63 v_set);
64 if (a_x >= 0 && a_x > G->a_size - (int)sizeof(int))
65 xerror("glp_asnprob_hall: a_x = %d; invalid offset\n", a_x);
66 if (glp_check_asnprob(G, v_set))
67 return -1;
68 /* determine the number of vertices in sets R and S and renumber
69 vertices in S which correspond to columns of the matrix; skip
70 all isolated vertices */
71 num = xcalloc(1+G->nv, sizeof(int));
72 n1 = n2 = 0;
73 for (i = 1; i <= G->nv; i++)
74 { v = G->v[i];
75 if (v->in == NULL && v->out != NULL)
76 n1++, num[i] = 0; /* vertex in R */
77 else if (v->in != NULL && v->out == NULL)
78 n2++, num[i] = n2; /* vertex in S */
79 else
80 { xassert(v->in == NULL && v->out == NULL);
81 num[i] = -1; /* isolated vertex */
82 }
83 }
84 /* the matrix must be square, thus, if it has more columns than
85 rows, extra rows will be just empty, and vice versa */
86 n = (n1 >= n2 ? n1 : n2);
87 /* allocate working arrays */
88 icn = xcalloc(1+G->na, sizeof(int));
89 ip = xcalloc(1+n, sizeof(int));
90 lenr = xcalloc(1+n, sizeof(int));
91 iperm = xcalloc(1+n, sizeof(int));
92 pr = xcalloc(1+n, sizeof(int));
93 arp = xcalloc(1+n, sizeof(int));
94 cv = xcalloc(1+n, sizeof(int));
95 out = xcalloc(1+n, sizeof(int));
96 /* build the adjacency matrix of the bipartite graph in row-wise
97 format (rows are vertices in R, columns are vertices in S) */
98 k = 0, loc = 1;
99 for (i = 1; i <= G->nv; i++)
100 { if (num[i] != 0) continue;
101 /* vertex i in R */
102 ip[++k] = loc;
103 v = G->v[i];
104 for (a = v->out; a != NULL; a = a->t_next)
105 { xassert(num[a->head->i] != 0);
106 icn[loc++] = num[a->head->i];
107 }
108 lenr[k] = loc - ip[k];
109 }
110 xassert(loc-1 == G->na);
111 /* make all extra rows empty (all extra columns are empty due to
112 the row-wise format used) */
113 for (k++; k <= n; k++)
114 ip[k] = loc, lenr[k] = 0;
115 /* find a row permutation that maximizes the number of non-zeros
116 on the main diagonal */
117 card = mc21a(n, icn, ip, lenr, iperm, pr, arp, cv, out);
118 #if 1 /* 18/II-2010 */
119 /* FIXED: if card = n, arp remains clobbered on exit */
120 for (i = 1; i <= n; i++)
121 arp[i] = 0;
122 for (i = 1; i <= card; i++)
123 { k = iperm[i];
124 xassert(1 <= k && k <= n);
125 xassert(arp[k] == 0);
126 arp[k] = i;
127 }
128 #endif
129 /* store solution, if necessary */
130 if (a_x < 0) goto skip;
131 k = 0;
132 for (i = 1; i <= G->nv; i++)
133 { if (num[i] != 0) continue;
134 /* vertex i in R */
135 k++;
136 v = G->v[i];
137 for (a = v->out; a != NULL; a = a->t_next)
138 { /* arp[k] is the number of matched column or zero */
139 if (arp[k] == num[a->head->i])
140 { xassert(arp[k] != 0);
141 xij = 1;
142 }
143 else
144 xij = 0;
145 memcpy((char *)a->data + a_x, &xij, sizeof(int));
146 }
147 }
148 skip: /* free working arrays */
149 xfree(num);
150 xfree(icn);
151 xfree(ip);
152 xfree(lenr);
153 xfree(iperm);
154 xfree(pr);
155 xfree(arp);
156 xfree(cv);
157 xfree(out);
158 return card;
159 }
160
161 /* eof */
162