5 #include "heuristical.h"
6 #include "html_dumper.h"
10 #include "pbqp_edge_t.h"
11 #include "pbqp_node.h"
12 #include "pbqp_node_t.h"
15 static pbqp_edge **edge_bucket;
16 static pbqp_node **node_buckets[4];
17 static pbqp_node **reduced_bucket = NULL;
18 static int buckets_filled = 0;
20 static void insert_into_edge_bucket(pbqp_edge *edge)
22 unsigned bucket_len = ARR_LEN(edge_bucket);
24 if (edge->bucket_index < bucket_len && edge_bucket[edge->bucket_index]
26 /* Edge is already inserted. */
29 edge->bucket_index = bucket_len;
30 ARR_APP1(pbqp_edge *, edge_bucket, edge);
33 static void init_buckets(void)
37 edge_bucket = NEW_ARR_F(pbqp_edge *, 0);
38 reduced_bucket = NEW_ARR_F(pbqp_node *, 0);
40 for (i = 0; i < 4; ++i) {
41 node_buckets[i] = NEW_ARR_F(pbqp_node *, 0);
45 static void free_buckets(void)
49 for (i = 0; i < 4; ++i) {
50 DEL_ARR_F(node_buckets[i]);
51 node_buckets[i] = NULL;
54 DEL_ARR_F(edge_bucket);
57 DEL_ARR_F(reduced_bucket);
58 reduced_bucket = NULL;
63 static void fill_node_buckets(pbqp *pbqp)
69 node_len = pbqp->num_nodes;
71 for (node_index = 0; node_index < node_len; ++node_index) {
73 pbqp_node *node = get_node(pbqp, node_index);
77 arity = ARR_LEN(node->edges);
79 /* We have only one bucket for nodes with arity >= 3. */
84 node->bucket_index = ARR_LEN(node_buckets[arity]);
86 ARR_APP1(pbqp_node *, node_buckets[arity], node);
92 static void normalize_towards_source(pbqp *pbqp, pbqp_edge *edge)
106 src_node = edge->src;
107 tgt_node = edge->tgt;
111 src_vec = src_node->costs;
112 tgt_vec = tgt_node->costs;
116 src_len = src_vec->len;
117 tgt_len = tgt_vec->len;
124 /* Normalize towards source node. */
125 for (src_index = 0; src_index < src_len; ++src_index) {
126 num min = pbqp_matrix_get_row_min(mat, src_index, tgt_vec);
129 if (src_vec->entries[src_index].data == INF_COSTS) {
130 pbqp_matrix_set_row_value(mat, src_index, 0);
132 pbqp_matrix_sub_row_value(mat, src_index, tgt_vec, min);
134 src_vec->entries[src_index].data = pbqp_add(
135 src_vec->entries[src_index].data, min);
137 if (min == INF_COSTS) {
139 unsigned edge_len = ARR_LEN(src_node->edges);
141 for (edge_index = 0; edge_index < edge_len; ++edge_index) {
142 pbqp_edge *edge_candidate = src_node->edges[edge_index];
143 if (edge_candidate != edge) {
144 insert_into_edge_bucket(edge_candidate);
152 static void normalize_towards_target(pbqp *pbqp, pbqp_edge *edge)
166 src_node = edge->src;
167 tgt_node = edge->tgt;
171 src_vec = src_node->costs;
172 tgt_vec = tgt_node->costs;
176 src_len = src_vec->len;
177 tgt_len = tgt_vec->len;
184 for (tgt_index = 0; tgt_index < tgt_len; ++tgt_index) {
185 num min = pbqp_matrix_get_col_min(mat, tgt_index, src_vec);
188 if (tgt_vec->entries[tgt_index].data == INF_COSTS) {
189 pbqp_matrix_set_col_value(mat, tgt_index, 0);
191 pbqp_matrix_sub_col_value(mat, tgt_index, src_vec, min);
193 tgt_vec->entries[tgt_index].data = pbqp_add(
194 tgt_vec->entries[tgt_index].data, min);
196 if (min == INF_COSTS) {
198 unsigned edge_len = ARR_LEN(tgt_node->edges);
200 for (edge_index = 0; edge_index < edge_len; ++edge_index) {
201 pbqp_edge *edge_candidate = tgt_node->edges[edge_index];
202 if (edge_candidate != edge) {
203 insert_into_edge_bucket(edge_candidate);
211 static void reorder_node(pbqp_node *node)
215 unsigned old_bucket_len;
216 unsigned old_bucket_index;
217 pbqp_node **old_bucket;
220 if (!buckets_filled) return;
224 arity = ARR_LEN(node->edges);
226 /* Same bucket as before */
227 if (arity > 2) return;
229 /* Assume node lost one incident edge. */
230 old_arity = arity + 1;
231 old_bucket = node_buckets[old_arity];
232 old_bucket_len = ARR_LEN(old_bucket);
233 old_bucket_index = node->bucket_index;
235 if (old_bucket_len <= old_bucket_index || old_bucket[old_bucket_index]
237 unsigned bucket_len = ARR_LEN(node_buckets[arity]);
239 /* Old arity is new arity, so we have nothing to do. */
240 assert(old_bucket_index < bucket_len);
241 assert(node_buckets[arity][old_bucket_index] == node);
245 assert(old_bucket[old_bucket_index] == node);
247 /* Delete node from old bucket... */
248 other = old_bucket[old_bucket_len - 1];
249 other->bucket_index = old_bucket_index;
250 old_bucket[old_bucket_index] = other;
251 ARR_SHRINKLEN(node_buckets[old_arity], old_bucket_len - 1);
253 /* ..and add to new one. */
254 node->bucket_index = ARR_LEN(node_buckets[arity]);
255 ARR_APP1(pbqp_node*, node_buckets[arity], node);
258 static void simplify_edge(pbqp *pbqp, pbqp_edge *edge)
271 src_node = edge->src;
272 tgt_node = edge->tgt;
276 /* If edge are already deleted, we have nothing to do. */
277 if (!is_connected(src_node, edge) || !is_connected(tgt_node, edge))
280 if (pbqp->dump_file) {
282 sprintf(txt, "Simplification of Edge n%d-n%d", src_node->index, tgt_node->index);
283 dump_section(pbqp->dump_file, 3, txt);
286 src_vec = src_node->costs;
287 tgt_vec = tgt_node->costs;
291 src_len = src_vec->len;
292 tgt_len = tgt_vec->len;
299 if (pbqp->dump_file) {
300 fputs("Input:<br>\n", pbqp->dump_file);
301 dump_simplifyedge(pbqp, edge);
304 normalize_towards_source(pbqp, edge);
305 normalize_towards_target(pbqp, edge);
307 if (pbqp->dump_file) {
308 fputs("<br>\nOutput:<br>\n", pbqp->dump_file);
309 dump_simplifyedge(pbqp, edge);
312 if (pbqp_matrix_is_zero(mat, src_vec, tgt_vec)) {
313 if (pbqp->dump_file) {
314 fputs("edge has been eliminated<br>\n", pbqp->dump_file);
318 reorder_node(src_node);
319 reorder_node(tgt_node);
323 void solve_pbqp_heuristical(pbqp *pbqp)
330 if (pbqp->dump_file) {
331 pbqp_dump_input(pbqp);
332 dump_section(pbqp->dump_file, 1, "2. Simplification of Cost Matrices");
335 node_len = pbqp->num_nodes;
339 /* First simplify all edges. */
340 for (node_index = 0; node_index < node_len; ++node_index) {
342 pbqp_node *node = get_node(pbqp, node_index);
349 edge_len = ARR_LEN(edges);
351 for (edge_index = 0; edge_index < edge_len; ++edge_index) {
352 pbqp_edge *edge = edges[edge_index];
354 /* Simplify only once per edge. */
355 if (node != edge->src) continue;
357 simplify_edge(pbqp, edge);
361 /* Put node into bucket representing their arity. */
362 fill_node_buckets(pbqp);
365 if (ARR_LEN(edge_bucket) > 0) {
367 } else if (ARR_LEN(node_buckets[1]) > 0) {
369 } else if (ARR_LEN(node_buckets[2]) > 0) {
371 } else if (ARR_LEN(node_buckets[3]) > 0) {
378 if (pbqp->dump_file) {
379 dump_section(pbqp->dump_file, 1, "4. Determine Solution/Minimum");
380 dump_section(pbqp->dump_file, 2, "4.1. Trivial Solution");
383 /* Solve trivial nodes and calculate solution. */
384 node_len = ARR_LEN(node_buckets[0]);
385 for (node_index = 0; node_index < node_len; ++node_index) {
386 pbqp_node *node = node_buckets[0][node_index];
389 node->solution = vector_get_min_index(node->costs);
390 pbqp->solution = pbqp_add(pbqp->solution,
391 node->costs->entries[node->solution].data);
392 if (pbqp->dump_file) {
393 fprintf(pbqp->dump_file, "node n%d is set to %d<br>\n", node->index, node->solution);
394 dump_node(pbqp, node);
398 if (pbqp->dump_file) {
399 dump_section(pbqp->dump_file, 2, "Minimum");
400 fprintf(pbqp->dump_file, "Minimum is equal to %d.", pbqp->solution);
401 dump_section(pbqp->dump_file, 2, "Back Propagation");
404 /* Solve reduced nodes. */
405 node_len = ARR_LEN(reduced_bucket);
406 for (node_index = node_len; node_index > 0; --node_index) {
407 pbqp_node *node = reduced_bucket[node_index - 1];
410 switch (ARR_LEN(node->edges)) {
412 back_propagate_RI(pbqp, node);
415 back_propagate_RII(pbqp, node);
418 panic("Only nodes with degree one or two should be in this bucket");
426 void apply_edge(pbqp *pbqp)
428 unsigned bucket_len = ARR_LEN(edge_bucket);
429 pbqp_edge *edge = edge_bucket[bucket_len - 1];
431 ARR_SHRINKLEN(edge_bucket, (int)bucket_len - 1);
433 simplify_edge(pbqp, edge);
436 void apply_RI(pbqp *pbqp)
438 pbqp_node **bucket = node_buckets[1];
439 unsigned bucket_len = ARR_LEN(bucket);
440 pbqp_node *node = bucket[bucket_len - 1];
441 pbqp_edge *edge = node->edges[0];
442 pbqp_matrix *mat = edge->costs;
443 int is_src = edge->src == node;
444 pbqp_node *other_node;
447 other_node = edge->tgt;
449 other_node = edge->src;
452 if (pbqp->dump_file) {
454 sprintf(txt, "RI-Reduction of Node n%d", node->index);
455 dump_section(pbqp->dump_file, 2, txt);
456 pbqp_dump_graph(pbqp);
457 fputs("<br>\nBefore reduction:<br>\n", pbqp->dump_file);
458 dump_node(pbqp, node);
459 dump_node(pbqp, other_node);
460 dump_edge(pbqp, edge);
464 pbqp_matrix_add_to_all_cols(mat, node->costs);
465 normalize_towards_target(pbqp, edge);
467 pbqp_matrix_add_to_all_rows(mat, node->costs);
468 normalize_towards_source(pbqp, edge);
470 disconnect_edge(other_node, edge);
472 if (pbqp->dump_file) {
473 fputs("<br>\nAfter reduction:<br>\n", pbqp->dump_file);
474 dump_node(pbqp, other_node);
477 /* Remove node from bucket... */
478 ARR_SHRINKLEN(bucket, (int)bucket_len - 1);
479 reorder_node(other_node);
481 /* ...and add it to back propagation list. */
482 node->bucket_index = ARR_LEN(reduced_bucket);
483 ARR_APP1(pbqp_node *, reduced_bucket, node);
486 void apply_RII(pbqp *pbqp)
488 pbqp_node **bucket = node_buckets[2];
489 unsigned bucket_len = ARR_LEN(bucket);
490 pbqp_node *node = bucket[bucket_len - 1];
491 pbqp_edge *src_edge = node->edges[0];
492 pbqp_edge *tgt_edge = node->edges[1];
493 int src_is_src = src_edge->src == node;
494 int tgt_is_src = tgt_edge->src == node;
495 pbqp_matrix *src_mat;
496 pbqp_matrix *tgt_mat;
513 src_node = src_edge->tgt;
515 src_node = src_edge->src;
519 tgt_node = tgt_edge->tgt;
521 tgt_node = tgt_edge->src;
524 /* Swap nodes if necessary. */
525 if (tgt_node->index < src_node->index) {
537 src_is_src = src_edge->src == node;
538 tgt_is_src = tgt_edge->src == node;
541 if (pbqp->dump_file) {
543 sprintf(txt, "RII-Reduction of Node n%d", node->index);
544 dump_section(pbqp->dump_file, 2, txt);
545 pbqp_dump_graph(pbqp);
546 fputs("<br>\nBefore reduction:<br>\n", pbqp->dump_file);
547 dump_node(pbqp, src_node);
548 dump_edge(pbqp, src_edge);
549 dump_node(pbqp, node);
550 dump_edge(pbqp, tgt_edge);
551 dump_node(pbqp, tgt_node);
554 src_mat = src_edge->costs;
555 tgt_mat = tgt_edge->costs;
557 src_vec = src_node->costs;
558 tgt_vec = tgt_node->costs;
559 node_vec = node->costs;
561 row_len = src_vec->len;
562 col_len = tgt_vec->len;
563 node_len = node_vec->len;
565 mat = pbqp_matrix_alloc(pbqp, row_len, col_len);
567 for (row_index = 0; row_index < row_len; ++row_index) {
568 for (col_index = 0; col_index < col_len; ++col_index) {
569 vec = vector_copy(pbqp, node_vec);
572 vector_add_matrix_col(vec, src_mat, row_index);
574 vector_add_matrix_row(vec, src_mat, row_index);
578 vector_add_matrix_col(vec, tgt_mat, col_index);
580 vector_add_matrix_row(vec, tgt_mat, col_index);
583 mat->entries[row_index * col_len + col_index] = vector_get_min(vec);
585 obstack_free(&pbqp->obstack, vec);
589 pbqp_edge *edge = get_edge(pbqp, src_node->index, tgt_node->index);
591 /* Disconnect node. */
592 disconnect_edge(src_node, src_edge);
593 disconnect_edge(tgt_node, tgt_edge);
595 /* Remove node from bucket... */
596 ARR_SHRINKLEN(bucket, (int)bucket_len - 1);
598 /* ...and add it to back propagation list. */
599 node->bucket_index = ARR_LEN(reduced_bucket);
600 ARR_APP1(pbqp_node *, reduced_bucket, node);
603 edge = alloc_edge(pbqp, src_node->index, tgt_node->index, mat);
605 pbqp_matrix_add(edge->costs, mat);
607 /* Free local matrix. */
608 obstack_free(&pbqp->obstack, mat);
610 reorder_node(src_node);
611 reorder_node(tgt_node);
614 if (pbqp->dump_file) {
615 fputs("<br>\nAfter reduction:<br>\n", pbqp->dump_file);
616 dump_edge(pbqp, edge);
619 /* Edge has changed so we simplify it. */
620 simplify_edge(pbqp, edge);
623 void apply_RN(pbqp *pbqp)
625 pbqp_node **bucket = node_buckets[3];
626 unsigned bucket_len = ARR_LEN(bucket);
627 pbqp_node *node = bucket[bucket_len - 1];
629 vector *node_vec = node->costs;
633 unsigned edge_len = ARR_LEN(node->edges);
635 unsigned node_len = node_vec->len;
636 unsigned min_index = 0;
642 if (pbqp->dump_file) {
644 sprintf(txt, "RN-Reduction of Node n%d", node->index);
645 dump_section(pbqp->dump_file, 2, txt);
646 pbqp_dump_graph(pbqp);
649 for (node_index = 0; node_index < node_len; ++node_index) {
652 for (edge_index = 0; edge_index < edge_len; ++edge_index) {
653 edge = node->edges[edge_index];
655 is_src = edge->src == node;
658 vec = vector_copy(pbqp, edge->tgt->costs);
659 vector_add_matrix_row(vec, mat, node_index);
661 vec = vector_copy(pbqp, edge->src->costs);
662 vector_add_matrix_col(vec, mat, node_index);
665 value = pbqp_add(value, vector_get_min(vec));
667 obstack_free(&pbqp->obstack, vec);
672 min_index = node_index;
676 if (pbqp->dump_file) {
677 fprintf(pbqp->dump_file, "node n%d is set to %d<br><br>\n",
678 node->index, min_index);
679 fprintf(pbqp->dump_file, "Minimal cost of RN reduction: %d<br>\n",
683 node->solution = min_index;
685 /* Now that we found the local minimum set all other costs to infinity. */
686 for (node_index = 0; node_index < node_len; ++node_index) {
687 if (node_index != min_index) {
688 node_vec->entries[node_index].data = INF_COSTS;
692 /* Add all incident edges to edge bucket, since they are now independent. */
693 for (edge_index = 0; edge_index < edge_len; ++edge_index) {
694 insert_into_edge_bucket(node->edges[edge_index]);
698 void back_propagate_RI(pbqp *pbqp, pbqp_node *node)
709 edge = node->edges[0];
711 is_src = edge->src == node;
717 vector_add_matrix_col(vec, mat, other->solution);
721 vector_add_matrix_row(vec, mat, other->solution);
724 node->solution = vector_get_min_index(vec);
725 if (pbqp->dump_file) {
726 fprintf(pbqp->dump_file, "node n%d is set to %d<br>\n", node->index, node->solution);
730 void back_propagate_RII(pbqp *pbqp, pbqp_node *node)
732 pbqp_edge *src_edge = node->edges[0];
733 pbqp_edge *tgt_edge = node->edges[1];
734 int src_is_src = src_edge->src == node;
735 int tgt_is_src = tgt_edge->src == node;
736 pbqp_matrix *src_mat;
737 pbqp_matrix *tgt_mat;
748 src_node = src_edge->tgt;
750 src_node = src_edge->src;
754 tgt_node = tgt_edge->tgt;
756 tgt_node = tgt_edge->src;
759 /* Swap nodes if necessary. */
760 if (tgt_node->index < src_node->index) {
772 src_is_src = src_edge->src == node;
773 tgt_is_src = tgt_edge->src == node;
776 src_mat = src_edge->costs;
777 tgt_mat = tgt_edge->costs;
779 node_vec = node->costs;
781 row_index = src_node->solution;
782 col_index = tgt_node->solution;
784 vec = vector_copy(pbqp, node_vec);
787 vector_add_matrix_col(vec, src_mat, row_index);
789 vector_add_matrix_row(vec, src_mat, row_index);
793 vector_add_matrix_col(vec, tgt_mat, col_index);
795 vector_add_matrix_row(vec, tgt_mat, col_index);
798 node->solution = vector_get_min_index(vec);
799 if (pbqp->dump_file) {
800 fprintf(pbqp->dump_file, "node n%d is set to %d<br>\n", node->index, node->solution);
803 obstack_free(&pbqp->obstack, vec);
806 int node_is_reduced(pbqp_node *node)
808 if (!reduced_bucket) return 0;
811 if (ARR_LEN(node->edges) == 0) return 1;
813 unsigned bucket_length = ARR_LEN(reduced_bucket);
814 unsigned bucket_index = node->bucket_index;
816 return bucket_index < bucket_length && reduced_bucket[bucket_index] == node;