- introduce some defines to enable/disable PBQP dumping etc.
[libfirm] / heuristical.c
1 #include "adt/array.h"
2 #include "assert.h"
3 #include "error.h"
4
5 #include "bucket.h"
6 #include "heuristical.h"
7 #if     KAPS_DUMP
8 #include "html_dumper.h"
9 #endif
10 #include "kaps.h"
11 #include "matrix.h"
12 #include "pbqp_edge.h"
13 #include "pbqp_edge_t.h"
14 #include "pbqp_node.h"
15 #include "pbqp_node_t.h"
16 #include "vector.h"
17
18 static pbqp_edge **edge_bucket;
19 static pbqp_node **node_buckets[4];
20 static pbqp_node **reduced_bucket = NULL;
21 static int         buckets_filled = 0;
22
23 #if KAPS_STATISTIC
24 static int dump = 0;
25 #endif
26
27 /* Forward declarations. */
28 static void apply_Brute_Force(pbqp *pbqp);
29
30 static void insert_into_edge_bucket(pbqp_edge *edge)
31 {
32         if (edge_bucket_contains(edge_bucket, edge)) {
33                 /* Edge is already inserted. */
34                 return;
35         }
36
37         edge_bucket_insert(&edge_bucket, edge);
38 }
39
40 static void init_buckets(void)
41 {
42         int i;
43
44         edge_bucket_init(&edge_bucket);
45         node_bucket_init(&reduced_bucket);
46
47         for (i = 0; i < 4; ++i) {
48                 node_bucket_init(&node_buckets[i]);
49         }
50 }
51
52 static void free_buckets(void)
53 {
54         int i;
55
56         for (i = 0; i < 4; ++i) {
57                 node_bucket_free(&node_buckets[i]);
58         }
59
60         edge_bucket_free(&edge_bucket);
61         node_bucket_free(&reduced_bucket);
62
63         buckets_filled = 0;
64 }
65
66 static void fill_node_buckets(pbqp *pbqp)
67 {
68         unsigned node_index;
69         unsigned node_len;
70
71         assert(pbqp);
72         node_len = pbqp->num_nodes;
73
74         for (node_index = 0; node_index < node_len; ++node_index) {
75                 unsigned   degree;
76                 pbqp_node *node = get_node(pbqp, node_index);
77
78                 if (!node) continue;
79
80                 degree = pbqp_node_get_degree(node);
81
82                 /* We have only one bucket for nodes with arity >= 3. */
83                 if (degree > 3) {
84                         degree = 3;
85                 }
86
87                 node_bucket_insert(&node_buckets[degree], node);
88         }
89
90         buckets_filled = 1;
91 }
92
93 static void normalize_towards_source(pbqp *pbqp, pbqp_edge *edge)
94 {
95         pbqp_matrix    *mat;
96         pbqp_node      *src_node;
97         pbqp_node      *tgt_node;
98         vector         *src_vec;
99         vector         *tgt_vec;
100         int             src_len;
101         int             tgt_len;
102         int             src_index;
103
104         assert(pbqp);
105         assert(edge);
106
107         src_node = edge->src;
108         tgt_node = edge->tgt;
109         assert(src_node);
110         assert(tgt_node);
111
112         src_vec = src_node->costs;
113         tgt_vec = tgt_node->costs;
114         assert(src_vec);
115         assert(tgt_vec);
116
117         src_len = src_vec->len;
118         tgt_len = tgt_vec->len;
119         assert(src_len > 0);
120         assert(tgt_len > 0);
121
122         mat = edge->costs;
123         assert(mat);
124
125         /* Normalize towards source node. */
126         for (src_index = 0; src_index < src_len; ++src_index) {
127                 num min = pbqp_matrix_get_row_min(mat, src_index, tgt_vec);
128
129                 if (min != 0) {
130                         if (src_vec->entries[src_index].data == INF_COSTS) {
131                                 pbqp_matrix_set_row_value(mat, src_index, 0);
132                         } else {
133                                 pbqp_matrix_sub_row_value(mat, src_index, tgt_vec, min);
134                         }
135                         src_vec->entries[src_index].data = pbqp_add(
136                                         src_vec->entries[src_index].data, min);
137
138                         if (min == INF_COSTS) {
139                                 unsigned edge_index;
140                                 unsigned edge_len = pbqp_node_get_degree(src_node);
141
142                                 for (edge_index = 0; edge_index < edge_len; ++edge_index) {
143                                         pbqp_edge *edge_candidate = src_node->edges[edge_index];
144                                         if (edge_candidate != edge) {
145                                                 insert_into_edge_bucket(edge_candidate);
146                                         }
147                                 }
148                         }
149                 }
150         }
151 }
152
153 static void normalize_towards_target(pbqp *pbqp, pbqp_edge *edge)
154 {
155         pbqp_matrix    *mat;
156         pbqp_node      *src_node;
157         pbqp_node      *tgt_node;
158         vector         *src_vec;
159         vector         *tgt_vec;
160         int             src_len;
161         int             tgt_len;
162         int             tgt_index;
163
164         assert(pbqp);
165         assert(edge);
166
167         src_node = edge->src;
168         tgt_node = edge->tgt;
169         assert(src_node);
170         assert(tgt_node);
171
172         src_vec = src_node->costs;
173         tgt_vec = tgt_node->costs;
174         assert(src_vec);
175         assert(tgt_vec);
176
177         src_len = src_vec->len;
178         tgt_len = tgt_vec->len;
179         assert(src_len > 0);
180         assert(tgt_len > 0);
181
182         mat = edge->costs;
183         assert(mat);
184
185         for (tgt_index = 0; tgt_index < tgt_len; ++tgt_index) {
186                 num min = pbqp_matrix_get_col_min(mat, tgt_index, src_vec);
187
188                 if (min != 0) {
189                         if (tgt_vec->entries[tgt_index].data == INF_COSTS) {
190                                 pbqp_matrix_set_col_value(mat, tgt_index, 0);
191                         } else {
192                                 pbqp_matrix_sub_col_value(mat, tgt_index, src_vec, min);
193                         }
194                         tgt_vec->entries[tgt_index].data = pbqp_add(
195                                         tgt_vec->entries[tgt_index].data, min);
196
197                         if (min == INF_COSTS) {
198                                 unsigned edge_index;
199                                 unsigned edge_len = pbqp_node_get_degree(tgt_node);
200
201                                 for (edge_index = 0; edge_index < edge_len; ++edge_index) {
202                                         pbqp_edge *edge_candidate = tgt_node->edges[edge_index];
203                                         if (edge_candidate != edge) {
204                                                 insert_into_edge_bucket(edge_candidate);
205                                         }
206                                 }
207                         }
208                 }
209         }
210 }
211
212 static void reorder_node(pbqp_node *node)
213 {
214         unsigned    degree     = pbqp_node_get_degree(node);
215         /* Assume node lost one incident edge. */
216         unsigned    old_degree = degree + 1;
217
218         if (!buckets_filled) return;
219
220         /* Same bucket as before */
221         if (degree > 2) return;
222
223         if (!node_bucket_contains(node_buckets[old_degree], node)) {
224                 /* Old arity is new arity, so we have nothing to do. */
225                 assert(node_bucket_contains(node_buckets[degree], node));
226                 return;
227         }
228
229         /* Delete node from old bucket... */
230         node_bucket_remove(&node_buckets[old_degree], node);
231
232         /* ..and add to new one. */
233         node_bucket_insert(&node_buckets[degree], node);
234 }
235
236 #if 0
237 static void check_melting_possibility(pbqp *pbqp, pbqp_edge *edge)
238 {
239         pbqp_matrix    *mat;
240         pbqp_node      *src_node;
241         pbqp_node      *tgt_node;
242         vector         *src_vec;
243         vector         *tgt_vec;
244         int             src_len;
245         int             tgt_len;
246         int             src_index;
247         int             tgt_index;
248
249         assert(pbqp);
250         assert(edge);
251
252         src_node = edge->src;
253         tgt_node = edge->tgt;
254         assert(src_node);
255         assert(tgt_node);
256
257         src_vec = src_node->costs;
258         tgt_vec = tgt_node->costs;
259         assert(src_vec);
260         assert(tgt_vec);
261
262         src_len = src_vec->len;
263         tgt_len = tgt_vec->len;
264         assert(src_len > 0);
265         assert(tgt_len > 0);
266
267         mat = edge->costs;
268         assert(mat);
269
270         if (src_len == 1 && tgt_len == 1) {
271                 //panic("Something is wrong");
272         }
273
274         int allRowsOk = 1;
275         for (src_index = 0; src_index < src_len; ++src_index) {
276                 int onlyOneZero = 0;
277                 if (src_vec->entries[src_index].data == INF_COSTS) {
278                         continue;
279                 }
280                 for (tgt_index = 0; tgt_index < tgt_len; ++tgt_index) {
281                         if (tgt_vec->entries[tgt_index].data == INF_COSTS) {
282                                 continue;
283                         }
284                         if (mat->entries[src_index * tgt_len + tgt_index] == 0) {
285                                 if (onlyOneZero) {
286                                         onlyOneZero = 0;
287                                         break;
288                                 } else {
289                                         onlyOneZero = 1;
290                                         continue;
291                                 }
292                         }
293                         if (mat->entries[src_index * tgt_len + tgt_index] == INF_COSTS) {
294                                 continue;
295                         }
296                         onlyOneZero = 0;
297                         break;
298                 }
299                 allRowsOk &= onlyOneZero;
300         }
301
302         int allColsOk = 1;
303         for (tgt_index = 0; tgt_index < tgt_len; ++tgt_index) {
304                 int onlyOneZero = 0;
305                 if (tgt_vec->entries[tgt_index].data == INF_COSTS) {
306                         continue;
307                 }
308                 for (src_index = 0; src_index < src_len; ++src_index) {
309                         if (src_vec->entries[src_index].data == INF_COSTS) {
310                                 continue;
311                         }
312                         if (mat->entries[src_index * tgt_len + tgt_index] == 0) {
313                                 if (onlyOneZero) {
314                                         onlyOneZero = 0;
315                                         break;
316                                 } else {
317                                         onlyOneZero = 1;
318                                         continue;
319                                 }
320                         }
321                         if (mat->entries[src_index * tgt_len + tgt_index] == INF_COSTS) {
322                                 continue;
323                         }
324                         onlyOneZero = 0;
325                         break;
326                 }
327                 allColsOk &= onlyOneZero;
328         }
329
330         if (allRowsOk && allColsOk) {
331                 panic("Hurray");
332         }
333 }
334 #endif
335
336 static void simplify_edge(pbqp *pbqp, pbqp_edge *edge)
337 {
338         pbqp_matrix    *mat;
339         pbqp_node      *src_node;
340         pbqp_node      *tgt_node;
341         vector         *src_vec;
342         vector         *tgt_vec;
343         int             src_len;
344         int             tgt_len;
345
346         assert(pbqp);
347         assert(edge);
348
349         src_node = edge->src;
350         tgt_node = edge->tgt;
351         assert(src_node);
352         assert(tgt_node);
353
354         /* If edge are already deleted, we have nothing to do. */
355         if (!is_connected(src_node, edge) || !is_connected(tgt_node, edge))
356                 return;
357
358 #if     KAPS_DUMP
359         if (pbqp->dump_file) {
360                 char txt[100];
361                 sprintf(txt, "Simplification of Edge n%d-n%d", src_node->index, tgt_node->index);
362                 dump_section(pbqp->dump_file, 3, txt);
363         }
364 #endif
365
366         src_vec = src_node->costs;
367         tgt_vec = tgt_node->costs;
368         assert(src_vec);
369         assert(tgt_vec);
370
371         src_len = src_vec->len;
372         tgt_len = tgt_vec->len;
373         assert(src_len > 0);
374         assert(tgt_len > 0);
375
376         mat = edge->costs;
377         assert(mat);
378
379 #if     KAPS_DUMP
380         if (pbqp->dump_file) {
381                 fputs("Input:<br>\n", pbqp->dump_file);
382                 dump_simplifyedge(pbqp, edge);
383         }
384 #endif
385
386         normalize_towards_source(pbqp, edge);
387         normalize_towards_target(pbqp, edge);
388
389 #if     KAPS_DUMP
390         if (pbqp->dump_file) {
391                 fputs("<br>\nOutput:<br>\n", pbqp->dump_file);
392                 dump_simplifyedge(pbqp, edge);
393         }
394 #endif
395
396         if (pbqp_matrix_is_zero(mat, src_vec, tgt_vec)) {
397 #if     KAPS_DUMP
398                 if (pbqp->dump_file) {
399                         fputs("edge has been eliminated<br>\n", pbqp->dump_file);
400                 }
401 #endif
402
403 #if KAPS_STATISTIC
404                 if (dump == 0) {
405                         pbqp->num_edges++;
406                 }
407 #endif
408
409                 delete_edge(edge);
410                 reorder_node(src_node);
411                 reorder_node(tgt_node);
412         }
413 }
414
415 static void initial_simplify_edges(pbqp *pbqp)
416 {
417         unsigned node_index;
418         unsigned node_len;
419
420         assert(pbqp);
421
422 #if     KAPS_DUMP
423         if (pbqp->dump_file) {
424                 pbqp_dump_input(pbqp);
425                 dump_section(pbqp->dump_file, 1, "2. Simplification of Cost Matrices");
426         }
427 #endif
428
429         node_len = pbqp->num_nodes;
430
431         init_buckets();
432
433         /* First simplify all edges. */
434         for (node_index = 0; node_index < node_len; ++node_index) {
435                 unsigned    edge_index;
436                 pbqp_node  *node = get_node(pbqp, node_index);
437                 pbqp_edge **edges;
438                 unsigned    edge_len;
439
440                 if (!node) continue;
441
442                 edges = node->edges;
443                 edge_len = pbqp_node_get_degree(node);
444
445                 for (edge_index = 0; edge_index < edge_len; ++edge_index) {
446                         pbqp_edge *edge = edges[edge_index];
447
448                         /* Simplify only once per edge. */
449                         if (node != edge->src) continue;
450
451                         simplify_edge(pbqp, edge);
452                 }
453         }
454 }
455
456 static num determine_solution(pbqp *pbqp)
457 {
458         unsigned node_index;
459         unsigned node_len;
460         num      solution   = 0;
461 #if     KAPS_DUMP
462         FILE     *file;
463 #endif
464
465         assert(pbqp);
466
467 #if     KAPS_DUMP
468         file = pbqp->dump_file;
469
470         if (file) {
471                 dump_section(file, 1, "4. Determine Solution/Minimum");
472                 dump_section(file, 2, "4.1. Trivial Solution");
473         }
474 #endif
475
476         /* Solve trivial nodes and calculate solution. */
477         node_len = node_bucket_get_length(node_buckets[0]);
478
479 #if KAPS_STATISTIC
480         if (dump == 0) {
481                 pbqp->num_r0 = node_len;
482         }
483 #endif
484
485         for (node_index = 0; node_index < node_len; ++node_index) {
486                 pbqp_node *node = node_buckets[0][node_index];
487                 assert(node);
488
489                 node->solution = vector_get_min_index(node->costs);
490                 solution       = pbqp_add(solution,
491                                 node->costs->entries[node->solution].data);
492
493 #if     KAPS_DUMP
494                 if (file) {
495                         fprintf(file, "node n%d is set to %d<br>\n", node->index, node->solution);
496                         dump_node(file, node);
497                 }
498 #endif
499         }
500
501 #if     KAPS_DUMP
502         if (file) {
503                 dump_section(file, 2, "Minimum");
504                 fprintf(file, "Minimum is equal to %lld.", solution);
505         }
506 #endif
507
508         return solution;
509 }
510
511 static void back_propagate(pbqp *pbqp)
512 {
513         unsigned node_index;
514         unsigned node_len   = node_bucket_get_length(reduced_bucket);
515
516         assert(pbqp);
517
518 #if     KAPS_DUMP
519         if (pbqp->dump_file) {
520                 dump_section(pbqp->dump_file, 2, "Back Propagation");
521         }
522 #endif
523
524         for (node_index = node_len; node_index > 0; --node_index) {
525                 pbqp_node *node = reduced_bucket[node_index - 1];
526
527                 switch (pbqp_node_get_degree(node)) {
528                         case 1:
529                                 back_propagate_RI(pbqp, node);
530                                 break;
531                         case 2:
532                                 back_propagate_RII(pbqp, node);
533                                 break;
534                         default:
535                                 panic("Only nodes with degree one or two should be in this bucket");
536                                 break;
537                 }
538         }
539 }
540
541 static void apply_heuristic_reductions(pbqp *pbqp)
542 {
543         for (;;) {
544                 if (edge_bucket_get_length(edge_bucket) > 0) {
545                         apply_edge(pbqp);
546                 } else if (node_bucket_get_length(node_buckets[1]) > 0) {
547                         apply_RI(pbqp);
548                 } else if (node_bucket_get_length(node_buckets[2]) > 0) {
549                         apply_RII(pbqp);
550                 } else if (node_bucket_get_length(node_buckets[3]) > 0) {
551                         apply_RN(pbqp);
552                 } else {
553                         return;
554                 }
555         }
556 }
557
558 void solve_pbqp_heuristical(pbqp *pbqp)
559 {
560         /* Reduce nodes degree ... */
561         initial_simplify_edges(pbqp);
562
563         /* ... and put node into bucket representing their degree. */
564         fill_node_buckets(pbqp);
565
566 #if KAPS_STATISTIC
567         FILE *fh = fopen("solutions.pb", "a");
568         fprintf(fh, "Solution");
569         fclose(fh);
570 #endif
571
572         apply_heuristic_reductions(pbqp);
573
574         pbqp->solution = determine_solution(pbqp);
575
576 #if KAPS_STATISTIC
577         fh = fopen("solutions.pb", "a");
578         fprintf(fh, ": %lld RE:%u R0:%u R1:%u R2:%u RN/BF:%u\n", pbqp->solution,
579                                 pbqp->num_edges, pbqp->num_r0, pbqp->num_r1, pbqp->num_r2,
580                                 pbqp->num_rn);
581         fclose(fh);
582 #endif
583
584         /* Solve reduced nodes. */
585         back_propagate(pbqp);
586
587         free_buckets();
588 }
589
590 void apply_edge(pbqp *pbqp)
591 {
592         pbqp_edge *edge = edge_bucket_pop(&edge_bucket);
593
594         simplify_edge(pbqp, edge);
595 }
596
597 void apply_RI(pbqp *pbqp)
598 {
599         pbqp_node   *node       = node_bucket_pop(&node_buckets[1]);
600         pbqp_edge   *edge       = node->edges[0];
601         pbqp_matrix *mat        = edge->costs;
602         int          is_src     = edge->src == node;
603         pbqp_node   *other_node;
604
605         assert(pbqp_node_get_degree(node) == 1);
606
607         if (is_src) {
608                 other_node = edge->tgt;
609         } else {
610                 other_node = edge->src;
611         }
612
613 #if     KAPS_DUMP
614         if (pbqp->dump_file) {
615                 char     txt[100];
616                 sprintf(txt, "RI-Reduction of Node n%d", node->index);
617                 dump_section(pbqp->dump_file, 2, txt);
618                 pbqp_dump_graph(pbqp);
619                 fputs("<br>\nBefore reduction:<br>\n", pbqp->dump_file);
620                 dump_node(pbqp->dump_file, node);
621                 dump_node(pbqp->dump_file, other_node);
622                 dump_edge(pbqp->dump_file, edge);
623         }
624 #endif
625
626         if (is_src) {
627                 pbqp_matrix_add_to_all_cols(mat, node->costs);
628                 normalize_towards_target(pbqp, edge);
629         } else {
630                 pbqp_matrix_add_to_all_rows(mat, node->costs);
631                 normalize_towards_source(pbqp, edge);
632         }
633         disconnect_edge(other_node, edge);
634
635 #if     KAPS_DUMP
636         if (pbqp->dump_file) {
637                 fputs("<br>\nAfter reduction:<br>\n", pbqp->dump_file);
638                 dump_node(pbqp->dump_file, other_node);
639         }
640 #endif
641
642         reorder_node(other_node);
643
644 #if KAPS_STATISTIC
645         if (dump == 0) {
646                 pbqp->num_r1++;
647         }
648 #endif
649
650         /* Add node to back propagation list. */
651         node_bucket_insert(&reduced_bucket, node);
652 }
653
654 void apply_RII(pbqp *pbqp)
655 {
656         pbqp_node   *node       = node_bucket_pop(&node_buckets[2]);
657         pbqp_edge   *src_edge   = node->edges[0];
658         pbqp_edge   *tgt_edge   = node->edges[1];
659         int          src_is_src = src_edge->src == node;
660         int          tgt_is_src = tgt_edge->src == node;
661         pbqp_matrix *src_mat;
662         pbqp_matrix *tgt_mat;
663         pbqp_node   *src_node;
664         pbqp_node   *tgt_node;
665         pbqp_matrix *mat;
666         vector      *vec;
667         vector      *node_vec;
668         vector      *src_vec;
669         vector      *tgt_vec;
670         unsigned     col_index;
671         unsigned     col_len;
672         unsigned     row_index;
673         unsigned     row_len;
674         unsigned     node_len;
675
676         assert(pbqp);
677         assert(pbqp_node_get_degree(node) == 2);
678
679         if (src_is_src) {
680                 src_node = src_edge->tgt;
681         } else {
682                 src_node = src_edge->src;
683         }
684
685         if (tgt_is_src) {
686                 tgt_node = tgt_edge->tgt;
687         } else {
688                 tgt_node = tgt_edge->src;
689         }
690
691         /* Swap nodes if necessary. */
692         if (tgt_node->index < src_node->index) {
693                 pbqp_node *tmp_node;
694                 pbqp_edge *tmp_edge;
695
696                 tmp_node = src_node;
697                 src_node = tgt_node;
698                 tgt_node = tmp_node;
699
700                 tmp_edge = src_edge;
701                 src_edge = tgt_edge;
702                 tgt_edge = tmp_edge;
703
704                 src_is_src = src_edge->src == node;
705                 tgt_is_src = tgt_edge->src == node;
706         }
707
708 #if     KAPS_DUMP
709         if (pbqp->dump_file) {
710                 char     txt[100];
711                 sprintf(txt, "RII-Reduction of Node n%d", node->index);
712                 dump_section(pbqp->dump_file, 2, txt);
713                 pbqp_dump_graph(pbqp);
714                 fputs("<br>\nBefore reduction:<br>\n", pbqp->dump_file);
715                 dump_node(pbqp->dump_file, src_node);
716                 dump_edge(pbqp->dump_file, src_edge);
717                 dump_node(pbqp->dump_file, node);
718                 dump_edge(pbqp->dump_file, tgt_edge);
719                 dump_node(pbqp->dump_file, tgt_node);
720         }
721 #endif
722
723         src_mat = src_edge->costs;
724         tgt_mat = tgt_edge->costs;
725
726         src_vec  = src_node->costs;
727         tgt_vec  = tgt_node->costs;
728         node_vec = node->costs;
729
730         row_len  = src_vec->len;
731         col_len  = tgt_vec->len;
732         node_len = node_vec->len;
733
734         mat = pbqp_matrix_alloc(pbqp, row_len, col_len);
735
736         for (row_index = 0; row_index < row_len; ++row_index) {
737                 for (col_index = 0; col_index < col_len; ++col_index) {
738                         vec = vector_copy(pbqp, node_vec);
739
740                         if (src_is_src) {
741                                 vector_add_matrix_col(vec, src_mat, row_index);
742                         } else {
743                                 vector_add_matrix_row(vec, src_mat, row_index);
744                         }
745
746                         if (tgt_is_src) {
747                                 vector_add_matrix_col(vec, tgt_mat, col_index);
748                         } else {
749                                 vector_add_matrix_row(vec, tgt_mat, col_index);
750                         }
751
752                         mat->entries[row_index * col_len + col_index] = vector_get_min(vec);
753
754                         obstack_free(&pbqp->obstack, vec);
755                 }
756         }
757
758         pbqp_edge *edge = get_edge(pbqp, src_node->index, tgt_node->index);
759
760         /* Disconnect node. */
761         disconnect_edge(src_node, src_edge);
762         disconnect_edge(tgt_node, tgt_edge);
763
764 #if KAPS_STATISTIC
765         if (dump == 0) {
766                 pbqp->num_r2++;
767         }
768 #endif
769
770         /* Add node to back propagation list. */
771         node_bucket_insert(&reduced_bucket, node);
772
773         if (edge == NULL) {
774                 edge = alloc_edge(pbqp, src_node->index, tgt_node->index, mat);
775         } else {
776                 pbqp_matrix_add(edge->costs, mat);
777
778                 /* Free local matrix. */
779                 obstack_free(&pbqp->obstack, mat);
780
781                 reorder_node(src_node);
782                 reorder_node(tgt_node);
783         }
784
785 #if     KAPS_DUMP
786         if (pbqp->dump_file) {
787                 fputs("<br>\nAfter reduction:<br>\n", pbqp->dump_file);
788                 dump_edge(pbqp->dump_file, edge);
789         }
790 #endif
791
792         /* Edge has changed so we simplify it. */
793         simplify_edge(pbqp, edge);
794 }
795
796 static void select_alternative(pbqp_node *node, unsigned selected_index)
797 {
798         unsigned  edge_index;
799         unsigned  node_index;
800         unsigned  node_len;
801         vector   *node_vec;
802         unsigned  max_degree = pbqp_node_get_degree(node);
803
804         assert(node);
805         node->solution = selected_index;
806         node_vec = node->costs;
807         node_len = node_vec->len;
808         assert(selected_index < node_len);
809
810         /* Set all other costs to infinity. */
811         for (node_index = 0; node_index < node_len; ++node_index) {
812                 if (node_index != selected_index) {
813                         node_vec->entries[node_index].data = INF_COSTS;
814                 }
815         }
816
817         /* Add all incident edges to edge bucket, since they are now independent. */
818         for (edge_index = 0; edge_index < max_degree; ++edge_index) {
819                 insert_into_edge_bucket(node->edges[edge_index]);
820         }
821 }
822
823 static pbqp_node *get_node_with_max_degree(void)
824 {
825         pbqp_node  **bucket       = node_buckets[3];
826         unsigned     bucket_len   = node_bucket_get_length(bucket);
827         unsigned     bucket_index;
828         unsigned     max_degree   = 0;
829         pbqp_node   *result       = NULL;
830
831         for (bucket_index = 0; bucket_index < bucket_len; ++bucket_index) {
832                 pbqp_node *candidate = bucket[bucket_index];
833                 unsigned   degree    = pbqp_node_get_degree(candidate);
834
835                 if (degree > max_degree) {
836                         result = candidate;
837                         max_degree = degree;
838                 }
839         }
840
841         return result;
842 }
843
844 static unsigned get_local_minimal_alternative(pbqp *pbqp, pbqp_node *node)
845 {
846         pbqp_edge   *edge;
847         vector      *node_vec;
848         vector      *vec;
849         pbqp_matrix *mat;
850         unsigned     edge_index;
851         unsigned     max_degree   = 0;
852         unsigned     node_index;
853         unsigned     node_len;
854         unsigned     min_index    = 0;
855         num          min          = INF_COSTS;
856         int          is_src;
857
858         assert(pbqp);
859         assert(node);
860         node_vec = node->costs;
861         node_len = node_vec->len;
862
863         for (node_index = 0; node_index < node_len; ++node_index) {
864                 num value = node_vec->entries[node_index].data;
865
866                 for (edge_index = 0; edge_index < max_degree; ++edge_index) {
867                         edge   = node->edges[edge_index];
868                         mat    = edge->costs;
869                         is_src = edge->src == node;
870
871                         if (is_src) {
872                                 vec = vector_copy(pbqp, edge->tgt->costs);
873                                 vector_add_matrix_row(vec, mat, node_index);
874                         } else {
875                                 vec = vector_copy(pbqp, edge->src->costs);
876                                 vector_add_matrix_col(vec, mat, node_index);
877                         }
878
879                         value = pbqp_add(value, vector_get_min(vec));
880
881                         obstack_free(&pbqp->obstack, vec);
882                 }
883
884                 if (value < min) {
885                         min = value;
886                         min_index = node_index;
887                 }
888         }
889
890         return min_index;
891 }
892
893 void apply_RN(pbqp *pbqp)
894 {
895         pbqp_node   *node         = NULL;
896         unsigned     min_index    = 0;
897
898         assert(pbqp);
899
900         /* We want to reduce a node with maximum degree. */
901         node = get_node_with_max_degree();
902         assert(node);
903         assert(pbqp_node_get_degree(node) > 2);
904
905 #if     KAPS_DUMP
906         if (pbqp->dump_file) {
907                 char     txt[100];
908                 sprintf(txt, "RN-Reduction of Node n%d", node->index);
909                 dump_section(pbqp->dump_file, 2, txt);
910                 pbqp_dump_graph(pbqp);
911         }
912 #endif
913
914         min_index = get_local_minimal_alternative(pbqp, node);
915
916 #if     KAPS_DUMP
917         if (pbqp->dump_file) {
918                 fprintf(pbqp->dump_file, "node n%d is set to %d<br><br>\n",
919                                         node->index, min_index);
920         }
921 #endif
922
923 #if KAPS_STATISTIC
924         if (dump == 0) {
925                 FILE *fh = fopen("solutions.pb", "a");
926                 fprintf(fh, "[%u]", min_index);
927                 fclose(fh);
928                 pbqp->num_rn++;
929         }
930 #endif
931
932         /* Now that we found the local minimum set all other costs to infinity. */
933         select_alternative(node, min_index);
934 }
935
936 static void apply_brute_force_reductions(pbqp *pbqp)
937 {
938         for (;;) {
939                 if (edge_bucket_get_length(edge_bucket) > 0) {
940                         apply_edge(pbqp);
941                 } else if (node_bucket_get_length(node_buckets[1]) > 0) {
942                         apply_RI(pbqp);
943                 } else if (node_bucket_get_length(node_buckets[2]) > 0) {
944                         apply_RII(pbqp);
945                 } else if (node_bucket_get_length(node_buckets[3]) > 0) {
946                         apply_Brute_Force(pbqp);
947                 } else {
948                         return;
949                 }
950         }
951 }
952
953 static unsigned get_minimal_alternative(pbqp *pbqp, pbqp_node *node)
954 {
955         vector      *node_vec;
956         unsigned     node_index;
957         unsigned     node_len;
958         unsigned     min_index    = 0;
959         num          min          = INF_COSTS;
960         unsigned     bucket_index;
961
962         assert(pbqp);
963         assert(node);
964         node_vec     = node->costs;
965         node_len     = node_vec->len;
966         bucket_index = node->bucket_index;
967
968         for (node_index = 0; node_index < node_len; ++node_index) {
969                 pbqp_node_bucket bucket_deg3;
970                 num              value;
971                 unsigned         bucket_0_length;
972                 unsigned         bucket_red_length;
973
974                 char *tmp = obstack_finish(&pbqp->obstack);
975
976                 node_bucket_init(&bucket_deg3);
977
978                 /* Some node buckets and the edge bucket should be empty. */
979                 assert(node_bucket_get_length(node_buckets[1]) == 0);
980                 assert(node_bucket_get_length(node_buckets[2]) == 0);
981                 assert(edge_bucket_get_length(edge_bucket)     == 0);
982
983                 /* char *tmp = obstack_finish(&pbqp->obstack); */
984
985                 /* Save current PBQP state. */
986                 node_bucket_copy(&bucket_deg3, node_buckets[3]);
987                 node_bucket_shrink(&node_buckets[3], 0);
988                 node_bucket_deep_copy(pbqp, &node_buckets[3], bucket_deg3);
989                 node_bucket_update(pbqp, node_buckets[3]);
990                 bucket_0_length   = node_bucket_get_length(node_buckets[0]);
991                 bucket_red_length = node_bucket_get_length(reduced_bucket);
992
993                 /* Select alternative and solve PBQP recursively. */
994                 select_alternative(node_buckets[3][bucket_index], node_index);
995                 apply_brute_force_reductions(pbqp);
996
997                 value = determine_solution(pbqp);
998
999                 if (value < min) {
1000                         min = value;
1001                         min_index = node_index;
1002                 }
1003
1004                 /* Some node buckets and the edge bucket should still be empty. */
1005                 assert(node_bucket_get_length(node_buckets[1]) == 0);
1006                 assert(node_bucket_get_length(node_buckets[2]) == 0);
1007                 assert(edge_bucket_get_length(edge_bucket)     == 0);
1008
1009                 /* Clear modified buckets... */
1010                 node_bucket_shrink(&node_buckets[3], 0);
1011
1012                 /* ... and restore old PBQP state. */
1013                 node_bucket_shrink(&node_buckets[0], bucket_0_length);
1014                 node_bucket_shrink(&reduced_bucket, bucket_red_length);
1015                 node_bucket_copy(&node_buckets[3], bucket_deg3);
1016                 node_bucket_update(pbqp, node_buckets[3]);
1017
1018                 /* Free copies. */
1019                 /* obstack_free(&pbqp->obstack, tmp); */
1020                 node_bucket_free(&bucket_deg3);
1021
1022                 obstack_free(&pbqp->obstack, tmp);
1023         }
1024
1025         return min_index;
1026 }
1027
1028 void apply_Brute_Force(pbqp *pbqp)
1029 {
1030         pbqp_node   *node         = NULL;
1031         unsigned     min_index    = 0;
1032
1033         assert(pbqp);
1034
1035         /* We want to reduce a node with maximum degree. */
1036         node = get_node_with_max_degree();
1037         assert(node);
1038         assert(pbqp_node_get_degree(node) > 2);
1039
1040 #if     KAPS_DUMP
1041         if (pbqp->dump_file) {
1042                 char     txt[100];
1043                 sprintf(txt, "BF-Reduction of Node n%d", node->index);
1044                 dump_section(pbqp->dump_file, 2, txt);
1045                 pbqp_dump_graph(pbqp);
1046         }
1047 #endif
1048
1049 #if KAPS_STATISTIC
1050         dump++;
1051 #endif
1052
1053         min_index = get_minimal_alternative(pbqp, node);
1054         node = pbqp->nodes[node->index];
1055
1056 #if     KAPS_DUMP
1057         if (pbqp->dump_file) {
1058                 fprintf(pbqp->dump_file, "node n%d is set to %d<br><br>\n",
1059                                         node->index, min_index);
1060         }
1061 #endif
1062
1063 #if KAPS_STATISTIC
1064         dump--;
1065         if (dump == 0) {
1066                 FILE *fh = fopen("solutions.pb", "a");
1067                 fprintf(fh, "[%u]", min_index);
1068                 fclose(fh);
1069                 pbqp->num_bf++;
1070         }
1071 #endif
1072
1073         /* Now that we found the minimum set all other costs to infinity. */
1074         select_alternative(node, min_index);
1075 }
1076
1077 void solve_pbqp_brute_force(pbqp *pbqp)
1078 {
1079         /* Reduce nodes degree ... */
1080         initial_simplify_edges(pbqp);
1081
1082         /* ... and put node into bucket representing their degree. */
1083         fill_node_buckets(pbqp);
1084
1085 #if KAPS_STATISTIC
1086         FILE *fh = fopen("solutions.pb", "a");
1087         fprintf(fh, "Solution");
1088         fclose(fh);
1089 #endif
1090
1091         apply_brute_force_reductions(pbqp);
1092
1093         pbqp->solution = determine_solution(pbqp);
1094
1095 #if KAPS_STATISTIC
1096         fh = fopen("solutions.pb", "a");
1097         fprintf(fh, ": %lld RE:%u R0:%u R1:%u R2:%u RN/BF:%u\n", pbqp->solution,
1098                         pbqp->num_edges, pbqp->num_r0, pbqp->num_r1, pbqp->num_r2,
1099                         pbqp->num_bf);
1100         fclose(fh);
1101 #endif
1102
1103         /* Solve reduced nodes. */
1104         back_propagate(pbqp);
1105
1106         free_buckets();
1107 }
1108
1109 void back_propagate_RI(pbqp *pbqp, pbqp_node *node)
1110 {
1111         pbqp_edge   *edge;
1112         pbqp_node   *other;
1113         pbqp_matrix *mat;
1114         vector      *vec;
1115         int          is_src;
1116
1117         assert(pbqp);
1118         assert(node);
1119
1120         edge = node->edges[0];
1121         mat = edge->costs;
1122         is_src = edge->src == node;
1123         vec = node->costs;
1124
1125         if (is_src) {
1126                 other = edge->tgt;
1127                 assert(other);
1128
1129                 /* Update pointer for brute force solver. */
1130                 other = pbqp->nodes[other->index];
1131
1132                 node->solution = pbqp_matrix_get_col_min_index(mat, other->solution, vec);
1133         } else {
1134                 other = edge->src;
1135                 assert(other);
1136
1137                 /* Update pointer for brute force solver. */
1138                 other = pbqp->nodes[other->index];
1139
1140                 node->solution = pbqp_matrix_get_row_min_index(mat, other->solution, vec);
1141         }
1142
1143 #if     KAPS_DUMP
1144         if (pbqp->dump_file) {
1145                 fprintf(pbqp->dump_file, "node n%d is set to %d<br>\n", node->index, node->solution);
1146         }
1147 #endif
1148 }
1149
1150 void back_propagate_RII(pbqp *pbqp, pbqp_node *node)
1151 {
1152         pbqp_edge   *src_edge   = node->edges[0];
1153         pbqp_edge   *tgt_edge   = node->edges[1];
1154         int          src_is_src = src_edge->src == node;
1155         int          tgt_is_src = tgt_edge->src == node;
1156         pbqp_matrix *src_mat;
1157         pbqp_matrix *tgt_mat;
1158         pbqp_node   *src_node;
1159         pbqp_node   *tgt_node;
1160         vector      *vec;
1161         vector      *node_vec;
1162         unsigned     col_index;
1163         unsigned     row_index;
1164
1165         assert(pbqp);
1166
1167         if (src_is_src) {
1168                 src_node = src_edge->tgt;
1169         } else {
1170                 src_node = src_edge->src;
1171         }
1172
1173         if (tgt_is_src) {
1174                 tgt_node = tgt_edge->tgt;
1175         } else {
1176                 tgt_node = tgt_edge->src;
1177         }
1178
1179         /* Swap nodes if necessary. */
1180         if (tgt_node->index < src_node->index) {
1181                 pbqp_node *tmp_node;
1182                 pbqp_edge *tmp_edge;
1183
1184                 tmp_node = src_node;
1185                 src_node = tgt_node;
1186                 tgt_node = tmp_node;
1187
1188                 tmp_edge = src_edge;
1189                 src_edge = tgt_edge;
1190                 tgt_edge = tmp_edge;
1191
1192                 src_is_src = src_edge->src == node;
1193                 tgt_is_src = tgt_edge->src == node;
1194         }
1195
1196         /* Update pointer for brute force solver. */
1197         src_node = pbqp->nodes[src_node->index];
1198         tgt_node = pbqp->nodes[tgt_node->index];
1199
1200         src_mat = src_edge->costs;
1201         tgt_mat = tgt_edge->costs;
1202
1203         node_vec = node->costs;
1204
1205         row_index = src_node->solution;
1206         col_index = tgt_node->solution;
1207
1208         vec = vector_copy(pbqp, node_vec);
1209
1210         if (src_is_src) {
1211                 vector_add_matrix_col(vec, src_mat, row_index);
1212         } else {
1213                 vector_add_matrix_row(vec, src_mat, row_index);
1214         }
1215
1216         if (tgt_is_src) {
1217                 vector_add_matrix_col(vec, tgt_mat, col_index);
1218         } else {
1219                 vector_add_matrix_row(vec, tgt_mat, col_index);
1220         }
1221
1222         node->solution = vector_get_min_index(vec);
1223
1224 #if     KAPS_DUMP
1225         if (pbqp->dump_file) {
1226                 fprintf(pbqp->dump_file, "node n%d is set to %d<br>\n", node->index, node->solution);
1227         }
1228 #endif
1229
1230         obstack_free(&pbqp->obstack, vec);
1231 }
1232
1233 int node_is_reduced(pbqp_node *node)
1234 {
1235         if (!reduced_bucket) return 0;
1236
1237         if (pbqp_node_get_degree(node) == 0) return 1;
1238
1239         return node_bucket_contains(reduced_bucket, node);
1240 }