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