max/edmonds-matching-algorithm
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity

refactor: store attributes in own vectors to improve data locality

Browse source
This commit is contained in:
Maximilian Keßler 2023-11-06 13:24:25 +01:00
parent 391e0761a9
commit 65aa75bc46
Signed by: max
GPG key ID: BCC5A619923C0BA5
4 changed files with 203 additions and 174 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

178
src/edmonds.cpp
View file

@ -2,58 +2,59 @@
#include <iostream>
#include <stack>
#include <tuple>
#include "graph_attributes.h"
using namespace ED;
namespace Edmonds {
void check_integrity(Graph const & graph)
void check_integrity(GraphAttributes const & attrs)
{
for(NodeId id = 0; id < graph.num_nodes(); ++id)
for(NodeId id = 0; id < attrs.num_nodes(); ++id)
{
// Check that μ encodes a valid matching
NodeId matched = graph.matched_neighbor(id);
NodeId matched = attrs.matched_neighbor(id);
if(matched != id)
{
assert(graph.matched_neighbor(matched) == id);
assert(attrs.matched_neighbor(matched) == id);
}
if (graph.is_out_of_forest(id))
if (attrs.is_out_of_forest(id))
{
assert(graph.phi(id) == id);
assert(graph.rho(id) == id);
assert(attrs.phi(id) == id);
assert(attrs.rho(id) == id);
}
else
{
// check for a path to the root, i.e. ρ(node)
NodeId cur_node = id;
while(cur_node != graph.rho(cur_node))
while(cur_node != attrs.rho(cur_node))
{
// If the condition was true, then cur_node is outer, part of a blossom
// and we want to follow its path
// therefore, we check that both φ and μ are not the identity on this node
// and point to vertices that have the same rho
assert(graph.matched_neighbor(cur_node) != cur_node);
assert(graph.phi(cur_node) != cur_node);
assert(graph.rho(graph.matched_neighbor(cur_node)) == graph.rho(cur_node));
assert(graph.rho(graph.phi(cur_node)) == graph.rho(cur_node));
assert(attrs.matched_neighbor(cur_node) != cur_node);
assert(attrs.phi(cur_node) != cur_node);
assert(attrs.rho(attrs.matched_neighbor(cur_node)) == attrs.rho(cur_node));
assert(attrs.rho(attrs.phi(cur_node)) == attrs.rho(cur_node));
// now, walk along the matched edge
cur_node = graph.matched_neighbor(cur_node);
cur_node = attrs.matched_neighbor(cur_node);
// now we want to walk along φ, this will again
// - not be the identity
// - result in a node that has the same rho
assert(graph.phi(cur_node) != cur_node);
assert(graph.rho(graph.phi(cur_node)) == graph.rho(cur_node));
assert(attrs.phi(cur_node) != cur_node);
assert(attrs.rho(attrs.phi(cur_node)) == attrs.rho(cur_node));
cur_node = graph.matched_neighbor(graph.phi(cur_node));
cur_node = attrs.matched_neighbor(attrs.phi(cur_node));
}
}
if (not graph.is_outer(id))
if (not attrs.is_outer(id))
{
assert(graph.rho(id) == id);
assert(attrs.rho(id) == id);
}
}
}
@ -66,66 +67,64 @@ void check_integrity(Graph const & graph)
* blossom forest on the graph when this method is called.
* **/
__attribute__((noinline))
std::vector<NodeId> path_to_forest_root(Graph const & graph, NodeId id)
std::vector<NodeId> path_to_forest_root(GraphAttributes const & attrs, NodeId id)
{
std::vector<NodeId> retval;
retval.push_back(id);
while (graph.matched_neighbor(id) != id)
while (attrs.matched_neighbor(id) != id)
{
id = graph.matched_neighbor(id);
id = attrs.matched_neighbor(id);
retval.push_back(id);
// Note that it is guaranteed that this does not produce a loop:
// We are traversing the path to a root of the forest,
// but we know that each root is exposed by M, so after traversing
// the matching edge, we cannot have reached a root.
id = graph.phi(id);
id = attrs.phi(id);
retval.push_back(id);
}
return retval;
}
void collect_exposed_vertices(Graph & graph, std::stack<NodeId> & container)
void collect_exposed_vertices(GraphAttributes & attrs, std::stack<NodeId> & container)
{
std::stack<NodeId>().swap(container);
for(NodeId id = 0; id < graph.num_nodes(); id++)
for(NodeId id = 0; id < attrs.num_nodes(); id++)
{
if (graph.matched_neighbor(id) == id)
if (attrs.matched_neighbor(id) == id)
{
container.push(id);
graph.node(id).scanned = true;
attrs.scanned_[id] = true;
}
}
}
__attribute__((noinline))
void augment(Graph & graph, std::vector<NodeId> const & x_path, std::vector<NodeId> const & y_path,
void augment(GraphAttributes & attrs, std::vector<NodeId> const & x_path, std::vector<NodeId> const & y_path,
std::stack<NodeId> & outer_unvisited_nodes)
{
//std::cout << "Augment" << std::endl;
// Paths are disjoint -> augment
graph.node(x_path.front()).matched_neighbor = y_path.front();
graph.node(y_path.front()).matched_neighbor = x_path.front();
attrs.mu_[x_path.front()] = y_path.front();
attrs.mu_[y_path.front()] = x_path.front();
// TODO: put this into own method?
for(size_t i = 1; i < x_path.size(); i += 2)
{
graph.node(x_path[i]).matched_neighbor = x_path[i+1];
graph.node(x_path[i+1]).matched_neighbor = x_path[i];
attrs.mu_[x_path[i]] = x_path[i+1];
attrs.mu_[x_path[i+1]] = x_path[i];
}
for(size_t i = 1; i < y_path.size(); i += 2)
{
graph.node(y_path[i]).matched_neighbor = y_path[i+1];
graph.node(y_path[i+1]).matched_neighbor = y_path[i];
attrs.mu_[y_path[i]] = y_path[i+1];
attrs.mu_[y_path[i+1]] = x_path[i+1];
}
// Note that since this is tail-recursion, this will not generate
// new stack frames in OPT mode
graph.reset_forest();
collect_exposed_vertices(graph, outer_unvisited_nodes);
attrs.reset_forest();
collect_exposed_vertices(attrs, outer_unvisited_nodes);
}
__attribute__((noinline))
std::tuple<NodeId, size_type, size_type> find_blossom_root_id(Graph const & graph, std::vector<NodeId> const & x_path, std::vector<NodeId> const & y_path)
std::tuple<NodeId, size_type, size_type> find_blossom_root_id(GraphAttributes const & attrs, std::vector<NodeId> const & x_path, std::vector<NodeId> const & y_path)
{
size_t distance_from_x = x_path.size() - 1;
size_t distance_from_y = y_path.size() - 1;
@ -136,7 +135,7 @@ std::tuple<NodeId, size_type, size_type> find_blossom_root_id(Graph const & grap
--distance_from_y;
}
// found first vertex of x_path \cap y_path
while (graph.rho(x_path[distance_from_x]) != x_path[distance_from_x])
while (attrs.rho(x_path[distance_from_x]) != x_path[distance_from_x])
{
++distance_from_x;
++distance_from_y;
@ -146,7 +145,7 @@ std::tuple<NodeId, size_type, size_type> find_blossom_root_id(Graph const & grap
}
__attribute__((noinline))
void update_phi_along_blossom_paths(Graph & graph, std::vector<NodeId> const & x_path, std::vector<NodeId> const & y_path,
void update_phi_along_blossom_paths(GraphAttributes & attrs, std::vector<NodeId> const & x_path, std::vector<NodeId> const & y_path,
std::tuple<NodeId, size_type, size_type> const & blossom_root)
{
auto const [blossom_root_id, distance_from_x, distance_from_y] = blossom_root;
@ -154,46 +153,46 @@ void update_phi_along_blossom_paths(Graph & graph, std::vector<NodeId> const & x
// Update φ along the paths to encode the ear decomposition
for (size_t i = 1; i <= distance_from_x; i += 2)
{
if (graph.rho(graph.phi(x_path[i])) != blossom_root_id)
if (attrs.rho(attrs.phi(x_path[i])) != blossom_root_id)
{
graph.node(graph.phi(x_path[i])).phi = x_path[i];
attrs.phi_[attrs.phi(x_path[i])] = x_path[i];
}
}
for (size_t i = 1; i <= distance_from_y; i += 2)
{
if (graph.rho(graph.phi(y_path[i])) != blossom_root_id)
if (attrs.rho(attrs.phi(y_path[i])) != blossom_root_id)
{
graph.node(graph.phi(y_path[i])).phi = y_path[i];
attrs.phi_[attrs.phi(y_path[i])] = y_path[i];
}
}
// Link x and y
if (graph.rho(x_path.front()) != blossom_root_id)
if (attrs.rho(x_path.front()) != blossom_root_id)
{
graph.node(x_path.front()).phi = y_path.front();
attrs.phi_[x_path.front()] = y_path.front();
}
if (graph.rho(y_path.front()) != blossom_root_id)
if (attrs.rho(y_path.front()) != blossom_root_id)
{
graph.node(y_path.front()).phi = x_path.front();
attrs.phi_[y_path.front()] = x_path.front();
}
}
__attribute__((noinline))
void contract_rho(Graph & graph, NodeId blossom_root_id)
void contract_rho(GraphAttributes & attrs, NodeId blossom_root_id)
{
// Iterating over whole graph.
for (NodeId node_id = 0; node_id < graph.num_nodes(); ++node_id)
// Iterating over whole attrs.
for (NodeId node_id = 0; node_id < attrs.num_nodes(); ++node_id)
{
if (graph.rho(graph.rho(node_id)) == blossom_root_id)
if (attrs.rho(attrs.rho(node_id)) == blossom_root_id)
{
graph.node(node_id).rho = blossom_root_id;
attrs.rho_[node_id] = blossom_root_id;
}
}
}
__attribute__((noinline))
void update_rho(Graph & graph, std::vector<NodeId> const & x_path, std::vector<NodeId> const & y_path,
void update_rho(GraphAttributes & attrs, std::vector<NodeId> const & x_path, std::vector<NodeId> const & y_path,
std::tuple<NodeId, size_type, size_type> const & blossom_root_description,
std::stack<NodeId> & outer_unvisited_nodes)
{
@ -206,38 +205,38 @@ void update_rho(Graph & graph, std::vector<NodeId> const & x_path, std::vector<N
auto const [blossom_root_id, distance_from_x, distance_from_y] = blossom_root_description;
for (size_t i = 0; i <= distance_from_x; ++i)
{
graph.node(x_path[i]).rho = blossom_root_id;
if (not graph.node(x_path[i]).scanned)
attrs.rho_[x_path[i]] = blossom_root_id;
if (not attrs.scanned_[x_path[i]])
{
outer_unvisited_nodes.push(x_path[i]);
}
}
for (size_t i = 0; i <= distance_from_y; ++i)
{
graph.node(y_path[i]).rho = blossom_root_id;
if (not graph.node(y_path[i]).scanned)
attrs.rho_[y_path[i]] = blossom_root_id;
if (not attrs.scanned_[y_path[i]])
{
outer_unvisited_nodes.push(y_path[i]);
}
}
contract_rho(graph, blossom_root_id);
contract_rho(attrs, blossom_root_id);
}
__attribute__((noinline))
void contract_blossom(Graph & graph, std::vector<NodeId> const & x_path, std::vector<NodeId> const & y_path,
void contract_blossom(GraphAttributes & attrs, std::vector<NodeId> const & x_path, std::vector<NodeId> const & y_path,
std::stack<NodeId> & outer_unvisited_nodes)
{
//std::cout << "Contract blossom" << std::endl;
std::tuple<NodeId, size_type, size_type> const blossom_root_description = find_blossom_root_id(graph, x_path, y_path);
update_phi_along_blossom_paths(graph, x_path, y_path, blossom_root_description);
std::tuple<NodeId, size_type, size_type> const blossom_root_description = find_blossom_root_id(attrs, x_path, y_path);
update_phi_along_blossom_paths(attrs, x_path, y_path, blossom_root_description);
//check_integrity(graph);
update_rho(graph, x_path, y_path, blossom_root_description, outer_unvisited_nodes);
//check_integrity(attrs);
update_rho(attrs, x_path, y_path, blossom_root_description, outer_unvisited_nodes);
}
void maximum_matching_from_initial_matching(Graph & graph)
void maximum_matching_from_initial_matching(Graph const & graph, GraphAttributes & attrs)
{
graph.reset_forest();
attrs.reset_forest();
// Over the course of the algorithm, this will maintain all outer vertices
// that have not been scanned yet.
// Note that at the beginning, this is exactly the exposed edges.
@ -246,7 +245,7 @@ void maximum_matching_from_initial_matching(Graph & graph)
// When this stack runs out, then we know that all vertices marked 'scanned' have already been processed,
// but also all vertices not marked 'scanned' are not outer vertices, so we can in fact terminate.
std::stack<NodeId> outer_unvisited_nodes;
collect_exposed_vertices(graph, outer_unvisited_nodes);
collect_exposed_vertices(attrs, outer_unvisited_nodes);
while(not outer_unvisited_nodes.empty())
{
NodeId const id = outer_unvisited_nodes.top();
@ -255,48 +254,48 @@ void maximum_matching_from_initial_matching(Graph & graph)
{
//check_integrity(graph);
//std::cout << "Check passed" << std::endl;
if (graph.is_out_of_forest(neighbor_id))
if (attrs.is_out_of_forest(neighbor_id))
{
//std::cout << "Grow forest" << std::endl;
// Grow Forest
graph.node(neighbor_id).phi = id;
assert(graph.matched_neighbor(neighbor_id) != neighbor_id);
outer_unvisited_nodes.push(graph.matched_neighbor(neighbor_id));
attrs.phi_[neighbor_id] = id;
assert(attrs.matched_neighbor(neighbor_id) != neighbor_id);
outer_unvisited_nodes.push(attrs.matched_neighbor(neighbor_id));
}
else if (graph.is_outer(neighbor_id) and graph.rho(id) != graph.rho(neighbor_id))
else if (attrs.is_outer(neighbor_id) and attrs.rho(id) != attrs.rho(neighbor_id))
{
std::vector<NodeId> x_path = path_to_forest_root(graph, id);
std::vector<NodeId> y_path = path_to_forest_root(graph, neighbor_id);
std::vector<NodeId> x_path = path_to_forest_root(attrs, id);
std::vector<NodeId> y_path = path_to_forest_root(attrs, neighbor_id);
if (x_path.back() != y_path.back())
{
// paths are disjoint -> can augment
augment(graph, x_path, y_path, outer_unvisited_nodes);
augment(attrs, x_path, y_path, outer_unvisited_nodes);
break;
}
else
{
// Paths are not disjoint -> contract the new blossom
contract_blossom(graph, x_path, y_path, outer_unvisited_nodes);
contract_blossom(attrs, x_path, y_path, outer_unvisited_nodes);
}
}
}
graph.node(id).scanned = true;
attrs.scanned_[id] = true;
}
};
void find_greedy_matching(Graph & graph)
void find_greedy_matching(Graph const & graph, GraphAttributes & attrs)
{
graph.reset_matching();
attrs.reset_matching();
for(NodeId node_id = 0; node_id < graph.num_nodes(); ++node_id)
{
if (graph.matched_neighbor(node_id) == node_id) {
if (attrs.matched_neighbor(node_id) == node_id) {
for(NodeId const neighbor_id : graph.node(node_id).neighbors())
{
if(graph.matched_neighbor(neighbor_id) == neighbor_id)
if(attrs.matched_neighbor(neighbor_id) == neighbor_id)
{
graph.node(neighbor_id).matched_neighbor = node_id;
graph.node(node_id).matched_neighbor = neighbor_id;
attrs.mu_[neighbor_id] = node_id;
attrs.mu_[node_id] = neighbor_id;
break;
}
}
@ -305,17 +304,18 @@ void find_greedy_matching(Graph & graph)
}
Graph maximum_matching(Graph & graph) {
graph.reset_forest();
find_greedy_matching(graph);
check_integrity(graph);
maximum_matching_from_initial_matching(graph);
GraphAttributes attrs(graph.num_nodes());
attrs.reset_forest();
find_greedy_matching(graph, attrs);
check_integrity(attrs);
maximum_matching_from_initial_matching(graph, attrs);
ED::Graph matching = ED::Graph(graph.num_nodes());
for (NodeId id = 0; id < graph.num_nodes(); ++id)
{
if (graph.matched_neighbor(id) > id)
if (attrs.matched_neighbor(id) > id)
{
matching.add_edge(id, graph.matched_neighbor(id));
matching.add_edge(id, attrs.matched_neighbor(id));
}
}
return matching;

21
src/graph.cpp
View file

@ -176,26 +176,5 @@ std::ostream & operator<<(std::ostream & output, Graph const & graph)
return output;
}
void Graph::reset_forest()
{
NodeId cur_id = 0;
for(auto & node : _nodes) {
node.phi = cur_id;
node.rho = cur_id;
node.scanned = false;
// Note that we do not change the matching itself here
++cur_id;
}
}
void Graph::reset_matching()
{
NodeId cur_id = 0;
for(auto & node : _nodes)
{
node.matched_neighbor = cur_id;
++cur_id;
}
}
} // namespace ED

62
src/graph.hpp
View file

@ -91,12 +91,6 @@ public:
/** @return The array of ids of the neighbors of this node. **/
std::vector<NodeId> const & neighbors() const;
public:
NodeId matched_neighbor {invalid_node_id};
NodeId phi {invalid_node_id};
NodeId rho {invalid_node_id};
bool scanned;
private:
// This allows each Graph to access private members of this class,
// in our case the add_neighbor function
@ -170,21 +164,6 @@ public:
**/
friend std::ostream & operator<<(std::ostream & str, Graph const & graph);
NodeId matched_neighbor(NodeId const id) const;
NodeId phi(NodeId const id) const;
NodeId rho(NodeId const id) const;
bool is_outer(NodeId const id) const;
bool is_inner(NodeId const id) const;
bool is_out_of_forest(NodeId const id) const;
void reset_forest();
void reset_matching();
private:
std::vector<Node> _nodes;
size_type _num_edges;
@ -239,47 +218,6 @@ Node & Graph::node(NodeId const id)
//END: Inline section
inline
NodeId Graph::matched_neighbor(NodeId const id) const
{
assert(id <= num_nodes());
return _nodes[id].matched_neighbor;
}
inline
NodeId Graph::phi(const NodeId id) const
{
assert(id <= num_nodes());
return _nodes[id].phi;
}
inline
NodeId Graph::rho(const NodeId id) const
{
assert(id <= num_nodes());
return _nodes[id].rho;
}
inline
bool Graph::is_outer(NodeId const id) const {
return matched_neighbor(id) == id or \
phi(matched_neighbor(id)) != matched_neighbor(id);
}
inline
bool Graph::is_inner(NodeId const id) const
{
return phi(id) != id and \
phi(matched_neighbor(id)) == matched_neighbor(id);
}
inline
bool Graph::is_out_of_forest(const ED::NodeId id) const
{
return matched_neighbor(id) != id and \
phi(id) == id and \
phi(matched_neighbor(id)) == matched_neighbor(id);
}
} // namespace ED

112
src/graph_attributes.h Normal file
View file

@ -0,0 +1,112 @@
#ifndef GRAPH_ATTRIBUTES_H
#define GRAPH_ATTRIBUTES_H
#include <numeric>
#include "graph.hpp"
namespace ED
{
struct GraphAttributes
{
explicit GraphAttributes(NodeId num_nodes);
[[nodiscard]] NodeId num_nodes() const;
[[nodiscard]] NodeId matched_neighbor(NodeId id) const;
[[nodiscard]] NodeId phi(NodeId id) const;
[[nodiscard]] NodeId rho(NodeId id) const;
[[nodiscard]] bool is_outer(NodeId id) const;
[[nodiscard]] bool is_inner(NodeId id) const;
[[nodiscard]] bool is_out_of_forest(NodeId id) const;
void reset_forest();
void reset_matching();
std::vector<NodeId> phi_;
std::vector<NodeId> rho_;
std::vector<NodeId> mu_;
std::vector<bool> scanned_;
};
inline
GraphAttributes::GraphAttributes(const ED::NodeId num_nodes):
phi_(num_nodes),
rho_(num_nodes),
mu_(num_nodes),
scanned_(num_nodes)
{
}
inline
NodeId GraphAttributes::num_nodes() const
{
assert(phi_.size() == rho_.size());
assert(phi_.size() == mu_.size());
assert(phi_.size() == scanned_.size());
return phi_.size();
}
inline
NodeId GraphAttributes::matched_neighbor(NodeId const id) const
{
assert(id <= num_nodes());
return mu_[id];
}
inline
NodeId GraphAttributes::phi(const NodeId id) const
{
assert(id <= num_nodes());
return phi_[id];
}
inline
NodeId GraphAttributes::rho(const NodeId id) const
{
assert(id <= num_nodes());
return rho_[id];
}
inline
bool GraphAttributes::is_outer(NodeId const id) const {
return matched_neighbor(id) == id or \
phi(matched_neighbor(id)) != matched_neighbor(id);
}
inline
bool GraphAttributes::is_inner(NodeId const id) const
{
return phi(id) != id and \
phi(matched_neighbor(id)) == matched_neighbor(id);
}
inline
bool GraphAttributes::is_out_of_forest(const ED::NodeId id) const
{
return matched_neighbor(id) != id and \
phi(id) == id and \
phi(matched_neighbor(id)) == matched_neighbor(id);
}
inline
void GraphAttributes::reset_forest()
{
std::iota(phi_.begin(), phi_.end(), 0);
std::iota(rho_.begin(), rho_.end(), 0);
std::fill(scanned_.begin(), scanned_.end(), false);
}
inline
void GraphAttributes::reset_matching()
{
std::iota(mu_.begin(), mu_.end(), 0);
}
}
#endif //GRAPH_ATTRIBUTES_H