-
Notifications
You must be signed in to change notification settings - Fork 20
/
hopcroft-karp.lisp
247 lines (236 loc) · 10.9 KB
/
hopcroft-karp.lisp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
(defpackage :cp/hopcroft-karp
(:use :cl)
(:export #:bipartite-graph #:make-bgraph #:bipartite-graph-p #:coerce-to-bgraph
#:bgraph-matching1 #:bgraph-matching2 #:bgraph-size1 #:bgraph-size2 #:bgraph-graph1
#:bgraph-add-edge! #:bgraph-build-matching! #:bgraph-decompose)
(:documentation "Provides Hopcroft-Karp algorithm for maximum bipartite matching.
Time complexity: O(E sqrt(V))"))
(in-package :cp/hopcroft-karp)
(defconstant +graph-inf-distance+ #xffffffff)
(defstruct (bipartite-graph
(:constructor make-bgraph
(size1
size2
&aux
(graph1 (make-array size1 :element-type 'list :initial-element nil))
(matching1 (make-array size1 :element-type 'fixnum :initial-element -1))
(matching2 (make-array size2 :element-type 'fixnum :initial-element -1))))
(:conc-name bgraph-)
(:copier nil))
(size1 0 :type (unsigned-byte 32))
(size2 0 :type (unsigned-byte 32))
(graph1 nil :type (simple-array list (*)))
(matching1 nil :type (simple-array fixnum (*)))
(matching2 nil :type (simple-array fixnum (*))))
(declaim (inline bgraph-add-edge!))
(defun bgraph-add-edge! (bgraph vertex1 vertex2)
(push vertex2 (aref (bgraph-graph1 bgraph) vertex1))
bgraph)
(defun %fill-levels (bgraph levels1 levels2 queue)
"Does BFS and fills LEVELS."
(declare (optimize (speed 3) (safety 0))
((simple-array (unsigned-byte 32) (*)) levels1 levels2 queue))
(let ((graph1 (bgraph-graph1 bgraph))
(matching1 (bgraph-matching1 bgraph))
(matching2 (bgraph-matching2 bgraph))
(q-front 0)
(q-end 0)
(found nil))
(declare ((mod #.array-dimension-limit) q-front q-end))
(labels ((enqueue (obj)
(setf (aref queue q-end) obj)
(incf q-end))
(dequeue ()
(prog1 (aref queue q-front)
(incf q-front))))
(declare (inline enqueue dequeue))
(fill levels1 +graph-inf-distance+)
(fill levels2 +graph-inf-distance+)
(dotimes (i (bgraph-size1 bgraph))
(when (= -1 (aref matching1 i))
(setf (aref levels1 i) 0)
(enqueue i)))
(loop until (= q-front q-end)
for vertex = (dequeue)
do (dolist (next (aref graph1 vertex))
(when (= +graph-inf-distance+ (aref levels2 next))
(setf (aref levels2 next) (+ 1 (aref levels1 vertex)))
(let ((partner (aref matching2 next)))
(when (= -1 partner)
(setq found t)
(return))
(setf (aref levels1 partner) (+ 1 (aref levels2 next)))
(enqueue partner))))))
found))
(defun %find-matching (bgraph src levels1 levels2)
"Does DFS and makes matching greedily on the residual network."
(declare (optimize (speed 3) (safety 0))
((mod #.array-dimension-limit) src)
((simple-array (unsigned-byte 32) (*)) levels1 levels2))
(let ((matching1 (bgraph-matching1 bgraph))
(matching2 (bgraph-matching2 bgraph))
(graph1 (bgraph-graph1 bgraph)))
(labels ((dfs (v)
(declare ((mod #.array-dimension-limit) v))
(dolist (next (aref graph1 v))
(when (= (aref levels2 next) (+ 1 (aref levels1 v)))
(setf (aref levels2 next) +graph-inf-distance+) ; mark visited
(let ((partner (aref matching2 next)))
(when (or (= -1 partner) (dfs partner))
(setf (aref matching1 v) next
(aref matching2 next) v
(aref levels1 v) +graph-inf-distance+ ; mark visited
)
(return-from dfs t)))))
(setf (aref levels1 v) +graph-inf-distance+) ; mark visited
nil ; not matched
))
(dfs src))))
(declaim (ftype (function * (values (unsigned-byte 32) &optional)) bgraph-build-matching!))
(defun bgraph-build-matching! (bgraph)
"Makes a maximum bipartite matching and returns two vectors: correspondence
from group 1 to group 2, and correspondence from group 2 to group 1. At an
unmatched vertex, -1 is stored.
NOTE: Pay attention to the stack size!"
(declare (optimize (speed 3)))
(let* ((size1 (bgraph-size1 bgraph))
(size2 (bgraph-size2 bgraph))
(matching1 (bgraph-matching1 bgraph))
(levels1 (make-array size1 :element-type '(unsigned-byte 32)))
(levels2 (make-array size2 :element-type '(unsigned-byte 32)))
(queue (make-array (+ size1 size2) :element-type '(unsigned-byte 32)))
(count 0))
(declare ((mod #.array-dimension-limit) count))
(loop while (%fill-levels bgraph levels1 levels2 queue)
do (dotimes (v size1)
(when (and (= -1 (aref matching1 v))
(%find-matching bgraph v levels1 levels2))
(incf count))))
count))
;; not tested
(defun coerce-to-bgraph (graph)
"Converts adjacency lists representation of undirected graph to
BIPARTITE-GRAPH.
GRAPH := vector of adjacency lists"
(declare (vector graph))
(let* ((n (length graph))
(visited (make-array n :element-type 'bit :initial-element 0))
(colors (make-array n :element-type 'bit :initial-element 0))
(nums (make-array n :element-type 'fixnum))
(size0 0)
(size1 0))
(declare ((mod #.array-dimension-limit) size0 size1))
(labels ((dfs (vertex color)
(cond ((zerop (aref visited vertex))
(setf (aref visited vertex) 1
(aref colors vertex) color)
(if (zerop color)
(setf (aref nums vertex) size0
size0 (+ size0 1))
(setf (aref nums vertex) size1
size1 (+ size1 1)))
(if (= color 1)
(dolist (neighbor (aref graph vertex))
(dfs neighbor 0))
(dolist (neighbor (aref graph vertex))
(dfs neighbor 1))))
((/= color (aref colors vertex))
(error "Not bipartite.")))))
(dotimes (i n)
(when (zerop (aref visited i))
(dfs i 1)))
(let ((bgraph (make-bgraph size0 size1)))
(dotimes (i n)
(when (zerop (aref colors i))
(let ((i-num (aref nums i)))
(dolist (j (aref graph i))
(let ((j-num (aref nums j)))
(bgraph-add-edge! bgraph i-num j-num))))))
bgraph))))
;; not tested
(declaim (ftype (function * (values (simple-array (integer 0 #.most-positive-fixnum) (*))
(simple-array (integer 0 #.most-positive-fixnum) (*))
&optional))
bgraph-decompose!))
(defun bgraph-decompose (bgraph)
"Decomposes a residual network to strongly connected components by Tarjan's
algorithm. BGRAPH-BUILD-MATCHING! must be called beforehand."
(declare (optimize (speed 3)))
(let* ((size1 (bgraph-size1 bgraph))
(size2 (bgraph-size2 bgraph))
(total-size (+ size1 size2 2))
(source (+ size1 size2))
(sink (+ size1 size2 1))
(graph1 (bgraph-graph1 bgraph))
(matching1 (bgraph-matching1 bgraph))
(matching2 (bgraph-matching2 bgraph))
(ord 0) ; in-order
(ords (make-array total-size :element-type 'fixnum :initial-element -1))
(lowlinks (make-array total-size :element-type 'fixnum))
(components (make-array total-size
:element-type '(integer 0 #.most-positive-fixnum)))
(comp-index 0) ; index number of component
(sizes (make-array total-size
:element-type '(integer 0 #.most-positive-fixnum)
:initial-element 0))
(stack (make-array total-size :element-type '(integer 0 #.most-positive-fixnum)))
(end 0) ; stack pointer
(in-stack (make-array total-size :element-type 'bit :initial-element 0)))
(declare ((mod #.array-dimension-limit) ord end comp-index source sink total-size))
(labels ((%push (v)
(setf (aref stack end) v
(aref in-stack v) 1)
(incf end))
(%pop ()
(decf end)
(let ((v (aref stack end)))
(setf (aref in-stack v) 0)
v))
(frob (v next)
(cond ((= -1 (aref ords next))
(visit next)
(setf (aref lowlinks v)
(min (aref lowlinks v) (aref lowlinks next))))
((= 1 (aref in-stack next))
(setf (aref lowlinks v)
(min (aref lowlinks v) (aref ords next))))))
(visit (v)
(setf (aref ords v) ord
(aref lowlinks v) ord)
(incf ord)
(%push v)
(cond ((= v source)
(loop for next below size1
when (= -1 (aref matching1 next))
do (frob v next)))
((= v sink)
(loop for next below size2
unless (= -1 (aref matching2 next))
do (frob v (+ size1 next))))
((and (< v size1) (= -1 (aref matching1 v)))
(dolist (next (aref graph1 v))
(declare ((mod #.array-dimension-limit) next))
(frob v (+ next size1))))
((and (< v size1) (/= -1 (aref matching1 v)))
(frob v source)
(dolist (next (aref graph1 v))
(declare ((mod #.array-dimension-limit) next))
(unless (= next (aref matching1 v))
(frob v (+ next size1)))))
((= -1 (aref matching2 (- v size1)))
(frob v sink))
(t
;; (assert (/= -1 (aref matching2 (- v size1))))
(frob v (aref matching2 (- v size1)))))
(when (= (aref lowlinks v) (aref ords v))
(loop for size of-type (mod #.array-dimension-limit) from 1
for w = (%pop)
do (setf (aref components w) comp-index)
until (= v w)
finally (setf (aref sizes comp-index) size)
(incf comp-index)))))
(dotimes (v total-size)
(when (= -1 (aref ords v))
(visit v)))
(values (subseq components 0 size1)
(subseq components size1 (+ size1 size2))))))