From c7ca55a2348dce3f8c05ae92fb3f344de8c1ebfb Mon Sep 17 00:00:00 2001 From: Krho Date: Mon, 23 Mar 2026 11:38:20 +0100 Subject: [PATCH 1/2] feat: implementation with hierholzer and naive matching for eulerian cycle --- graph-test.py | 36 ++++++++++ graph.py | 123 +++++++++++++++++++++++++++++++--- instances/farfalle.txt | 12 ++++ instances/line.txt | 6 ++ instances/triangle.txt | 7 ++ requirements.txt | 4 +- solver/chinesePostman-test.py | 51 ++++++++++++++ solver/chinesePostman.py | 29 ++++++++ solver/dijkstra-test.py | 29 ++++++++ solver/dijkstra.py | 15 +++++ solver/hierholzer-test.py | 48 +++++++++++++ solver/hierholzer.py | 80 ++++++++++++++++++++++ solver/naiveMatching-test.py | 53 +++++++++++++++ solver/naiveMatching.py | 52 ++++++++++++++ solver/solution.md | 57 ++++++++++++++++ 15 files changed, 592 insertions(+), 10 deletions(-) create mode 100644 graph-test.py create mode 100644 instances/farfalle.txt create mode 100644 instances/line.txt create mode 100644 instances/triangle.txt create mode 100644 solver/chinesePostman-test.py create mode 100644 solver/chinesePostman.py create mode 100644 solver/dijkstra-test.py create mode 100644 solver/dijkstra.py create mode 100644 solver/hierholzer-test.py create mode 100644 solver/hierholzer.py create mode 100644 solver/naiveMatching-test.py create mode 100644 solver/naiveMatching.py create mode 100644 solver/solution.md diff --git a/graph-test.py b/graph-test.py new file mode 100644 index 0000000..d5a0711 --- /dev/null +++ b/graph-test.py @@ -0,0 +1,36 @@ +from graph import Graph +from input import parse_file + + +# Refacto : have a fixture file +def read_graph(in_file: str) -> Graph: + vertices, edges = parse_file(in_file) + print(f"#E={len(edges)}, #V={len(vertices)}") + return Graph(vertices, edges) + + +triangle = read_graph("instances/triangle.txt") +line = read_graph("instances/line.txt") +paris = read_graph("instances/paris_map.txt") + +## End refacto needed + + +def odd(): + assert triangle.is_eulerian() + assert not line.is_eulerian() + assert triangle.odd_vertices() == [] + assert line.odd_vertices() == [0, 2] + assert not paris.is_eulerian() + + +def next(): + assert triangle.next(0, 1) == [1, 2] + assert triangle.next(1, 2) == [2, 0] # Next + assert triangle.next(0, 0) == [2, 2] # Reverse + assert line.next(1, 2) == [-1, None] # End of the line + + +def tests(): + odd() + next() diff --git a/graph.py b/graph.py index 5ce4d4a..e6c4d14 100644 --- a/graph.py +++ b/graph.py @@ -1,10 +1,15 @@ from __future__ import annotations +from typing import Union + import matplotlib.pyplot as plt from matplotlib.collections import LineCollection +from scipy.sparse import csr_array Coordinates = tuple[float, float] -Edge = tuple[int, int, int, Coordinates, Coordinates] +True_Edge = tuple[int, int, int, Coordinates, Coordinates] +Virtual_Edge = tuple[int, int, int, Coordinates, Coordinates, list[True_Edge]] +Edge = Union[None, True_Edge, Virtual_Edge] class Graph: @@ -25,19 +30,21 @@ def __init__( """ self.vertices = vertices self.edges = edges + self.real_edges = edges + self.adjency_computed = False def plot(self): """ Plot the graph. """ - weights = list(set(edge[2] for edge in self.edges)) + weights = list(set(edge[2] for edge in self.edges if edge is not None)) colors = plt.cm.get_cmap("viridis", len(weights)) _, ax = plt.subplots() for i, weight in enumerate(weights): lines = [ [edge[-2][::-1], edge[-1][::-1]] for edge in self.edges - if edge[2] == weight + if edge is not None and edge[2] == weight ] ax.add_collection( LineCollection( @@ -60,10 +67,108 @@ def display_paths(*, paths: list[list[Edge]]): ax = figure.add_subplot() for path_index, path in enumerate(paths): for edge_index, edge in enumerate(path): - ax.annotate( - str(edge_index), - xytext=edge[3], - xy=edge[4], - arrowprops=dict(arrowstyle="->", color=colors(path_index)), - ) + if edge is None: + continue + if edge is True_Edge: + ax.annotate( + str(edge_index), + xytext=edge[3], + xy=edge[4], + width=1, + arrowprops=dict(arrowstyle="->", color=colors(path_index)), + ) + if edge is Virtual_Edge: + # TODO : more elegant way to display edges travelled twice + ax.annotate( + f"virtual {edge_index}", + xytext=edge[3], + xy=edge[4], + width=2, + arrowprops=dict(arrowstyle="->", color=colors(path_index)), + ) plt.show() + + def compute_adjency(self): + vertex_adjency = [0 for vertex in self.vertices] + for edge in self.edges: + if edge is not None: + [v1, v2] = edge[:2] + vertex_adjency[v1] += 1 + vertex_adjency[v2] += 1 + self.vertex_adjency = vertex_adjency + self.vertices_odd = [ + i for i in range(len(self.vertices)) if vertex_adjency[i] % 2 == 1 + ] + self.adjency_computed = True + self.is_eulerian_cache = len(self.vertices_odd) == 0 + + def is_eulerian(self): + if not self.adjency_computed: + self.compute_adjency() + return self.is_eulerian_cache + + def odd_vertices(self): + if not self.adjency_computed: + self.compute_adjency() + return self.vertices_odd + + def next(self, input_edge_index: int, vertex_index: int): + for i, edge in enumerate(self.edges): + if i == input_edge_index or edge is None: + continue # Actually not necessary for Hierholzer since we remove the edge before computing the next + if edge[0] == vertex_index: + return [i, edge[1]] + if edge[1] == vertex_index: + return [i, edge[0]] + return [-1, None] + + # TODO instead of None, add a boolean to mark the edge as visited + # no need to clone this way + # prevents all the "if edge is None" in the code + def remove_edge(self, index): + if not self.adjency_computed: + self.compute_adjency() + edge = self.edges[index] + self.vertex_adjency[edge[0]] -= 1 + self.vertex_adjency[edge[1]] -= 1 + self.edges[index] = None + + def add_virtual_edge( + self, v1_index: int, v2_index: int, weight: int, path: list[True_Edge] + ): + [v1, v2] = [self.vertices[v1_index], self.vertices[v2_index]] + new_edge: Virtual_Edge = (v1_index, v2_index, weight, v1, v2, path) + self.edges.append(new_edge) + self.adjency_computed = False + + def has_unvisited_edges(self): + if not self.adjency_computed: + return len(self.edges) > 0 + for e in self.edges: + if e is not None: + return True + return False + + def is_hub(self, vertex): + if not self.adjency_computed: + self.compute_adjency() + return self.vertex_adjency[vertex] > 0 + + # Probably a standard method somewhere + def matrix(self): + n = len(self.vertices) + matrix = [[0 for i in range(n)] for i in range(n)] + for edge in self.edges: + if edge is not None: + [i, j, weight] = [edge[0], edge[1], edge[2]] + matrix[i][j] = weight + matrix[j][i] = weight + return csr_array(matrix) + + def clone(self): + g = Graph( + [vertex for vertex in self.vertices], [edge for edge in self.real_edges] + ) + g.real_edges = g.edges + g.edges = [edge for edge in self.edges] + return g diff --git a/instances/farfalle.txt b/instances/farfalle.txt new file mode 100644 index 0000000..70c88cb --- /dev/null +++ b/instances/farfalle.txt @@ -0,0 +1,12 @@ +5 6 +0 0 +0 1 +1 1 +1 2 +2 2 +0 1 2 +1 2 1 +2 0 1 +2 3 1 +3 4 1 +4 2 1 diff --git a/instances/line.txt b/instances/line.txt new file mode 100644 index 0000000..e4bbda8 --- /dev/null +++ b/instances/line.txt @@ -0,0 +1,6 @@ +3 2 +0.5 0 +0.5 1 +0.5 2 +0 1 1 +1 2 1 diff --git a/instances/triangle.txt b/instances/triangle.txt new file mode 100644 index 0000000..271717a --- /dev/null +++ b/instances/triangle.txt @@ -0,0 +1,7 @@ +3 3 +0 0 +0 1 +1 1 +0 1 2 +1 2 1 +2 0 1 diff --git a/requirements.txt b/requirements.txt index 4b43f7e..fee5a8a 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1 +1,3 @@ -matplotlib \ No newline at end of file +matplotlib +pytest +scipy diff --git a/solver/chinesePostman-test.py b/solver/chinesePostman-test.py new file mode 100644 index 0000000..f2b74d0 --- /dev/null +++ b/solver/chinesePostman-test.py @@ -0,0 +1,51 @@ +import datetime + +from chinesePostman import chinese_postman_solver + +from graph import Graph +from input import parse_file + + +# Refacto : have a fixture file +def read_graph(in_file: str) -> Graph: + vertices, edges = parse_file(in_file) + print(f"#E={len(edges)}, #V={len(vertices)}") + return Graph(vertices, edges) + + +line = read_graph("instances/line.txt") +triangle = read_graph("instances/triangle.txt") +hard_to_choose = read_graph("instances/hard_to_choose.txt") +islands = read_graph("instances/islands.txt") +paris = read_graph("instances/paris_map.txt") + +# End Refacto + + +def integration(): + path = chinese_postman_solver(line) + assert path[0] == 4 # weight + assert len(path[1]) == 3 # two real one virtual + path = chinese_postman_solver(triangle) + assert path[0] == 4 + assert len(path[1]) == 3 + + +def performance(graph): + begin = datetime.datetime.now() + path = chinese_postman_solver(graph) + duration = datetime.datetime.now() - begin + print(f"found path of weight {path[0]} in {duration}s") + + +def performances(): + performance(hard_to_choose) + """ both these instances fail with indexes bound when reconstructing Dijkstra paths + performance(islands) + performance(paris) + """ + + +def tests(): + integration() + performances() diff --git a/solver/chinesePostman.py b/solver/chinesePostman.py new file mode 100644 index 0000000..349253b --- /dev/null +++ b/solver/chinesePostman.py @@ -0,0 +1,29 @@ +from hierholzer import hierholzer +from naiveMatching import naive_matching + +from graph import Graph + + +def chinese_postman_solver_generic(g: Graph, eulerian_cycle, matching): + if not g.is_eulerian(): + print("no Eulerian cycle in the graph ; adding missing edges", flush=True) + g = matching(g) + print("computing Eulerian cycle on enriched graph", flush=True) + else: + print("Graph is Eulerian, computing cycle directly") + cycle = eulerian_cycle(g) + return cycle_to_path(cycle) + + +def cycle_to_path(cycle): + # TODO + # if there are virtual edges, remove the one of maximum weight + # present the path in a readable way + return cycle + + +def chinese_postman_solver(g: Graph): + # TODO + # Refacto to inject smallest path strategy as well (here choice to use Dijkstra is hidden) + # Implement and use Edmond's blossom algorithm instead of naive matching : will increase computation time but will give optimal solution + return chinese_postman_solver_generic(g, hierholzer, naive_matching) diff --git a/solver/dijkstra-test.py b/solver/dijkstra-test.py new file mode 100644 index 0000000..d057624 --- /dev/null +++ b/solver/dijkstra-test.py @@ -0,0 +1,29 @@ +import datetime + +from dijkstra import compute_dijkstra + +from graph import Graph +from input import parse_file + + +# Refacto : have a fixture file +def read_graph(in_file: str) -> Graph: + vertices, edges = parse_file(in_file) + print(f"#E={len(edges)}, #V={len(vertices)}") + return Graph(vertices, edges) + + +paris = read_graph("instances/paris_map.txt") + +# End Refacto + + +def perf(): + begin = datetime.datetime.now() + compute_dijkstra(paris) + duration = datetime.datetime.now() - begin + print(f"computed dijkstra in {duration}") + + +def tests(): + perf() diff --git a/solver/dijkstra.py b/solver/dijkstra.py new file mode 100644 index 0000000..294b26c --- /dev/null +++ b/solver/dijkstra.py @@ -0,0 +1,15 @@ +from scipy.sparse.csgraph import dijkstra + +from graph import Graph + + +# TODO make it a class so it can be a service called by every method used for perfect-matching +def compute_dijkstra(g: Graph): + graph = g.matrix() + dist_matrix, predecessors = dijkstra( + csgraph=graph, + directed=False, + return_predecessors=True, + indices=g.odd_vertices(), + ) + return [dist_matrix, predecessors] diff --git a/solver/hierholzer-test.py b/solver/hierholzer-test.py new file mode 100644 index 0000000..4633ba5 --- /dev/null +++ b/solver/hierholzer-test.py @@ -0,0 +1,48 @@ +from hierholzer import cycle, hierholzer + +from graph import Graph +from input import parse_file + + +# Refacto : have a fixture file +def read_graph(in_file: str) -> Graph: + vertices, edges = parse_file(in_file) + print(f"#E={len(edges)}, #V={len(vertices)}") + return Graph(vertices, edges) + + +triangle = read_graph("instances/triangle.txt") +line = read_graph("instances/line.txt") +farfalle = read_graph("instances/farfalle.txt") + +# End todo refacto + + +def cycle_tests(): + assert cycle(triangle.clone(), 0, 0) == [ + 4, + [triangle.edges[0], triangle.edges[2], triangle.edges[1]], + ] + + +def hierholzer_tests(): + assert hierholzer(triangle.clone()) == [ + 4, + [triangle.edges[0], triangle.edges[2], triangle.edges[1]], + ] + assert hierholzer(farfalle.clone()) == [ + 7, + [ + farfalle.edges[0], + farfalle.edges[2], + farfalle.edges[3], + farfalle.edges[5], + farfalle.edges[4], + farfalle.edges[1], + ], + ] + + +def tests(): + cycle_tests() + hierholzer_tests() diff --git a/solver/hierholzer.py b/solver/hierholzer.py new file mode 100644 index 0000000..6f7a813 --- /dev/null +++ b/solver/hierholzer.py @@ -0,0 +1,80 @@ +from typing import Union + +from graph import Graph + + +# TODO : refacto to make it simpler +def cycle(g: Graph, initial_vertex: int, initial_edge_index: int): + i = initial_edge_index + edge = g.edges[i] + if edge is None: + return [0, []] + vertex = initial_vertex + weight = 0 + path = [] + # tour until either dead end or coming back to starting point + while ( + vertex is not None + and edge is not None + and (vertex != initial_vertex or len(path) <= 0) + ): + weight += edge[2] + path.append(edge) + g.remove_edge(i) + [j, next_vertex] = g.next(i, vertex) + i = j + edge = g.edges[i] + vertex = next_vertex + return [weight, path] + + +# TODO : linear complexity, find a way to have constant complexity if too much time spent here +def insert(host: list, to_add: list): + if len(to_add) == 0: + return host + if len(host) == 0: + return to_add + vertex = to_add[0][0] + for [i, edge] in enumerate(host): + if edge[1] == vertex: + return host[:i] + to_add + host[i:] + # found no place to insert + lisible_to_add = [edge[:2] for edge in to_add] + lisible_host = [edge[:2] for edge in host] + raise Exception(f"Impossible to insert {lisible_to_add} into {lisible_host}") + + +# TODO improve the complexity with better data structure +def next_start(g: Graph, path=[]): + if len(path) == 0: + if len(g.edges) < 1 or g.edges[0] is None: + raise Exception("Impossible to compute Eulerian cycle") + return [0, 0] + for p in path: + vertex = p[1] + if g.is_hub(vertex): + for [i, edge] in enumerate(g.edges): + if edge is None: + continue + if edge[0] == vertex: + return [vertex, i] + raise Exception("Impossible to compute Eulerian cycle") + + +# TODO : if majority of time spent on this part of code +# improve data structure to have constant time for each operation +def hierholzer(input_graph: Graph) -> Union[list, None]: + g = input_graph.clone() + if len(g.edges) < 1 or g.edges[0] is None: + raise Exception("Impossible to compute Eulerian cycle") + [vertex, edge] = next_start(g) + [weight, path] = cycle(g, vertex, edge) + nb_tours = 1 + while g.has_unvisited_edges() and nb_tours < len(g.edges): # safety to break + [vertex, edge] = next_start(g, path) + current_tour = cycle(g, vertex, edge) + nb_tours += 1 + path = insert(path, current_tour[1]) + weight += current_tour[0] + print(f"{nb_tours} tours", flush=True) + return [weight, path] diff --git a/solver/naiveMatching-test.py b/solver/naiveMatching-test.py new file mode 100644 index 0000000..a865483 --- /dev/null +++ b/solver/naiveMatching-test.py @@ -0,0 +1,53 @@ +import datetime + +from naiveMatching import naive_matching, vertices_to_path + +from graph import Graph +from input import parse_file + + +# Refacto : have a fixture file +def read_graph(in_file: str) -> Graph: + vertices, edges = parse_file(in_file) + print(f"#E={len(edges)}, #V={len(vertices)}") + return Graph(vertices, edges) + + +line = read_graph("instances/line.txt") +hard_to_choose = read_graph("instances/hard_to_choose.txt") +paris = read_graph("instances/paris_map.txt") + + +def to_path(): + assert vertices_to_path(line, [0, 1, 2]) == [line.edges[0], line.edges[1]] + assert vertices_to_path(line, [2, 1, 0]) == [line.edges[1], line.edges[0]] + assert vertices_to_path(line, [0, 2]) == [] + + +def naive(): + naive_matching(line) + assert line.is_eulerian() + # we only test the added edge + assert len(line.edges) == 3 + expected_edge = ( + 0, + 2, + 2, + (0.5, 0), + (0.5, 2), + [line.edges[1], line.edges[0]], + ) + assert line.edges[2] == expected_edge + + +def perf(graph): + begin = datetime.datetime.now() + naive_matching(graph) + duration = datetime.datetime.now() - begin + print(f"found matching in {duration}") + + +def tests(): + to_path() + naive() + perf(hard_to_choose) diff --git a/solver/naiveMatching.py b/solver/naiveMatching.py new file mode 100644 index 0000000..be73fb8 --- /dev/null +++ b/solver/naiveMatching.py @@ -0,0 +1,52 @@ +import math + +from dijkstra import compute_dijkstra + +from graph import Graph + + +def vertices_to_path(g: Graph, vertices): + edges = [] + for i in range(len(vertices) - 1): + # TODO complexity is linear when it should be constant because the data structures don't work well with each other + for edge in g.real_edges: + if edge is None: + continue + # TODO create an explicit method for this + if (edge[0] == vertices[i] and edge[1] == vertices[i + 1]) or ( + edge[0] == vertices[i + 1] and edge[1] == vertices[i] + ): + edges.append(edge) + break + return edges + + +# TODO : refacto separate recursive and global method +def compute_path(g: Graph, i: int, j: int, predecessors, vertices=[]): + vertices.append(j) + if i == j: + return vertices_to_path(g, vertices) + previous = predecessors[i][j] + if previous < 0: + return None + return compute_path(g, i, previous, predecessors, vertices) + + +# simplest way to make graph eulerian : pair each odd vertex with another one +def naive_matching(g: Graph): + odds = g.odd_vertices() + print("odd vertices computed") + if len(odds) % 2 != 0: + raise Exception("Something is wrong") + [dist_matrix, predecessors] = compute_dijkstra(g) + print("Dijkstra computed") + for k in range(len(odds) // 2): + [i, j] = odds[2 * k : 2 * k + 2] + weight = math.floor(dist_matrix[i][j]) + path = compute_path(g, i, j, predecessors) + g.add_virtual_edge(i, j, weight, path) + print( + f"added {len(g.edges) - len(g.real_edges)} to a graph with {len(odds)} odd vertices", + flush=True, + ) + return g diff --git a/solver/solution.md b/solver/solution.md new file mode 100644 index 0000000..974f125 --- /dev/null +++ b/solver/solution.md @@ -0,0 +1,57 @@ +# Problem + +## References +- https://en.wikipedia.org/wiki/Chinese_postman_problem +-- Kwan, Mei-ko (1960), "Graphic programming using odd or even points", Acta Mathematica Sinica (in Chinese), 10: 263–266, MR 0162630. Translated in Chinese Mathematics 1: 273–277, 1962. +- https://en.wikipedia.org/wiki/Eulerian_path#Hierholzer's_algorithm +-- Fleischner, Herbert (1991), "X.1 Algorithms for Eulerian Trails", Eulerian Graphs and Related Topics: Part 1, Volume 2, Annals of Discrete Mathematics, vol. 50, Elsevier, pp. X.1–13, ISBN 978-0-444-89110-5. +- https://en.wikipedia.org/wiki/Blossom_algorithm +-- Edmonds, Jack (1965). "Paths, trees, and flowers". Can. J. Math. 17: 449–467. doi:10.4153/CJM-1965-045-4. + +## Global algorithm + +Finding the path of smallest weight that covers all edges of an undirected graph is a **Chinese postman problem**. It has an exact solution with polynomial complexity 🎉 + +When a solution exists where each edge is travelled exactly once, we have an Eulerian graph, and the solution is called an Eulerian cycle. On an Eulerian graph, each node has an even number of vertices. + +When there is no Eulerian graph, the best solution is to add virtual vertices to make the graph Eulerian. Each virtual vertex is an existing path in the graph. + +To minimize the additional weight : ++ We only add virtual vertices for edges which need them + + only betwen edges with even vertices + + half as many vertices as even edges ++ minimum weight for each vertex : this is a shortest path problem, solvable with a Dijkstra algorithm ++ minimum total weight of added vertices : minimum weight perfect matching problem, solvable with an Edmond's blossom algorithm + +Once we have an Eulerian graph, finding a Eulerian cycle can be done with Hierholzer's algorithm. + +Since we want a path and not a cycle, we can remove the virtual vertex of max weight from the solution, its edges becoming the start and the finish of the path. + +# Implementation choices + +I wanted to have code that runs: since Dijkstra is already available in scipy and the global Chinese postman algorithm is simple, I chose to +- **implement Hierholzer's algorithm** : it was core to the solution, and other approachs such as depth first travel looked like they would sacrifice a lot in term of weight without being simpler to write from scratch +- not implement Edmond's blossom algorithm since it's complex and a naive approach, even if real bad weight-wise, would be really easy to write +- keep the existing graph data structure as much as possible to build on what already exists +- added simple instances for tests + +# Potential next steps + +## Debug +Paris and islands instancs are failing at the naive_matching step + +## Performance +- measure time spent in each algorithm to taylor data structure specifically for it ; try with an adjency matrix + +## Solution quality +- Use Edmond's blossom algorithm instead of naive matching to minimize weight + +## Code quality +- Have a way to call `edge.vertex[0]` instead of `edge[0]` and `edge.weight()` instead of `edge[2]` +- better handling of errors / None value / empty result +- stronger typing (once performances are stabilized) +- naming consistency (sometimes `vertex` is an `int`, sometimes it's a `list[Coordinates]`) +- organise solver in sub folders for each problem (perfect matching, shortest path, eulerian path, and chinese postman) +- consistency between case (camel, kebab) +- consistency between functional and object oriented style +- create and use a decorator to measure and log time spent in each method From 3a6795c4133a6afe4ef3781340764b534d4c4c2f Mon Sep 17 00:00:00 2001 From: Krho Date: Mon, 23 Mar 2026 11:38:47 +0100 Subject: [PATCH 2/2] fix: updated .gitignore --- .gitignore | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/.gitignore b/.gitignore index f9915ab..5c5e489 100644 --- a/.gitignore +++ b/.gitignore @@ -1,2 +1,3 @@ __pycache__ -in/ \ No newline at end of file +in/ +.pytest_cache/