1 #include <terraces/constraints.hpp>
2
3 #include "io_utils.hpp"
4 #include "trees_impl.hpp"
5 #include "union_find.hpp"
6
7 namespace terraces {
8
operator <<(std::ostream & s,const constraint & c)9 std::ostream& operator<<(std::ostream& s, const constraint& c) {
10 s << "lca(" << c.left << "," << c.shared << ") < lca(" << c.shared << "," << c.right << ")";
11 return s;
12 }
13
operator <<(std::ostream & stream,utils::named_output<constraints,name_map> output)14 std::ostream& operator<<(std::ostream& stream, utils::named_output<constraints, name_map> output) {
15 auto c = output.entry;
16 auto& n = *output.names;
17 stream << "lca(" << n[c.left] << "," << n[c.shared] << ") < lca(" << n[c.shared] << ","
18 << n[c.right] << ")";
19 return stream;
20 }
21
compute_constraints(const std::vector<tree> & trees)22 constraints compute_constraints(const std::vector<tree>& trees) {
23 constraints result;
24 auto num_nodes =
25 (*std::max_element(trees.begin(), trees.end(), [](const tree& a, const tree& b) {
26 return a.size() < b.size();
27 })).size();
28 std::vector<std::pair<index, index>> outermost_nodes(num_nodes, {none, none});
29
30 for (auto& t : trees) {
31 // collect outermost nodes for each subtree (these have lca i)
32 foreach_postorder(t, [&](index i) {
33 auto node = t[i];
34 if (is_leaf(node)) {
35 outermost_nodes[i].first = i;
36 outermost_nodes[i].second = i;
37 } else {
38 outermost_nodes[i].first = outermost_nodes[node.lchild()].first;
39 outermost_nodes[i].second = outermost_nodes[node.rchild()].second;
40 }
41 });
42
43 // extract constraints for each edge
44 foreach_preorder(t, [&](index i) {
45 auto node = t[i];
46 if (!is_leaf(node)) {
47 auto lchild = node.lchild();
48 auto rchild = node.rchild();
49 auto lnode = t[lchild];
50 auto rnode = t[rchild];
51 // taxon of leftmost descendant of i
52 auto i1 = t[outermost_nodes[i].first].taxon();
53 // taxon of rightmost descendant of lchild of i
54 auto i2 = t[outermost_nodes[lchild].second].taxon();
55 // taxon of leftmost descendant of rchild of i
56 auto i3 = t[outermost_nodes[rchild].first].taxon();
57 // taxon of rightmost descendant of i
58 auto i4 = t[outermost_nodes[i].second].taxon();
59
60 // if the left edge is an inner edge
61 if (!is_leaf(lnode)) {
62 result.emplace_back(i2, i1, i4);
63 }
64 // if the right edge is an inner edge
65 if (!is_leaf(rnode)) {
66 result.emplace_back(i3, i4, i1);
67 }
68 }
69 });
70 }
71
72 return result;
73 }
74
deduplicate_constraints(constraints & in_c)75 index deduplicate_constraints(constraints& in_c) {
76 for (auto& c : in_c) {
77 c = {std::min(c.left, c.shared), std::max(c.left, c.shared), c.right};
78 }
79 std::sort(in_c.begin(), in_c.end(), [](constraint a, constraint b) {
80 return std::tie(a.left, a.shared, a.right) < std::tie(b.left, b.shared, b.right);
81 });
82 auto it = std::unique(in_c.begin(), in_c.end());
83 index count = (index)std::distance(it, in_c.end());
84 in_c.erase(it, in_c.end());
85 return count;
86 }
87
88 } // namespace terraces
89