1#include "lb.h"
2#include "bm.h"
3#include <assert.h>
4#include <string.h>
5#ifdef LB_DEBUG
6#include <stdio.h>
7#endif
8
9void
10lb_init(struct lb * lb, size_t n) {
11 assert(n <= LB_SIZE);
12 assert(n > 0);
13 memset(lb, 0, sizeof(*lb));
14 lb->size = n;
15 lb->dirty = true;
16 for (size_t x = 0; x < lb->size; x++) {
17 lb->x_root[x] = x;
18 lb->x_has_equivs[x] = false;
19 BM_SET(lb->x_equivs[x], x);
20 }
21}
22
23void
24lb_init_equivalent(struct lb * lb, size_t x0, size_t x1) {
25 assert(x0 < x1);
26 assert(x1 < lb->size);
27 assert(lb->x_root[x0] == x0);
28 assert(!lb->x_has_equivs[x1]);
29 assert(!lb->marked);
30 if (lb->x_root[x1] == x0) return;
31 assert(lb->x_root[x1] == x1);
32
33 lb->x_root[x1] = x0;
34 BM_SET(lb->x_equivs[x0], x1);
35 lb->x_has_equivs[x0] = true;
36}
37
38static void
39lb_mark_positive_half(struct lb * lb, bool h, size_t x, size_t y) {
40 assert(lb->halfs[h].pair[x] == 0);
41 for (size_t i = 0; i < lb->size; i++) {
42 if (i != x) {
43 BM_SET(lb->halfs[h].bs[i], y);
44 }
45 if (i != y) {
46 BM_SET(lb->halfs[h].bs[x], i);
47 }
48 }
49 BM_SET(lb->halfs[h].paired, x);
50 lb->halfs[h].pair[x] = y;
51 lb->dirty = true;
52}
53
54void
55lb_mark_positive(struct lb * lb, size_t x, size_t y) {
56 assert(x < lb->size);
57 assert(y < lb->size);
58
59 if (lb->x_has_equivs[x]) {
60 size_t x_original = x;
61 while (BM_ISSET(lb->halfs[0].paired, x)) {
62 x++;
63 assert(x < lb->size && lb->x_root[x] == x_original);
64 }
65 }
66
67 assert(!BM_ISSET(lb->halfs[0].bs[x], y));
68 assert(!BM_ISSET(lb->halfs[1].bs[y], x));
69 if (BM_ISSET(lb->halfs[1].paired, y)) {
70 assert(BM_ISSET(lb->halfs[1].paired, y));
71 assert(lb->halfs[0].pair[x] == y);
72 assert(lb->halfs[1].pair[y] == x);
73 return;
74 }
75 assert(!BM_ISSET(lb->halfs[1].paired, y));
76
77 for (size_t i = lb->x_root[x]; i < lb->size; i++) {
78 if (!BM_ISSET(lb->x_equivs[x], i)) {
79 continue;
80 }
81 if (BM_ISSET(lb->halfs[0].paired, i)) {
82 continue;
83 }
84 lb_mark_positive_half(lb, 0, i, y);
85 lb_mark_positive_half(lb, 1, y, i);
86 return;
87 }
88 assert(0);
89}
90
91static void
92lb_mark_negative_half(struct lb * lb, bool h, size_t x, size_t y) {
93 if (BM_ISSET(lb->halfs[h].paired, x)) {
94 assert(lb->halfs[h].pair[x] != y);
95 }
96 if (!BM_ISSET(lb->halfs[h].bs[x], y)) {
97 BM_SET(lb->halfs[h].bs[x], y);
98 lb->dirty = true;
99 }
100}
101
102void
103lb_mark_negative(struct lb * lb, size_t x, size_t y) {
104 assert(x < lb->size);
105 assert(y < lb->size);
106 for (size_t i = lb->x_root[x]; i < lb->size; i++) {
107 if (!BM_ISSET(lb->x_equivs[x], i)) {
108 continue;
109 }
110 lb_mark_negative_half(lb, 0, i, y);
111 lb_mark_negative_half(lb, 1, y, i);
112 if (!lb->x_has_equivs[x]) {
113 break;
114 }
115 }
116}
117
118static void
119lb_deduce_unique(struct lb * lb, bool h) {
120 // Look for cases where there is only 1 possibility left
121 for (size_t i = 0; i < lb->size; i++) {
122 if (BM_ISSET(lb->halfs[h].paired, i)) {
123 continue;
124 }
125 size_t p = BM_POPCOUNT(lb->halfs[h].bs[i]);
126 assert(p < lb->size);
127 if (p == lb->size - 1) {
128 size_t z = BM_FFZ(lb->halfs[h].bs[i]);
129 assert(z < lb->size);
130 lb_mark_positive_half(lb, h, i, z);
131 lb_mark_positive_half(lb, !h, z, i);
132 }
133 }
134}
135
136static void
137lb_deduce_clique(struct lb * lb, bool h) {
138 // Look for self-contained subsets
139 size_t l = 0;
140 for (size_t i = 0; i < lb->size; i++) {
141 if (BM_ISSET(lb->halfs[h].paired, i)) {
142 continue;
143 }
144 size_t p = lb->size - BM_POPCOUNT(lb->halfs[h].bs[i]);
145 if (p <= lb->size / 2) {
146 uint64_t js[LB_BM_SIZE];
147 memset(js, 0xFF, sizeof(js));
148 for (size_t j = 0; j < lb->size; j++) {
149 if (!BM_ISSET(lb->halfs[h].bs[i], j)) {
150 continue;
151 }
152 BM_ANDEQ(js, lb->halfs[!h].bs[j]);
153 }
154 size_t q = BM_POPCOUNT(js);
155 assert(q <= p); // Hill's Criteria
156 if (q == p) {
157 // Found a clique
158 for (size_t x = 0; x < lb->size; x++) {
159 if (BM_ISSET(js, x)) {
160 continue;
161 }
162 for (size_t y = 0; y < lb->size; y++) {
163 if (BM_ISSET(lb->halfs[h].bs[i], y)) {
164 continue;
165 }
166 lb_mark_negative_half(lb, h, x, y);
167 lb_mark_negative_half(lb, !h, y, x);
168 }
169 }
170 return;
171 }
172 }
173 }
174}
175
176static void
177lb_deduce_equivalents(struct lb * lb) {
178 for (size_t x = 0; x < lb->size; x++) {
179 if (!lb->x_has_equivs[x]) {
180 continue;
181 }
182
183 uint64_t bs_union[LB_BM_SIZE];
184 size_t union_count = 0;
185 for (size_t i = lb->x_root[x]; i < lb->size; i++) {
186 if (!BM_ISSET(lb->x_equivs[x], i)) {
187 continue;
188 }
189 if (BM_ISSET(lb->halfs[0].paired, i)) {
190 continue;
191 }
192 if (union_count == 0) {
193 memcpy(bs_union, lb->halfs[0].bs[i], sizeof(bs_union));
194 } else {
195 BM_ANDEQ(bs_union, lb->halfs[0].bs[i]);
196 }
197 union_count++;
198 }
199 if (union_count < 2) {
200 // Would be handled by lb_deduce_unique
201 continue;
202 }
203 size_t c = lb->size - BM_POPCOUNT(bs_union);
204 assert(c >= union_count);
205 if (c == union_count) {
206 for (size_t i = lb->x_root[x]; i < lb->size; i++) {
207 if (BM_ISSET(bs_union, i)) {
208 continue;
209 }
210 if (BM_ISSET(lb->halfs[0].paired, i)) {
211 continue;
212 }
213 lb_mark_positive(lb, i, BM_FFZ(lb->halfs[0].bs[i]));
214 }
215
216 }
217 }
218}
219
220void
221lb_deduce(struct lb * lb) {
222 while (lb->dirty) {
223 lb->dirty = false;
224
225 lb_deduce_unique(lb, 0);
226 if (lb->dirty) continue;
227 lb_deduce_unique(lb, 1);
228 if (lb->dirty) continue;
229
230 lb_deduce_clique(lb, 0);
231 if (lb->dirty) continue;
232 lb_deduce_clique(lb, 1);
233 if (lb->dirty) continue;
234
235 lb_deduce_equivalents(lb);
236 if (lb->dirty) continue;
237 }
238}
239
240#ifdef LB_DEBUG
241static void
242lb_selftest_print(struct lb * lb, const char * name) {
243 if (name != NULL) {
244 printf(">>> %s\n", name);
245 lb_deduce(lb);
246 }
247
248 printf(" ");
249 for (size_t i = 0; i < lb->size; i++) {
250 printf("%c", 'A' + (char) (i % 26));
251 }
252 printf(" ");
253 for (size_t i = 0; i < lb->size; i++) {
254 printf("%c", '0' + (char) (i % 10));
255 }
256 printf("\n");
257
258 for (size_t i = 0; i < lb->size; i++) {
259 for (int h = 0; h < 2; h++) {
260 if (h == 0) {
261 printf(" %c ", '0' + (char) (i % 10));
262 } else {
263 printf(" %c ", 'A' + (char) (i % 26));
264 }
265 for (size_t j = 0; j < lb->size; j++) {
266 if (BM_ISSET(lb->halfs[h].bs[i], j)) {
267 printf("x");
268 } else if (BM_ISSET(lb->halfs[h].paired, i)) {
269 assert(lb->halfs[h].pair[i] == j);
270 printf("@");
271 } else {
272 printf(".");
273 }
274 }
275 }
276 printf("\n");
277 }
278 printf("\n");
279}
280
281void
282lb_selftest() {
283 struct lb lb[1];
284 lb_init(lb, 6);
285 lb_selftest_print(lb, "init");
286
287 lb_mark_negative(lb, 0, 1);
288 lb_selftest_print(lb, "-A1");
289
290 lb_mark_positive(lb, 2, 3);
291 lb_selftest_print(lb, "+C3");
292
293 lb_mark_negative(lb, 0, 2);
294 lb_mark_negative(lb, 0, 4);
295 lb_mark_negative(lb, 0, 5);
296 lb_selftest_print(lb, "+A0");
297
298 lb_mark_positive(lb, 1, 1);
299 lb_mark_positive(lb, 3, 5);
300 lb_mark_negative(lb, 4, 4);
301 lb_selftest_print(lb, "Finish");
302
303 lb_init(lb, 6);
304 lb_mark_negative(lb, 2, 0);
305 lb_mark_negative(lb, 2, 1);
306 lb_mark_negative(lb, 3, 0);
307 lb_mark_negative(lb, 3, 1);
308 lb_mark_negative(lb, 4, 0);
309 lb_mark_negative(lb, 4, 1);
310 lb_mark_negative(lb, 5, 0);
311 lb_mark_negative(lb, 5, 1);
312 lb_selftest_print(lb, "Clique 2x4");
313
314 lb_init(lb, 6);
315 lb_mark_negative(lb, 3, 0);
316 lb_mark_negative(lb, 3, 1);
317 lb_mark_negative(lb, 3, 2);
318 lb_mark_negative(lb, 4, 0);
319 lb_mark_negative(lb, 4, 1);
320 lb_mark_negative(lb, 4, 2);
321 lb_mark_negative(lb, 5, 0);
322 lb_mark_negative(lb, 5, 1);
323 lb_mark_negative(lb, 5, 2);
324 lb_mark_negative(lb, 2, 0);
325 lb_mark_negative(lb, 2, 1);
326 lb_selftest_print(lb, "Clique 3x3");
327
328 lb_init(lb, 6);
329 for (size_t i = 2; i < lb->size; i++) {
330 for (size_t j = 0; j < i && j < 4; j++) {
331 lb_mark_negative(lb, j, i);
332 }
333 }
334 lb_selftest_print(lb, NULL);
335 lb_selftest_print(lb, "diagonal on 6");
336
337 lb_mark_negative(lb, 0, 0);
338 lb_mark_negative(lb, 4, 4);
339 lb_selftest_print(lb, "full solve on 6");
340
341 // 5x5 tests
342 lb_init(lb, 5);
343 lb_mark_negative(lb, 2, 0);
344 lb_mark_negative(lb, 2, 1);
345 lb_mark_negative(lb, 3, 0);
346 lb_mark_negative(lb, 3, 1);
347 lb_mark_negative(lb, 4, 0);
348 lb_mark_negative(lb, 4, 1);
349 lb_selftest_print(lb, "Clique 2x3");
350
351 lb_init(lb, 5);
352 for (size_t i = 2; i < lb->size; i++) {
353 for (size_t j = 0; j < i && j < 4; j++) {
354 lb_mark_negative(lb, j, i);
355 }
356 }
357 lb_selftest_print(lb, NULL);
358 lb_selftest_print(lb, "diagonal on 5");
359
360 lb_mark_negative(lb, 0, 0);
361 lb_selftest_print(lb, "full solve on 5");
362
363 // 6x6 with equivs
364 lb_init(lb, 6);
365 lb_init_equivalent(lb, 0, 1);
366 lb_init_equivalent(lb, 0, 2);
367 lb_init_equivalent(lb, 0, 3);
368 lb_init_equivalent(lb, 0, 4);
369 lb_selftest_print(lb, "Init");
370
371 lb_mark_negative(lb, 0, 1);
372 lb_selftest_print(lb, "Solved");
373
374 lb_init(lb, 6);
375 lb_init_equivalent(lb, 0, 1);
376 lb_init_equivalent(lb, 0, 2);
377 lb_init_equivalent(lb, 0, 3);
378 lb_init_equivalent(lb, 0, 4);
379 lb_init_equivalent(lb, 0, 5);
380 lb_selftest_print(lb, "Fully solved");
381
382 // 6x6, replicate CT 2 key issue
383 lb_init(lb, 6);
384 lb_init_equivalent(lb, 0, 1);
385 lb_selftest_print(lb, "Init");
386
387 lb_mark_negative(lb, 0, 2);
388 lb_mark_negative(lb, 0, 3);
389 lb_mark_negative(lb, 0, 4);
390 lb_mark_negative(lb, 0, 5);
391 lb_selftest_print(lb, "Step 2");
392
393 // 5x5, mark equivalents positive
394 lb_init(lb, 5);
395 lb_init_equivalent(lb, 0, 1);
396 lb_init_equivalent(lb, 0, 2);
397 lb_selftest_print(lb, "Init");
398
399 lb_mark_positive(lb, 0, 0);
400 lb_mark_positive(lb, 1, 1);
401 lb_selftest_print(lb, "Step 2");
402
403 lb_mark_positive(lb, 0, 2);
404 lb_selftest_print(lb, "Step 3");
405
406 return;
407
408 // 256x256 test
409 lb_init(lb, 256);
410 for (size_t i = 3; i < lb->size; i++) {
411 for (size_t j = 0; j < i && j < lb->size - 3; j++) {
412 lb_mark_negative(lb, j, i);
413 }
414 }
415 lb_deduce(lb);
416 //lb_selftest_print(lb, "diagonal on 256x256");
417}
418
419#endif
lb.c