From 78ffdb7df416d24cbc1f6ff6d93c68c921afc19c Mon Sep 17 00:00:00 2001 From: Peter Corke Date: Mon, 12 Oct 2026 18:35:37 +1000 Subject: [PATCH 1/3] feat: Graph.within, Graph.near and Graph.subgraph within(box) lists the vertices inside an n-dimensional box given as a 2 x n array of lower and upper corner, near(coord, radius) those within a distance of a point by the graph's metric with their distances, nearest first, and subgraph(vertices) is the graph induced by a set of vertices, a new graph of the same kind in which every edge with an end outside is dropped. subgraph reuses the flat description that pickling uses, so classes, names, coordinates, costs, user attributes, metric and the order of vertices and edges are kept, parallel edges and loops included, and the result shares nothing with the original. Co-Authored-By: Claude Sonnet 5.5 --- docs/source/analysis.rst | 13 ++ src/pgraph/PGraph.py | 156 +++++++++++++++++++ tests/test_subgraph.py | 317 +++++++++++++++++++++++++++++++++++++++ 3 files changed, 486 insertions(+) create mode 100644 tests/test_subgraph.py diff --git a/docs/source/analysis.rst b/docs/source/analysis.rst index 10a157e0..a83d9879 100644 --- a/docs/source/analysis.rst +++ b/docs/source/analysis.rst @@ -49,10 +49,20 @@ Structure and connectivity Geometry and consistency ^^^^^^^^^^^^^^^^^^^^^^^^^ +For a graph whose vertices have coordinates, +:meth:`~pgraph.PGraph.Graph.closest` finds the nearest vertex to a point, +:meth:`~pgraph.PGraph.Graph.near` the vertices within a distance of a point and +:meth:`~pgraph.PGraph.Graph.within` the vertices inside a box. Passing the +result of the last to :meth:`~pgraph.PGraph.Graph.subgraph` gives the part of +the graph in a region, without the edges that cross its border. + .. autosummary:: :nosignatures: ~pgraph.PGraph.Graph.closest + ~pgraph.PGraph.Graph.near + ~pgraph.PGraph.Graph.within + ~pgraph.PGraph.Graph.subgraph ~pgraph.PGraph.Graph.validate The matrices of a graph (:meth:`~pgraph.PGraph.Graph.adjacency`, @@ -154,6 +164,9 @@ Details .. automethod:: pgraph.PGraph.Graph.connectivity .. automethod:: pgraph.PGraph.Graph.average_degree .. automethod:: pgraph.PGraph.Graph.closest +.. automethod:: pgraph.PGraph.Graph.near +.. automethod:: pgraph.PGraph.Graph.within +.. automethod:: pgraph.PGraph.Graph.subgraph .. automethod:: pgraph.PGraph.Graph.validate .. autoproperty:: pgraph.PGraph.Vertex.name diff --git a/src/pgraph/PGraph.py b/src/pgraph/PGraph.py index c2e2994a..44a40ffd 100644 --- a/src/pgraph/PGraph.py +++ b/src/pgraph/PGraph.py @@ -890,6 +890,70 @@ def copy(self) -> Graph: """ return copy.deepcopy(self) + def subgraph(self, vertices: Iterable[Vertex | str]) -> Graph: + """ + Subgraph induced by a set of vertices + + :param vertices: vertices of this graph, or their names + :type vertices: iterable of Vertex subclass or str + :raises ValueError: a vertex is not in this graph + :raises KeyError: a name is not that of a vertex in this graph + :raises TypeError: an item is neither a vertex nor a string + :return: a new graph of the same kind + :rtype: Graph + + The new graph has the given vertices and every edge of this graph that + has *both* its end vertices among them. An edge that has one end + outside is removed, so the edges that cross the border of a region + are dropped. Parallel edges and loops are kept. + + Like :meth:`copy` the new graph is independent of this one: its + vertices and edges are new objects that keep their class, name, + coordinate, cost and any attributes that were added, and the graph + keeps its metric, heuristic and vertex and edge classes. They are in + the order of this graph, however ``vertices`` was ordered. Names are + kept, so the vertex ``g[name]`` has its counterpart ``h[name]``. + + .. runblock:: pycon + + >>> from pgraph import UGraph + >>> g = UGraph() + >>> v = [g.add_vertex(coord=[i, 0], name=f"v{i}") for i in range(4)] + >>> for i in range(3): + ... g.add_edge(v[i], v[i + 1]) + >>> h = g.subgraph([v[1], v[2]]) + >>> print(h) + >>> h.edges() + + It is a single pass over the graph. The vertices might be a component, + the vertices that satisfy a test, or those in a region of an embedded + graph from :meth:`within`. + + :seealso: :meth:`copy` :meth:`within` :meth:`component` + """ + wanted = set() + for v in vertices: + v = self._resolve_vertex(v, "vertex") + if v._graph is not self: + raise ValueError(f"vertex {v.name!r} is not in this graph") + wanted.add(id(v)) + + state = self.__getstate__() + flat = state["_flat"] + keep = [i for i, v in enumerate(self._vertexlist) if id(v) in wanted] + renumber = {old: new for new, old in enumerate(keep)} + flat["vertices"] = [flat["vertices"][i] for i in keep] + flat["edges"] = [ + (cls, renumber[i], renumber[j], attrs) + for cls, i, j, attrs in flat["edges"] + if i in renumber and j in renumber + ] + + h = type(self).__new__(type(self)) + h.__setstate__(copy.deepcopy(state)) + h._connectivitychange = True # the components are those of the new graph + return h + # pickling and copying #: the graph's own containers, which are saved as a flat description @@ -1844,6 +1908,98 @@ def closest(self, coord: ArrayLike) -> tuple[Vertex | None, float]: return min_which, min_dist + def within(self, box: ArrayLike) -> list[Vertex]: + """ + Vertices inside a box + + :param box: lower and upper corners of the box + :type box: ndarray(2,n) + :raises ValueError: ``box`` is not 2 x n, its upper corner is below its + lower corner in some dimension, or the coordinate of a vertex does + not have n elements + :return: the vertices inside the box, in the order of the graph + :rtype: list of Vertex subclass + + The first row of ``box`` is the lower corner, the second row the upper + corner, so in two dimensions it is ``[[xmin, ymin], [xmax, ymax]]``. + The box is closed: a vertex on its boundary is inside. Vertices + without a coordinate are skipped, as they are by :meth:`closest`. + + It is a single pass over the vertices. Together with + :meth:`subgraph` it selects the part of an embedded graph in a region + + .. runblock:: pycon + + >>> from pgraph import UGraph + >>> g = UGraph() + >>> for i in range(5): + ... g.add_vertex(coord=[i, i], name=f"v{i}") + >>> g.within([[1, 1], [3, 3]]) + >>> g.subgraph(g.within([[1, 1], [3, 3]])) + + :seealso: :meth:`near` :meth:`closest` :meth:`subgraph` + """ + box = np.asarray(box, dtype=float) + if box.ndim != 2 or box.shape[0] != 2: + raise ValueError(f"box must be 2 x n, not {box.shape}") + lower, upper = box + if np.any(upper < lower): + raise ValueError("the upper corner of the box is below the lower corner") + + found = [] + for vertex in self: + if vertex.coord is None: + continue + coord = np.asarray(vertex.coord, dtype=float) + if coord.shape != lower.shape: + raise ValueError( + f"vertex {vertex.name!r} has a coordinate of {coord.size} " + f"elements, the box has {lower.size}" + ) + if np.all(coord >= lower) and np.all(coord <= upper): + found.append(vertex) + return found + + def near(self, coord: ArrayLike, radius: float) -> list[tuple[Vertex, float]]: + """ + Vertices within a distance of a point + + :param coord: coordinates of a point + :type coord: ndarray(n) + :param radius: greatest distance, inclusive + :raises ValueError: ``radius`` is negative + :return: the vertices within ``radius`` of the point with their + distances, nearest first + :rtype: list of (Vertex subclass, float) + + Distance is computed according to the graph's metric, as it is by + :meth:`closest`, of which this is the generalization to every vertex + within a radius. Vertices without a coordinate are skipped, and + vertices at equal distance are in the order of the graph. + + .. runblock:: pycon + + >>> from pgraph import UGraph + >>> g = UGraph() + >>> for i in range(5): + ... g.add_vertex(coord=[i, 0], name=f"v{i}") + >>> [(v.name, d) for v, d in g.near([2.2, 0], 1.5)] + + :seealso: :meth:`closest` :meth:`within` :meth:`metric` + """ + if radius < 0: + raise ValueError(f"radius must not be negative, got {radius!r}") + + found = [] + for vertex in self: + if vertex.coord is None: + continue + d = float(self.metric(vertex.coord - coord)) + if d <= radius: + found.append((vertex, d)) + found.sort(key=lambda pair: pair[1]) + return found + def edges(self) -> list[Edge]: """ All edges in the graph, in the order they were added diff --git a/tests/test_subgraph.py b/tests/test_subgraph.py new file mode 100644 index 00000000..f9d5e931 --- /dev/null +++ b/tests/test_subgraph.py @@ -0,0 +1,317 @@ +""" +``Graph.within`` (the vertices in a box), ``Graph.near`` (the vertices within +a distance of a point) and ``Graph.subgraph`` (the graph induced by a set of +vertices). + +Expectations are written out by hand. +""" + +import pickle # nosec B403 - only used to round-trip graphs made in this file +import unittest + +import numpy as np + +from pgraph import DGraph, Edge, UGraph, Vertex + + +def line(cls=UGraph, n=5, **kwargs): + """Vertices v0 .. v(n-1) at (0, 0), (1, 0), ... with no edges""" + g = cls(dim=2, **kwargs) + return g, [g.add_vertex([float(i), 0.0], name=f"v{i}") for i in range(n)] + + +def names(vertices): + return [v.name for v in vertices] + + +class TestWithin(unittest.TestCase): + def test_closed_box(self): + g, _ = line() + # the corners and edges of the box are inside + self.assertEqual(names(g.within([[1, -1], [3, 1]])), ["v1", "v2", "v3"]) + self.assertEqual(names(g.within([[1, 0], [3, 0]])), ["v1", "v2", "v3"]) + + def test_nothing_inside(self): + g, _ = line() + self.assertEqual(g.within([[1.2, -1], [1.8, 1]]), []) + self.assertEqual(g.within([[0, 1], [4, 2]]), []) + + def test_a_box_of_no_size(self): + g, _ = line() + self.assertEqual(names(g.within([[2, 0], [2, 0]])), ["v2"]) + + def test_graph_order_not_distance_order(self): + g = UGraph(dim=2) + for name, x in [("far", 9.0), ("near", 1.0), ("mid", 5.0)]: + g.add_vertex([x, 0.0], name=name) + self.assertEqual(names(g.within([[0, 0], [10, 0]])), ["far", "near", "mid"]) + + def test_three_dimensions(self): + g = UGraph(dim=3) + for name, c in [("a", (0, 0, 0)), ("b", (1, 1, 1)), ("c", (1, 1, 5))]: + g.add_vertex(list(c), name=name) + self.assertEqual(names(g.within([[0, 0, 0], [2, 2, 2]])), ["a", "b"]) + self.assertEqual(names(g.within([[0.5, 0.5, 0.5], [2, 2, 6]])), ["b", "c"]) + + def test_vertices_without_a_coordinate_are_skipped(self): + g = UGraph() + g.add_vertex(name="nowhere") + g.add_vertex([1.0, 1.0], name="here") + self.assertEqual(names(g.within([[0, 0], [2, 2]])), ["here"]) + + def test_the_box_must_be_two_by_n(self): + g, _ = line() + for box in ([0, 0, 1, 1], [[0, 0], [1, 1], [2, 2]], [[0], [1], [2]]): + with self.subTest(box=box): + with self.assertRaises(ValueError): + g.within(box) + + def test_the_upper_corner_must_not_be_below_the_lower(self): + g, _ = line() + with self.assertRaises(ValueError): + g.within([[3, 0], [1, 0]]) + + def test_wrong_dimension(self): + g, _ = line() # 2-D coordinates + with self.assertRaises(ValueError): + g.within([[0, 0, 0], [1, 1, 1]]) + + def test_does_not_change_the_graph(self): + g, vs = line() + g.within([[0, 0], [2, 0]]) + self.assertEqual(g.n, 5) + self.assertEqual(names(g), ["v0", "v1", "v2", "v3", "v4"]) + + +class TestNear(unittest.TestCase): + def test_sorted_by_distance(self): + g, _ = line() + found = g.near([2.2, 0], 1.5) + self.assertEqual(names(v for v, _ in found), ["v2", "v3", "v1"]) + np.testing.assert_allclose([d for _, d in found], [0.2, 0.8, 1.2]) + + def test_the_radius_is_inclusive(self): + g, _ = line() + self.assertEqual(names(v for v, _ in g.near([0, 0], 2)), ["v0", "v1", "v2"]) + + def test_nothing_near(self): + g, _ = line() + self.assertEqual(g.near([10, 10], 1), []) + + def test_equal_distances_in_graph_order(self): + g = UGraph(dim=2) + for name, c in [("east", (1, 0)), ("west", (-1, 0)), ("north", (0, 1))]: + g.add_vertex(list(c), name=name) + found = g.near([0, 0], 1) + self.assertEqual(names(v for v, _ in found), ["east", "west", "north"]) + + def test_the_closest_vertex_is_the_first(self): + g, _ = line() + near = g.near([2.2, 0], 10) + closest, distance = g.closest([2.2, 0]) + self.assertIs(near[0][0], closest) + self.assertAlmostEqual(near[0][1], distance) + + def test_uses_the_metric_of_the_graph(self): + for metric, expected in (("L2", ["a", "b"]), ("L1", ["a"])): + with self.subTest(metric=metric): + g = UGraph(dim=2, metric=metric) + g.add_vertex([0.0, 0.0], name="a") + g.add_vertex([1.0, 1.0], name="b") # 1.41 by L2, 2 by L1 + self.assertEqual(names(v for v, _ in g.near([0, 0], 1.5)), expected) + + def test_vertices_without_a_coordinate_are_skipped(self): + g = UGraph() + g.add_vertex(name="nowhere") + g.add_vertex([0.0, 0.0], name="here") + self.assertEqual(names(v for v, _ in g.near([0, 0], 1)), ["here"]) + + def test_negative_radius(self): + g, _ = line() + with self.assertRaises(ValueError): + g.near([0, 0], -1) + + +class MyVertex(Vertex): + pass + + +class MyEdge(Edge): + pass + + +class _Subgraph: + cls: type + + def graph(self): + """v0-v1, v1=v2 (two edges, costs 9 and 4), a loop on v2, v2-v3, v3-v4""" + g, v = line(self.cls) + g.add_edge(v[0], v[1]) + g.add_edge(v[1], v[2], cost=9) + g.add_edge(v[1], v[2], cost=4) + g.add_edge(v[2], v[2]) + g.add_edge(v[2], v[3]) + g.add_edge(v[3], v[4]) + return g, v + + @staticmethod + def description(g): + return [(e.v1.name, e.v2.name, e.cost) for e in g.edges()] + + def test_the_edges_between_the_vertices_are_kept(self): + g, v = self.graph() + h = g.subgraph([v[1], v[2], v[3]]) + self.assertIsInstance(h, self.cls) + self.assertEqual(h.n, 3) + # the edges v0-v1 and v3-v4 cross the border: removed; the parallel + # edges and the loop stay, in the order they were added + self.assertEqual( + self.description(h), + [("v1", "v2", 9), ("v1", "v2", 4), ("v2", "v2", 0.0), ("v2", "v3", 1.0)], + ) + + def test_the_original_is_unchanged(self): + g, v = self.graph() + before = self.description(g) + g.subgraph([v[1], v[2]]) + self.assertEqual(g.n, 5) + self.assertEqual(self.description(g), before) + self.assertEqual(g.validate(), []) + + def test_the_result_is_a_valid_graph(self): + g, v = self.graph() + self.assertEqual(g.subgraph([v[1], v[2], v[3]]).validate(), []) + + def test_the_new_graph_is_independent(self): + g, v = self.graph() + h = g.subgraph([v[1], v[2]]) + self.assertIsNot(h["v1"], g["v1"]) + h["v1"].coord[0] = 99.0 + self.assertEqual(g["v1"].coord[0], 1.0) + h.edges()[0].cost = 123 + self.assertEqual(g.edges()[1].cost, 9) + h.add_vertex([5.0, 5.0], name="new") + self.assertEqual(g.n, 5) + self.assertNotIn("new", g) + + def test_names_are_kept_and_order_is_the_graphs(self): + g, v = self.graph() + h = g.subgraph([v[3], v[1], v[2]]) # asked for in another order + self.assertEqual(names(h), ["v1", "v2", "v3"]) + np.testing.assert_array_equal(h["v2"].coord, g["v2"].coord) + + def test_vertices_by_name(self): + g, v = self.graph() + h = g.subgraph(["v2", v[3]]) + self.assertEqual(names(h), ["v2", "v3"]) + self.assertEqual(h.ne, 2) # the loop and v2-v3 + + def test_a_cost_that_was_given_is_kept(self): + g, v = self.graph() + h = g.subgraph([v[1], v[2]]) + # 9 is not the distance between the vertices, which is 1 + self.assertEqual([e.cost for e in h.edges()][:2], [9, 4]) + + def test_nothing(self): + g, _ = self.graph() + h = g.subgraph([]) + self.assertIsInstance(h, self.cls) + self.assertEqual((h.n, h.ne, h.nc), (0, 0, 0)) + + def test_everything_is_a_copy(self): + g, v = self.graph() + h = g.subgraph(g) + self.assertEqual(names(h), names(g)) + self.assertEqual(self.description(h), self.description(g)) + + def test_components_are_those_of_the_new_graph(self): + g, v = self.graph() + self.assertEqual(g.nc, 1) # make the original's cached components stale + h = g.subgraph([v[0], v[1], v[3], v[4]]) # v0-v1 and v3-v4 + self.assertEqual(h.nc, 2) + self.assertEqual(h["v0"].label, h["v1"].label) + self.assertEqual(h["v3"].label, h["v4"].label) + self.assertNotEqual(h["v1"].label, h["v3"].label) + self.assertEqual(g.nc, 1) + + def test_a_vertex_from_another_graph(self): + g, _ = self.graph() + other, w = line(self.cls) + with self.assertRaises(ValueError): + g.subgraph([w[0]]) + + def test_a_vertex_that_was_removed(self): + g, v = self.graph() + g.remove_vertex(v[4]) + with self.assertRaises(ValueError): + g.subgraph([v[4]]) + + def test_unknown_name_and_wrong_type(self): + g, _ = self.graph() + with self.assertRaises(KeyError): + g.subgraph(["nobody"]) + with self.assertRaises(TypeError): + g.subgraph([3]) + + def test_classes_metric_and_attributes_are_kept(self): + g = self.cls(dim=2, vertex_class=MyVertex, edge_class=MyEdge, metric="L1") + a = g.add_vertex([0.0, 0.0], name="a") + b = g.add_vertex([1.0, 1.0], name="b") + g.add_vertex([5.0, 5.0], name="c") + e = g.add_edge(a, b) + e.data = {"road": ["x"]} + a.note = "start" + + h = g.subgraph([a, b]) + self.assertIsInstance(h["a"], MyVertex) + self.assertIsInstance(h.edges()[0], MyEdge) + self.assertEqual(h["a"].note, "start") + self.assertEqual(h.edges()[0].cost, 2.0) # L1 + self.assertEqual(h.metric(np.array([1.0, 1.0])), 2.0) + self.assertIs(h.vertex_class, MyVertex) + self.assertIs(h.edge_class, MyEdge) + # the data are copies, not shared with the original + self.assertEqual(h.edges()[0].data, {"road": ["x"]}) + h.edges()[0].data["road"].append("y") + self.assertEqual(e.data, {"road": ["x"]}) + + def test_with_within(self): + g, v = self.graph() + h = g.subgraph(g.within([[1, -1], [3, 1]])) + self.assertEqual(names(h), ["v1", "v2", "v3"]) + self.assertEqual(h.ne, 4) + + def test_can_be_pickled_and_modified(self): + g, v = self.graph() + h = g.subgraph([v[1], v[2], v[3]]) + h2 = pickle.loads(pickle.dumps(h)) # nosec B301 + self.assertEqual(self.description(h2), self.description(h)) + h.add_edge(h["v1"], h["v3"]) + self.assertEqual(h.ne, 5) + self.assertEqual(h.validate(), []) + + +class TestUndirectedSubgraph(_Subgraph, unittest.TestCase): + cls = UGraph + + def test_an_edge_is_found_from_either_end(self): + g, v = self.graph() + h = g.subgraph([v[2], v[3]]) + self.assertEqual(names(h["v2"].adjacent()), ["v2", "v2", "v3"]) + self.assertEqual(names(h["v3"].adjacent()), ["v2"]) + + +class TestDirectedSubgraph(_Subgraph, unittest.TestCase): + cls = DGraph + + def test_direction_is_kept(self): + g, v = self.graph() + h = g.subgraph([v[1], v[2], v[3]]) + self.assertEqual(names(h["v1"].successors()), ["v2", "v2"]) + self.assertEqual(h["v1"].predecessors(), []) # v0 -> v1 crossed the border + self.assertEqual(names(h["v3"].predecessors()), ["v2"]) + self.assertEqual(h["v3"].successors(), []) # v3 -> v4 crossed it too + + +if __name__ == "__main__": + unittest.main() From caec1ab321b960e0424cb4907267851bd82992bf Mon Sep 17 00:00:00 2001 From: Peter Corke Date: Tue, 13 Oct 2026 07:07:18 +1000 Subject: [PATCH 2/3] feat: coordinates(), build_kdtree(), subgraph mask, validate coordinates coordinates() is the matrix of vertex coordinates, one column per vertex in insertion order (NaN for a vertex with none), a view of a row-per-vertex array. within() and near() are now computed on it with numpy, for the built-in metrics, instead of vertex by vertex. build_kdtree() builds an optional scipy cKDTree that closest() and near() use for the built-in L1, L2 and SE2 (with its weight, using a periodic dimension for the angle). scipy is imported only inside the method and is not added to pyproject.toml. The index is dropped when a vertex is added or removed or a coord is assigned (Vertex.__setattr__), so a stale index is never used; a coordinate edited in place is documented as not seen. It is not copied or pickled. subgraph() also takes a boolean mask of length n. An integer vector is refused, since it could be indices. validate() now checks that coordinates are finite 1-D numeric vectors and, when dim was not given, all of one length, and with embedded=True that every vertex has one. Co-Authored-By: Claude Sonnet 5.5 --- docs/source/analysis.rst | 21 +- src/pgraph/PGraph.py | 392 ++++++++++++++++++++++++++++++++---- tests/test_spatial_index.py | 348 ++++++++++++++++++++++++++++++++ tests/test_subgraph.py | 58 ++++++ 4 files changed, 777 insertions(+), 42 deletions(-) create mode 100644 tests/test_spatial_index.py diff --git a/docs/source/analysis.rst b/docs/source/analysis.rst index a83d9879..fc1dd09f 100644 --- a/docs/source/analysis.rst +++ b/docs/source/analysis.rst @@ -50,18 +50,27 @@ Geometry and consistency ^^^^^^^^^^^^^^^^^^^^^^^^^ For a graph whose vertices have coordinates, -:meth:`~pgraph.PGraph.Graph.closest` finds the nearest vertex to a point, -:meth:`~pgraph.PGraph.Graph.near` the vertices within a distance of a point and -:meth:`~pgraph.PGraph.Graph.within` the vertices inside a box. Passing the -result of the last to :meth:`~pgraph.PGraph.Graph.subgraph` gives the part of -the graph in a region, without the edges that cross its border. +:meth:`~pgraph.PGraph.Graph.coordinates` is the matrix of them, one column per +vertex. :meth:`~pgraph.PGraph.Graph.closest` finds the nearest vertex to a +point, :meth:`~pgraph.PGraph.Graph.near` the vertices within a distance of a +point and :meth:`~pgraph.PGraph.Graph.within` the vertices inside a box. Passing +the result of the last to :meth:`~pgraph.PGraph.Graph.subgraph` gives the part of +the graph in a region, without the edges that cross its border, and a boolean +mask computed from the coordinates selects in the same way. + +These are a single pass over the vertices. For a large graph that is queried +often, :meth:`~pgraph.PGraph.Graph.build_kdtree` builds an index that +:meth:`~pgraph.PGraph.Graph.closest` and :meth:`~pgraph.PGraph.Graph.near` then +use, which needs scipy and is for the built-in metrics. .. autosummary:: :nosignatures: + ~pgraph.PGraph.Graph.coordinates ~pgraph.PGraph.Graph.closest ~pgraph.PGraph.Graph.near ~pgraph.PGraph.Graph.within + ~pgraph.PGraph.Graph.build_kdtree ~pgraph.PGraph.Graph.subgraph ~pgraph.PGraph.Graph.validate @@ -163,9 +172,11 @@ Details .. automethod:: pgraph.PGraph.Graph.samecomponent .. automethod:: pgraph.PGraph.Graph.connectivity .. automethod:: pgraph.PGraph.Graph.average_degree +.. automethod:: pgraph.PGraph.Graph.coordinates .. automethod:: pgraph.PGraph.Graph.closest .. automethod:: pgraph.PGraph.Graph.near .. automethod:: pgraph.PGraph.Graph.within +.. automethod:: pgraph.PGraph.Graph.build_kdtree .. automethod:: pgraph.PGraph.Graph.subgraph .. automethod:: pgraph.PGraph.Graph.validate diff --git a/src/pgraph/PGraph.py b/src/pgraph/PGraph.py index 44a40ffd..40786fff 100644 --- a/src/pgraph/PGraph.py +++ b/src/pgraph/PGraph.py @@ -12,6 +12,7 @@ from mpl_toolkits.mplot3d import proj3d import copy import heapq +import importlib from collections import deque from collections.abc import Iterable, Iterator import tempfile @@ -32,6 +33,16 @@ def _L2(v: NDArray) -> float: return float(np.linalg.norm(v)) +class _KDTree: + """A spatial index of a graph's vertices and how to query it (private)""" + + def __init__(self, tree: Any, index: NDArray, p: float, prepare: Any): + self.tree = tree # scipy.spatial.cKDTree + self.index = index # position in the graph's vertex list of each point + self.p = p # the Minkowski norm that the tree uses + self.prepare = prepare # point -> point for the tree, or None if out of range + + class SE2Metric: r""" Distance between planar poses, a weighted mixed norm @@ -335,6 +346,10 @@ def __init__( self._dim = dim self._ncomponents = 0 self._connectivitychange = False + # the optional spatial index of build_kdtree(), and whether the + # vertices have changed since it was built + self._kdtree: _KDTree | None = None + self._spatialchange = False if metric is None: self.metric = "L2" else: @@ -653,10 +668,12 @@ def edge_class(self) -> type[Edge]: """ return self._edge_class - def validate(self) -> list[str]: + def validate(self, embedded: bool = False) -> list[str]: """ Check that the graph's records agree with each other + :param embedded: also report a vertex that has no coordinate, defaults + to False :return: a description of each inconsistency found, empty if the graph is consistent :rtype: list of str @@ -682,6 +699,10 @@ def validate(self) -> list[str]: - every vertex points back at this graph, is an instance of :attr:`vertex_class`, is in the name table under its own name, and has a coordinate of the required length if ``dim`` was given + - every coordinate is a finite one-dimensional numeric vector, and if + ``dim`` was not given they all have the same length + - with ``embedded=True``, every vertex has a coordinate. A graph with + no coordinates is valid, but one that is half embedded rarely is - every entry of the name table is a vertex of this graph and is filed under its own name - each edge is in the right vertex lists, the right number of times: @@ -737,6 +758,7 @@ def validate(self) -> list[str]: f"name table entry {key!r} refers to vertex {v.name!r}, " "which has a different name" ) + length: int | None = None # of the coordinates, if dim was not given for v in vertices: if v._graph is not self: problems.append(f"vertex {v.name!r} does not point back at this graph") @@ -758,6 +780,32 @@ def validate(self) -> list[str]: f"vertex {v.name!r} has a coordinate of length " f"{len(v.coord)} but the graph requires dim={self._dim}" ) + if v.coord is None: + if embedded: + problems.append(f"vertex {v.name!r} has no coordinate") + else: + try: + coord = np.asarray(v.coord, dtype=float) + except (TypeError, ValueError): + problems.append(f"vertex {v.name!r} has a coordinate that is not numeric") + else: + if coord.ndim != 1: + problems.append( + f"vertex {v.name!r} has a coordinate of shape " + f"{coord.shape}, not a vector" + ) + elif not np.all(np.isfinite(coord)): + problems.append( + f"vertex {v.name!r} has a coordinate that is not finite" + ) + elif self._dim is None: + if length is None: + length = len(coord) + elif len(coord) != length: + problems.append( + f"vertex {v.name!r} has a coordinate of length " + f"{len(coord)}, but others have {length}" + ) if not directed and v._inlist: problems.append( f"vertex {v.name!r} has {len(v._inlist)} incoming edges " @@ -890,13 +938,15 @@ def copy(self) -> Graph: """ return copy.deepcopy(self) - def subgraph(self, vertices: Iterable[Vertex | str]) -> Graph: + def subgraph(self, vertices: Iterable[Vertex | str] | NDArray) -> Graph: """ Subgraph induced by a set of vertices - :param vertices: vertices of this graph, or their names - :type vertices: iterable of Vertex subclass or str - :raises ValueError: a vertex is not in this graph + :param vertices: vertices of this graph, or their names, or a boolean + mask with one element per vertex + :type vertices: iterable of Vertex subclass or str, or ndarray(n) of bool + :raises ValueError: a vertex is not in this graph, or a mask is not + the length of the number of vertices :raises KeyError: a name is not that of a vertex in this graph :raises TypeError: an item is neither a vertex nor a string :return: a new graph of the same kind @@ -929,18 +979,47 @@ def subgraph(self, vertices: Iterable[Vertex | str]) -> Graph: the vertices that satisfy a test, or those in a region of an embedded graph from :meth:`within`. + A boolean mask, such as the result of a comparison on + :meth:`coordinates`, selects the vertices by their position in the + graph, the ``True`` elements being kept:: + + h = g.subgraph(g.coordinates()[0] > 0) # the vertices with x > 0 + + The mask must be of boolean type (a list of ``bool`` also works) and of + length ``g.n``. A vector of 0 and 1 of integer type is not accepted, since + it could as well be a list of indices. + :seealso: :meth:`copy` :meth:`within` :meth:`component` """ - wanted = set() - for v in vertices: - v = self._resolve_vertex(v, "vertex") - if v._graph is not self: - raise ValueError(f"vertex {v.name!r} is not in this graph") - wanted.add(id(v)) + mask = None + if isinstance(vertices, np.ndarray) and vertices.dtype == np.bool_: + mask = vertices + elif ( + isinstance(vertices, (list, tuple)) + and len(vertices) > 0 + and all(isinstance(x, (bool, np.bool_)) for x in vertices) + ): + mask = np.array(vertices, dtype=bool) + + if mask is not None: + if mask.shape != (len(self._vertexlist),): + raise ValueError( + f"a mask must have one element per vertex, {len(self._vertexlist)}, " + f"not shape {mask.shape}" + ) + keep = [int(i) for i in np.flatnonzero(mask)] + else: + wanted = set() + items: Iterable[Any] = vertices + for item in items: + v = self._resolve_vertex(item, "vertex") + if v._graph is not self: + raise ValueError(f"vertex {v.name!r} is not in this graph") + wanted.add(id(v)) + keep = [i for i, v in enumerate(self._vertexlist) if id(v) in wanted] state = self.__getstate__() flat = state["_flat"] - keep = [i for i, v in enumerate(self._vertexlist) if id(v) in wanted] renumber = {old: new for new, old in enumerate(keep)} flat["vertices"] = [flat["vertices"][i] for i in keep] flat["edges"] = [ @@ -957,7 +1036,7 @@ def subgraph(self, vertices: Iterable[Vertex | str]) -> Graph: # pickling and copying #: the graph's own containers, which are saved as a flat description - _STRUCTURE = ("_vertexlist", "_vertexdict", "_edges") + _STRUCTURE = ("_vertexlist", "_vertexdict", "_edges", "_kdtree") #: the references a vertex holds, which are rebuilt rather than saved _VERTEX_LINKS = ("_edgelist", "_inlist", "_graph") @@ -1019,6 +1098,8 @@ def __setstate__(self, state: dict[str, Any]) -> None: if isinstance(value, tuple) and value[0] == "builtin": state[key] = _METRICS[value[1]] self.__dict__.update(state) + self._kdtree = None + self._spatialchange = False self._vertexlist = [] self._vertexdict = {} @@ -1158,6 +1239,7 @@ def add_vertex( print(f"New vertex {vertex.name}: {vertex.coord}") vertex._graph = self self._connectivitychange = True + self._spatialchange = True return vertex @classmethod @@ -1442,6 +1524,7 @@ def remove_vertex(self, vertex: Vertex) -> None: vertex._graph = None # it is not in any graph now vertex._label = None # and so not in a component self._connectivitychange = True + self._spatialchange = True def remove(self, x: Edge | Vertex) -> None: """ @@ -1870,6 +1953,197 @@ def __contains__(self, item: Vertex | str) -> bool: elif isinstance(item, Vertex): return item in self._vertexdict.values() + def coordinates(self) -> NDArray: + """ + Coordinates of all the vertices as a matrix + + :return: coordinates, one column per vertex in the order they were added + :rtype: ndarray(dim, n) + :raises ValueError: the coordinates of two vertices differ in length + + Column ``i`` is the coordinate of ``g[i]``, the same order as the + rows and columns of :meth:`adjacency`. A vertex that has no + coordinate has a column of NaN, so the matrix always has one column + per vertex. The number of rows is ``dim`` if the graph was given one, + otherwise the length of the first coordinate, and 0 if no vertex has + one. Elements are float. + + The result is a *view* of a new array that has one row per vertex, so + the coordinate of each vertex is contiguous in memory and ``C[:, i]`` + is a contiguous slice, while a row such as ``C[0]`` is strided. Use + :func:`numpy.ascontiguousarray` if you need the rows to be contiguous. + + .. note:: This is a snapshot, made each time: changing it does not + move the vertices, and moving a vertex does not change a matrix + that was already returned. + + With NumPy this selects vertices in a single operation, as + :meth:`within` and :meth:`near` do, or it can be given to a spatial + index of your own. + + .. runblock:: pycon + + >>> from pgraph import UGraph + >>> g = UGraph() + >>> for i in range(3): + ... g.add_vertex(coord=[i, 10 * i]) + >>> g.coordinates() + + :seealso: :meth:`within` :meth:`near` :meth:`closest` :meth:`adjacency` + """ + n = len(self._vertexlist) + dim = self._dim + if dim is None: + dim = next((len(v.coord) for v in self._vertexlist if v.coord is not None), 0) + out = np.full((n, dim), np.nan) + for i, vertex in enumerate(self._vertexlist): + if vertex.coord is None: + continue + if len(vertex.coord) != dim: + raise ValueError( + f"vertex {vertex.name!r} has a coordinate of length " + f"{len(vertex.coord)}, but the others have {dim}" + ) + out[i] = vertex.coord + return out.T + + def build_kdtree(self) -> None: + """ + Build a spatial index of the vertices, used by :meth:`closest` and :meth:`near` + + :raises ValueError: the metric is not one of the built-in ``"L1"``, + ``"L2"`` or ``"SE2"``, or ``"SE2"`` and the coordinates are not + 3 elements, or no vertex has a coordinate + :raises ImportError: scipy is not installed + + Finding the nearest vertex to a point is a pass over all the vertices, + which is the cost of a query on a large graph that is queried often. A + k-d tree answers in a time that grows more slowly, at the cost of + building it once. :meth:`closest` and :meth:`near` use the index if it + has been built, and fall back to the pass over the vertices if it has + not or if the graph has changed. + + The index is for the built-in metrics ``"L1"``, ``"L2"`` and + ``"SE2"`` (with any weight, see :class:`SE2Metric`), and results + are the same as without it, except that when two vertices are exactly + the same distance the one returned may differ. + + .. warning:: The index is dropped when a vertex is added or removed + or when ``vertex.coord`` is assigned (including ``+=``), and + :meth:`closest` is then correct but slow, until it is built again. + It cannot see a coordinate that is changed in place, + ``vertex.coord[0] = x``. Do not change vertices in place after + building, or build again. + + The index is not copied or pickled with the graph. + + .. note:: This needs `scipy `_. pgraph does not + list it as a dependency and imports it only when this method is + called, but `spatialmath-python`, which pgraph does depend on, + installs it. + + :seealso: :meth:`closest` :meth:`near` :meth:`coordinates` + """ + metric = self._metric + weight: float | None = None + if metric is _METRICS["L1"]: + p = 1 + elif metric is _METRICS["L2"]: + p = 2 + elif isinstance(metric, SE2Metric): + p, weight = 2, metric.weight + else: + raise ValueError( + "build_kdtree is for the built-in metrics 'L1', 'L2' and 'SE2', " + f"not {metric!r}" + ) + try: + spatial = importlib.import_module("scipy.spatial") + except ImportError as err: + raise ImportError( + "build_kdtree needs scipy, which pgraph does not itself require: " + "pip install scipy" + ) from err + + points = self.coordinates().T # one row per vertex + if points.shape[1] == 0: + raise ValueError("no vertex has a coordinate") + index = np.flatnonzero(np.isfinite(points).all(axis=1)) + if len(index) == 0: + raise ValueError("no vertex has a coordinate") + points = points[index] + + prepare: Any = None + if weight is None: + tree = spatial.cKDTree(points) + else: + if points.shape[1] != 3: + raise ValueError( + f"the SE2 metric needs 3-element coordinates, not {points.shape[1]}" + ) + # scipy has periodic dimensions for the angle, but a periodic box + # must be finite and the data in it, so put the plane inside a box + # much larger than the data and wrap the angle, scaled by its weight + span = max(float(np.ptp(points[:, :2], axis=0).max()), 1.0) + size = 1e3 * span + shift = points[:, :2].min(axis=0) - size / 4 + period = 2 * np.pi * weight + + def prepare(q: NDArray) -> NDArray | None: # noqa: F811 + xy = q[:2] - shift + if not np.all((xy >= 0) & (xy < size)): + return None # far outside the data: let the caller do it + if weight == 0: + return xy + return np.r_[xy, np.mod(weight * q[2], period)] + + if weight == 0: + tree = spatial.cKDTree(points[:, :2] - shift, boxsize=[size, size]) + else: + data = np.column_stack( + [points[:, :2] - shift, np.mod(weight * points[:, 2], period)] + ) + tree = spatial.cKDTree(data, boxsize=[size, size, period]) + + self._kdtree = _KDTree(tree, index, p, prepare) + self._spatialchange = False + + def _index(self) -> _KDTree | None: + """The spatial index, if it is built and still valid (private method)""" + if self._kdtree is not None and self._spatialchange: + self._kdtree = None # vertices were added or removed: no longer right + return self._kdtree + + def _index_point(self, index: _KDTree, coord: ArrayLike) -> NDArray | None: + """The point to query the index with, or None if it cannot be (private)""" + q = np.asarray(coord, dtype=float) + if q.ndim != 1 or not np.all(np.isfinite(q)): + return None + if index.prepare is None: + return q if q.shape[0] == index.tree.m else None + if q.shape[0] != 3: + return None + return index.prepare(q) + + def _metric_columns(self, diff: NDArray) -> NDArray | None: + """ + The metric of each column of a matrix, or None (private method) + + Only for the built-in metrics, which can be applied to all columns in + one operation; the metric of a user is a function of one vector. + """ + metric = self._metric + if metric is _METRICS["L1"]: + return np.asarray(np.abs(diff).sum(axis=0)) + if metric is _METRICS["L2"]: + return np.asarray(np.sqrt((diff**2).sum(axis=0))) + if isinstance(metric, SE2Metric) and diff.shape[0] == 3: + theta = (diff[2] + np.pi) % (2 * np.pi) - np.pi + return np.asarray( + np.sqrt(diff[0] ** 2 + diff[1] ** 2 + (metric.weight * theta) ** 2) + ) + return None + def closest(self, coord: ArrayLike) -> tuple[Vertex | None, float]: """ Vertex closest to point @@ -1884,6 +2158,9 @@ def closest(self, coord: ArrayLike) -> tuple[Vertex | None, float]: according to the graph's metric. Vertices without a coordinate (``coord`` is None) are skipped -- they have no position to compare. + This is a pass over the vertices, unless a spatial index has been + built with :meth:`build_kdtree` and is still valid. + .. runblock:: pycon >>> from pgraph import UGraph @@ -1893,8 +2170,16 @@ def closest(self, coord: ArrayLike) -> tuple[Vertex | None, float]: >>> vertex, d = g.closest([1, 1]) >>> print(vertex, d) - :seealso: :meth:`metric` + :seealso: :meth:`metric` :meth:`near` :meth:`build_kdtree` """ + index = self._index() + if index is not None: + q = self._index_point(index, coord) + if q is not None: + d, i = index.tree.query(q, p=index.p) + if np.isfinite(d): + return self._vertexlist[int(index.index[i])], float(d) + min_dist = np.inf min_which: Vertex | None = None @@ -1915,8 +2200,8 @@ def within(self, box: ArrayLike) -> list[Vertex]: :param box: lower and upper corners of the box :type box: ndarray(2,n) :raises ValueError: ``box`` is not 2 x n, its upper corner is below its - lower corner in some dimension, or the coordinate of a vertex does - not have n elements + lower corner in some dimension, or the coordinates of the vertices + do not have n elements :return: the vertices inside the box, in the order of the graph :rtype: list of Vertex subclass @@ -1925,8 +2210,9 @@ def within(self, box: ArrayLike) -> list[Vertex]: The box is closed: a vertex on its boundary is inside. Vertices without a coordinate are skipped, as they are by :meth:`closest`. - It is a single pass over the vertices. Together with - :meth:`subgraph` it selects the part of an embedded graph in a region + The test is made for all the vertices at once on :meth:`coordinates`. + Together with :meth:`subgraph` it selects the part of an embedded + graph in a region .. runblock:: pycon @@ -1937,7 +2223,7 @@ def within(self, box: ArrayLike) -> list[Vertex]: >>> g.within([[1, 1], [3, 3]]) >>> g.subgraph(g.within([[1, 1], [3, 3]])) - :seealso: :meth:`near` :meth:`closest` :meth:`subgraph` + :seealso: :meth:`near` :meth:`closest` :meth:`subgraph` :meth:`coordinates` """ box = np.asarray(box, dtype=float) if box.ndim != 2 or box.shape[0] != 2: @@ -1946,19 +2232,16 @@ def within(self, box: ArrayLike) -> list[Vertex]: if np.any(upper < lower): raise ValueError("the upper corner of the box is below the lower corner") - found = [] - for vertex in self: - if vertex.coord is None: - continue - coord = np.asarray(vertex.coord, dtype=float) - if coord.shape != lower.shape: - raise ValueError( - f"vertex {vertex.name!r} has a coordinate of {coord.size} " - f"elements, the box has {lower.size}" - ) - if np.all(coord >= lower) and np.all(coord <= upper): - found.append(vertex) - return found + C = self.coordinates() + if C.shape[0] == 0 or C.shape[1] == 0: + return [] # no vertex has a coordinate + if C.shape[0] != lower.size: + raise ValueError( + f"the vertices have coordinates of {C.shape[0]} elements, " + f"the box has {lower.size}" + ) + inside = np.all((C >= lower[:, None]) & (C <= upper[:, None]), axis=0) + return [self._vertexlist[i] for i in np.flatnonzero(inside)] def near(self, coord: ArrayLike, radius: float) -> list[tuple[Vertex, float]]: """ @@ -1977,6 +2260,11 @@ def near(self, coord: ArrayLike, radius: float) -> list[tuple[Vertex, float]]: within a radius. Vertices without a coordinate are skipped, and vertices at equal distance are in the order of the graph. + For the built-in metrics the distance to every vertex is found in one + operation on :meth:`coordinates`, or by the spatial index if + :meth:`build_kdtree` was called; a metric of your own is applied + to each vertex in turn. + .. runblock:: pycon >>> from pgraph import UGraph @@ -1985,18 +2273,38 @@ def near(self, coord: ArrayLike, radius: float) -> list[tuple[Vertex, float]]: ... g.add_vertex(coord=[i, 0], name=f"v{i}") >>> [(v.name, d) for v, d in g.near([2.2, 0], 1.5)] - :seealso: :meth:`closest` :meth:`within` :meth:`metric` + :seealso: :meth:`closest` :meth:`within` :meth:`metric` :meth:`build_kdtree` """ if radius < 0: raise ValueError(f"radius must not be negative, got {radius!r}") + candidates: list[Vertex] | None = None + index = self._index() + if index is not None: + q = self._index_point(index, coord) + if q is not None: + # the index finds the candidates, the metric gives the distance + hits = index.tree.query_ball_point(q, radius * (1 + 1e-9) + 1e-12, p=index.p) + candidates = [self._vertexlist[int(index.index[i])] for i in sorted(hits)] + + if candidates is None: + C = self.coordinates() + point = np.asarray(coord, dtype=float) + if C.shape[0] != 0 and point.shape == (C.shape[0],): + d = self._metric_columns(C - point[:, None]) + if d is not None: + hit = np.flatnonzero(d <= radius) + hit = hit[np.argsort(d[hit], kind="stable")] + return [(self._vertexlist[i], float(d[i])) for i in hit] + candidates = [v for v in self._vertexlist if v.coord is not None] + found = [] - for vertex in self: + for vertex in candidates: if vertex.coord is None: continue - d = float(self.metric(vertex.coord - coord)) - if d <= radius: - found.append((vertex, d)) + dist = float(self.metric(vertex.coord - coord)) + if dist <= radius: + found.append((vertex, dist)) found.sort(key=lambda pair: pair[1]) return found @@ -4194,6 +4502,16 @@ def __init__(self, coord: ArrayLike | None = None, name: str | None = None): self._graph: Graph | None = None # reference to owning graph # print('Vertex init', type(self)) + def __setattr__(self, name: str, value: Any) -> None: + # a vertex that moves makes its graph's spatial index (build_kdtree) + # stale; this sees ``v.coord = x`` and ``v.coord += x`` but not + # ``v.coord[0] = x``, which is why that is documented + if name == "coord": + graph = self.__dict__.get("_graph") + if graph is not None: + graph._spatialchange = True + object.__setattr__(self, name, value) + def __reduce__(self) -> Any: """ Refuse to pickle or copy a vertex on its own diff --git a/tests/test_spatial_index.py b/tests/test_spatial_index.py new file mode 100644 index 00000000..36dc5b0a --- /dev/null +++ b/tests/test_spatial_index.py @@ -0,0 +1,348 @@ +""" +``Graph.coordinates``, the optional spatial index of ``Graph.build_kdtree`` +that ``closest`` and ``near`` use, and what ``validate`` says about +coordinates. + +The index needs scipy, which pgraph does not depend on, so its tests are +skipped where it is absent; the test that scipy is only imported on demand +is not. Expectations are computed here by brute force, not by the code +under test. +""" + +import ast +import copy +import importlib +import importlib.util +import pathlib +import pickle # nosec B403 - only used to round-trip graphs made in this file +import unittest +import unittest.mock + +import numpy as np + +import pgraph +from pgraph import DGraph, SE2Metric, UGraph + +HAVE_SCIPY = importlib.util.find_spec("scipy") is not None + + +def points_graph(cls=UGraph, n=200, dim=2, seed=0, **kwargs): + rng = np.random.default_rng(seed) + g = cls(dim=dim, **kwargs) + P = rng.uniform(-10, 10, (n, dim)) + if dim == 3: + P[:, 2] = rng.uniform(-np.pi, np.pi, n) + for i, p in enumerate(P): + g.add_vertex(p, name=f"v{i}") + return g, P + + +def wrap(a): + return (a + np.pi) % (2 * np.pi) - np.pi + + +def se2(D, weight): + return np.sqrt(D[:, 0] ** 2 + D[:, 1] ** 2 + (weight * wrap(D[:, 2])) ** 2) + + +class TestCoordinates(unittest.TestCase): + def test_one_column_per_vertex_in_order(self): + g = UGraph() + g.add_vertex([1, 2], name="a") + g.add_vertex([3, 4], name="b") + g.add_vertex([5, 6], name="c") + C = g.coordinates() + self.assertEqual(C.shape, (2, 3)) + np.testing.assert_array_equal(C, [[1, 3, 5], [2, 4, 6]]) + np.testing.assert_array_equal(C[:, 1], g["b"].coord) + self.assertEqual(C.dtype, np.float64) + + def test_is_a_view_with_each_coordinate_contiguous(self): + g, _ = points_graph(n=10) + C = g.coordinates() + self.assertTrue(C.T.flags.c_contiguous) # a row per vertex underneath + self.assertIsNotNone(C.base) + self.assertTrue(C[:, 3].flags.c_contiguous) + + def test_vertex_without_coordinate_is_nan(self): + g = UGraph() + g.add_vertex([1.0, 2.0], name="a") + g.add_vertex(name="nowhere") + g.add_vertex([5.0, 6.0], name="c") + C = g.coordinates() + self.assertEqual(C.shape, (2, 3)) + self.assertTrue(np.isnan(C[:, 1]).all()) + np.testing.assert_array_equal(C[:, [0, 2]], [[1, 5], [2, 6]]) + + def test_no_coordinates_at_all(self): + g = UGraph() + g.add_vertex() + g.add_vertex() + self.assertEqual(g.coordinates().shape, (0, 2)) + self.assertEqual(UGraph().coordinates().shape, (0, 0)) + + def test_dim_of_the_graph_gives_the_rows(self): + g = UGraph(dim=3) + g.add_vertex() + self.assertEqual(g.coordinates().shape, (3, 1)) + + def test_different_lengths(self): + g = UGraph() + g.add_vertex([1.0, 2.0], name="a") + g.add_vertex([1.0, 2.0, 3.0], name="b") + with self.assertRaises(ValueError) as cm: + g.coordinates() + self.assertIn("'b'", str(cm.exception)) + + def test_a_snapshot(self): + g, _ = points_graph(n=5) + C = g.coordinates() + before = C.copy() + g["v0"].coord = np.r_[99.0, 99.0] + np.testing.assert_array_equal(C, before) # not changed by moving a vertex + C[:] = 0 + self.assertNotEqual(g["v1"].coord[0], 0.0) # nor the vertices by it + + def test_columns_follow_removal(self): + g = UGraph() + for i in range(3): + g.add_vertex([float(i), 0.0], name=f"v{i}") + g.remove_vertex(g["v1"]) + np.testing.assert_array_equal(g.coordinates(), [[0, 2], [0, 0]]) + + +@unittest.skipUnless(HAVE_SCIPY, "scipy is not installed") +class TestIndexAgreesWithTheLoop(unittest.TestCase): + def queries(self, dim, seed=1, n=60): + rng = np.random.default_rng(seed) + Q = rng.uniform(-12, 12, (n, dim)) + if dim == 3: + Q[:, 2] = rng.uniform(-4, 4, n) + return Q + + def check(self, cls, metric, expected_distances, dim): + g, P = points_graph(cls, dim=dim, metric=metric) + g.build_kdtree() + for q in self.queries(dim): + d = expected_distances(P - q) + vertex, dist = g.closest(q) + self.assertAlmostEqual(dist, d.min(), places=9) + self.assertAlmostEqual(d[int(vertex.name[1:])], d.min(), places=9) + for radius in (0.0, 1.0, 4.0): + found = g.near(q, radius) + self.assertEqual( + sorted(int(v.name[1:]) for v, _ in found), + sorted(np.flatnonzero(d <= radius).tolist()), + ) + np.testing.assert_allclose([x for _, x in found], np.sort(d[d <= radius]), atol=1e-9) + + def test_L2(self): + for cls in (UGraph, DGraph): + self.check(cls, "L2", lambda D: np.linalg.norm(D, axis=1), 2) + + def test_L1(self): + self.check(UGraph, "L1", lambda D: np.abs(D).sum(axis=1), 2) + + def test_L2_in_three_dimensions(self): + self.check(UGraph, "L2", lambda D: np.linalg.norm(D, axis=1), 3) + + def test_SE2_with_weights(self): + for weight in (0.0, 0.3, 1.0, 4.0): + with self.subTest(weight=weight): + self.check(UGraph, SE2Metric(weight), lambda D, w=weight: se2(D, w), 3) + + def test_the_named_SE2(self): + self.check(UGraph, "SE2", lambda D: se2(D, 1.0), 3) + + def test_a_query_across_the_angle_wrap(self): + g = UGraph(dim=3, metric=SE2Metric(2.0)) + g.add_vertex([0.0, 0.0, np.pi - 0.05], name="a") + g.add_vertex([0.0, 0.1, 0.0], name="b") + g.build_kdtree() + v, d = g.closest([0.0, 0.0, -np.pi + 0.05]) # 0.1 rad from a, weight 2 + self.assertEqual(v.name, "a") + self.assertAlmostEqual(d, 0.2) + + def test_a_query_far_outside_is_still_right(self): + g, P = points_graph(UGraph, dim=3, metric=SE2Metric(1.0)) + g.build_kdtree() + q = np.r_[1e9, -1e9, 0.0] + v, d = g.closest(q) + self.assertAlmostEqual(d, se2(P - q, 1.0).min(), delta=1e-3) + + def test_vertices_without_coordinates_are_skipped(self): + g = UGraph(dim=2) + g.add_vertex(name="nowhere") + g.add_vertex([5.0, 5.0], name="a") + g.add_vertex(name="also nowhere") + g.add_vertex([0.0, 0.0], name="b") + g.build_kdtree() + self.assertEqual(g.closest([4.0, 4.0])[0].name, "a") + self.assertEqual([v.name for v, _ in g.near([0, 0], 1)], ["b"]) + + +@unittest.skipUnless(HAVE_SCIPY, "scipy is not installed") +class TestIndexLifetime(unittest.TestCase): + def graph(self): + g = UGraph(dim=2) + g.add_vertex([0.0, 0.0], name="a") + g.add_vertex([10.0, 0.0], name="b") + g.build_kdtree() + return g + + def test_a_vertex_added_is_found(self): + g = self.graph() + self.assertEqual(g.closest([9.0, 0.0])[0].name, "b") + g.add_vertex([9.0, 0.0], name="new") + self.assertEqual(g.closest([9.0, 0.0])[0].name, "new") + self.assertEqual({v.name for v, _ in g.near([9.0, 0.0], 1.5)}, {"b", "new"}) + + def test_a_vertex_removed_is_not(self): + g = self.graph() + g.remove_vertex(g["b"]) + self.assertEqual(g.closest([9.0, 0.0])[0].name, "a") + + def test_assigning_a_coordinate(self): + g = self.graph() + g["a"].coord = np.r_[9.0, 0.0] + self.assertEqual(g.closest([9.0, 0.0]), (g["a"], 0.0)) + + def test_augmented_assignment_of_a_coordinate(self): + g = self.graph() + g["a"].coord += np.r_[9.0, 0.0] + self.assertEqual(g.closest([9.0, 0.0]), (g["a"], 0.0)) + + def test_adding_an_edge_keeps_it(self): + g = self.graph() + g.add_edge(g["a"], g["b"]) + self.assertIsNotNone(g._kdtree) # white box: edges do not move vertices + self.assertEqual(g.closest([9.0, 0.0])[0].name, "b") + + def test_a_coordinate_changed_in_place_is_not_seen(self): + # the documented limitation: the index is kept and gives the old answer + g = self.graph() + g["a"].coord[0] = 9.0 + self.assertEqual(g.closest([9.0, 0.0])[0].name, "b") # a is really nearer + g.build_kdtree() # as the documentation says to do + self.assertEqual(g.closest([9.0, 0.0])[0].name, "a") + + def test_not_in_copies_or_pickles(self): + g = self.graph() + for h in (g.copy(), copy.deepcopy(g), pickle.loads(pickle.dumps(g)), g.subgraph(g)): # nosec B301 + self.assertIsNone(h._kdtree) + self.assertEqual(h.closest([9.0, 0.0])[0].name, "b") + + def test_a_metric_of_your_own(self): + g = UGraph(dim=2, metric=lambda v: float(np.abs(v).max())) + g.add_vertex([0.0, 0.0]) + with self.assertRaises(ValueError): + g.build_kdtree() + self.assertEqual(g.closest([1.0, 0.5])[1], 1.0) # still works, by the loop + + def test_se2_needs_three_elements(self): + g = UGraph(dim=2, metric="SE2") + g.add_vertex([0.0, 0.0]) + with self.assertRaises(ValueError): + g.build_kdtree() + + def test_no_coordinates(self): + g = UGraph() + g.add_vertex() + with self.assertRaises(ValueError): + g.build_kdtree() + + +class TestScipyIsOnlyImportedOnDemand(unittest.TestCase): + def test_no_module_of_pgraph_imports_scipy_at_import_time(self): + # spatialmath, which pgraph depends on, imports scipy itself, so + # whether scipy is loaded cannot be tested; what pgraph does can + package = pathlib.Path(pgraph.__file__).parent + for path in package.glob("*.py"): + tree = ast.parse(path.read_text()) + for node in tree.body: # module level only + names = [] + if isinstance(node, ast.Import): + names = [a.name for a in node.names] + elif isinstance(node, ast.ImportFrom): + names = [node.module or ""] + for name in names: + self.assertFalse( + name.split(".")[0] == "scipy", f"{path.name} imports {name}" + ) + + def test_without_scipy_there_is_a_clear_error(self): + g = UGraph(dim=2) + g.add_vertex([0.0, 0.0]) + real = importlib.import_module + + def refuse(name, *args, **kwargs): + if name.startswith("scipy"): + raise ImportError("No module named 'scipy'") + return real(name, *args, **kwargs) + + with unittest.mock.patch("importlib.import_module", refuse): + with self.assertRaises(ImportError) as cm: + g.build_kdtree() + self.assertIn("scipy", str(cm.exception)) + self.assertIn("pip install scipy", str(cm.exception)) + self.assertEqual(g.closest([1.0, 1.0])[1], np.sqrt(2)) # the loop works + + +class TestValidateCoordinates(unittest.TestCase): + def test_a_good_graph(self): + g, _ = points_graph(n=5) + self.assertEqual(g.validate(), []) + self.assertEqual(g.validate(embedded=True), []) + + def test_a_graph_with_no_coordinates_is_valid(self): + g = UGraph() + g.add_vertex() + self.assertEqual(g.validate(), []) + + def test_presence_only_when_asked_for(self): + g = UGraph() + g.add_vertex([1.0, 2.0], name="a") + g.add_vertex(name="b") + self.assertEqual(g.validate(), []) + problems = g.validate(embedded=True) + self.assertEqual(len(problems), 1) + self.assertIn("'b'", problems[0]) + + def test_lengths_must_agree_when_dim_was_not_given(self): + g = UGraph() + g.add_vertex([1.0, 2.0], name="a") + g.add_vertex([1.0, 2.0, 3.0], name="b") + problems = g.validate() + self.assertEqual(len(problems), 1) + self.assertIn("'b'", problems[0]) + + def test_not_finite(self): + g = UGraph() + g.add_vertex([1.0, 2.0], name="a") + g["a"].coord = np.r_[np.nan, 2.0] + g.add_vertex([np.inf, 0.0], name="b") + problems = g.validate() + self.assertEqual(len(problems), 2) + self.assertTrue(all("not finite" in p for p in problems)) + + def test_not_a_vector(self): + g = UGraph() + g.add_vertex([1.0, 2.0], name="a") + g["a"].coord = np.zeros((2, 2)) + self.assertIn("not a vector", g.validate()[0]) + + def test_not_numeric(self): + g = UGraph() + g.add_vertex([1.0, 2.0], name="a") + g["a"].coord = np.array(["x", "y"]) + self.assertIn("not numeric", g.validate()[0]) + + def test_length_of_dim_still_checked(self): + g = UGraph(dim=2) + g.add_vertex([1.0, 2.0], name="a") + g["a"].coord = np.r_[1.0, 2.0, 3.0] + self.assertEqual(len(g.validate()), 1) + + +if __name__ == "__main__": + unittest.main() diff --git a/tests/test_subgraph.py b/tests/test_subgraph.py index f9d5e931..765b33d2 100644 --- a/tests/test_subgraph.py +++ b/tests/test_subgraph.py @@ -313,5 +313,63 @@ def test_direction_is_kept(self): self.assertEqual(h["v3"].successors(), []) # v3 -> v4 crossed it too +class TestSubgraphMask(unittest.TestCase): + def graph(self, cls=UGraph): + g, v = line(cls) + for i in range(4): + g.add_edge(v[i], v[i + 1]) + return g + + def test_numpy_boolean_mask(self): + for cls in (UGraph, DGraph): + with self.subTest(cls=cls.__name__): + g = self.graph(cls) + h = g.subgraph(g.coordinates()[0] >= 2) # v2, v3, v4 + self.assertEqual(names(h), ["v2", "v3", "v4"]) + self.assertEqual(h.ne, 2) + self.assertIsInstance(h, cls) + + def test_list_of_booleans(self): + g = self.graph() + h = g.subgraph([True, True, False, False, False]) + self.assertEqual(names(h), ["v0", "v1"]) + self.assertEqual(h.ne, 1) + + def test_the_same_as_vertices(self): + g = self.graph() + mask = np.array([False, True, True, True, False]) + by_mask = g.subgraph(mask) + by_vertices = g.subgraph([g[i] for i in np.flatnonzero(mask)]) + self.assertEqual(names(by_mask), names(by_vertices)) + self.assertEqual( + [(e.v1.name, e.v2.name) for e in by_mask.edges()], + [(e.v1.name, e.v2.name) for e in by_vertices.edges()], + ) + + def test_all_false(self): + g = self.graph() + h = g.subgraph(np.zeros(5, dtype=bool)) + self.assertEqual((h.n, h.ne), (0, 0)) + + def test_wrong_length(self): + g = self.graph() + for mask in (np.ones(4, dtype=bool), np.ones(6, dtype=bool), [True, False]): + with self.subTest(length=len(mask)): + with self.assertRaises(ValueError): + g.subgraph(mask) + + def test_integers_are_not_a_mask(self): + # [0, 1, 1, 0, 0] could as well be vertex indices, so it is refused + g = self.graph() + with self.assertRaises(TypeError): + g.subgraph(np.array([0, 1, 1, 0, 0])) + with self.assertRaises(TypeError): + g.subgraph([0, 1, 1, 0, 0]) + + def test_an_empty_list_is_no_vertices_not_an_empty_mask(self): + g = self.graph() + self.assertEqual(g.subgraph([]).n, 0) + + if __name__ == "__main__": unittest.main() From b761767f7856c82eeee87a6e076c778b2bf596ba Mon Sep 17 00:00:00 2001 From: Peter Corke Date: Tue, 13 Oct 2026 09:10:06 +1000 Subject: [PATCH 3/3] tests: cover near without an index, the guards, and the pickle markers near() with SE2Metric and with a metric of your own, within() with no coordinates, closest() skipping a vertex without one, the index refusing a graph whose coordinates are all NaN, a query that is not finite or the wrong length. The pickle round trips are routed through one helper marked for the scanners. An unreachable guard in near() is an assertion. Co-Authored-By: Claude Sonnet 5.5 --- src/pgraph/PGraph.py | 3 +- tests/test_spatial_index.py | 70 +++++++++++++++++++++++++++++++++++-- tests/test_subgraph.py | 9 +++-- 3 files changed, 76 insertions(+), 6 deletions(-) diff --git a/src/pgraph/PGraph.py b/src/pgraph/PGraph.py index 40786fff..232eb68c 100644 --- a/src/pgraph/PGraph.py +++ b/src/pgraph/PGraph.py @@ -2300,8 +2300,7 @@ def near(self, coord: ArrayLike, radius: float) -> list[tuple[Vertex, float]]: found = [] for vertex in candidates: - if vertex.coord is None: - continue + assert vertex.coord is not None # the candidates all have a position dist = float(self.metric(vertex.coord - coord)) if dist <= radius: found.append((vertex, dist)) diff --git a/tests/test_spatial_index.py b/tests/test_spatial_index.py index 36dc5b0a..b4d98167 100644 --- a/tests/test_spatial_index.py +++ b/tests/test_spatial_index.py @@ -14,7 +14,7 @@ import importlib import importlib.util import pathlib -import pickle # nosec B403 - only used to round-trip graphs made in this file +import pickle # nosec B403 nosemgrep - only used to round-trip graphs made in this file import unittest import unittest.mock @@ -26,6 +26,11 @@ HAVE_SCIPY = importlib.util.find_spec("scipy") is not None +def roundtrip(obj): + """Pickle and unpickle: these tests only read what they have just written""" + return pickle.loads(pickle.dumps(obj)) # nosec B301 nosemgrep + + def points_graph(cls=UGraph, n=200, dim=2, seed=0, **kwargs): rng = np.random.default_rng(seed) g = cls(dim=dim, **kwargs) @@ -228,7 +233,7 @@ def test_a_coordinate_changed_in_place_is_not_seen(self): def test_not_in_copies_or_pickles(self): g = self.graph() - for h in (g.copy(), copy.deepcopy(g), pickle.loads(pickle.dumps(g)), g.subgraph(g)): # nosec B301 + for h in (g.copy(), copy.deepcopy(g), roundtrip(g), g.subgraph(g)): self.assertIsNone(h._kdtree) self.assertEqual(h.closest([9.0, 0.0])[0].name, "b") @@ -252,6 +257,67 @@ def test_no_coordinates(self): g.build_kdtree() +class TestVectorizedAndLoopPaths(unittest.TestCase): + """near() without an index: built-in metrics in one operation, others in a loop""" + + def test_SE2_without_an_index(self): + for weight in (0.0, 0.5, 2.0): + with self.subTest(weight=weight): + g, P = points_graph(UGraph, dim=3, metric=SE2Metric(weight)) + for q in ([0.0, 0.0, 3.0], [4.0, -2.0, -3.0]): + d = se2(P - np.asarray(q), weight) + found = g.near(q, 5.0) + self.assertEqual( + sorted(int(v.name[1:]) for v, _ in found), + sorted(np.flatnonzero(d <= 5.0).tolist()), + ) + np.testing.assert_allclose([x for _, x in found], np.sort(d[d <= 5.0]), atol=1e-9) + + def test_metric_of_your_own(self): + g = UGraph(dim=2, metric=lambda v: float(np.abs(v).max())) # Chebyshev + for name, c in [("a", (0, 0)), ("b", (2, 2)), ("c", (3, 0)), ("d", (0, 5))]: + g.add_vertex(list(c), name=name) + g.add_vertex(name="nowhere") + found = g.near([0.0, 0.0], 3.0) + self.assertEqual([(v.name, d) for v, d in found], [("a", 0.0), ("b", 2.0), ("c", 3.0)]) + + def test_vertices_without_coordinates(self): + g = UGraph() + g.add_vertex(name="nowhere") + g.add_vertex(name="nowhere else") + self.assertEqual(g.within([[0, 0], [1, 1]]), []) + self.assertEqual(g.near([0, 0], 5), []) + + def test_closest_skips_vertices_without_coordinates(self): + g = UGraph(dim=2) + g.add_vertex(name="nowhere") + g.add_vertex([3.0, 4.0], name="a") + vertex, d = g.closest([0.0, 0.0]) # no index: the loop + self.assertEqual((vertex.name, d), ("a", 5.0)) + + +@unittest.skipUnless(HAVE_SCIPY, "scipy is not installed") +class TestIndexGuards(unittest.TestCase): + def test_a_query_that_is_not_finite_falls_back(self): + g, _ = points_graph(UGraph, n=20) + g.build_kdtree() + vertex, d = g.closest([np.nan, 0.0]) # the loop: nothing is nearer than nan + self.assertIsNone(vertex) + self.assertEqual(d, np.inf) + + def test_a_query_of_the_wrong_length_with_SE2(self): + g, _ = points_graph(UGraph, n=20, dim=3, metric="SE2") + g.build_kdtree() + with self.assertRaises(ValueError): # the metric itself says so + g.closest([0.0, 0.0]) + + def test_coordinates_that_are_none_finite(self): + g = UGraph(dim=2) + g.add_vertex([np.nan, np.nan]) + with self.assertRaises(ValueError): + g.build_kdtree() + + class TestScipyIsOnlyImportedOnDemand(unittest.TestCase): def test_no_module_of_pgraph_imports_scipy_at_import_time(self): # spatialmath, which pgraph depends on, imports scipy itself, so diff --git a/tests/test_subgraph.py b/tests/test_subgraph.py index 765b33d2..d8a9bd65 100644 --- a/tests/test_subgraph.py +++ b/tests/test_subgraph.py @@ -6,7 +6,7 @@ Expectations are written out by hand. """ -import pickle # nosec B403 - only used to round-trip graphs made in this file +import pickle # nosec B403 nosemgrep - only used to round-trip graphs made in this file import unittest import numpy as np @@ -14,6 +14,11 @@ from pgraph import DGraph, Edge, UGraph, Vertex +def roundtrip(obj): + """Pickle and unpickle: these tests only read what they have just written""" + return pickle.loads(pickle.dumps(obj)) # nosec B301 nosemgrep + + def line(cls=UGraph, n=5, **kwargs): """Vertices v0 .. v(n-1) at (0, 0), (1, 0), ... with no edges""" g = cls(dim=2, **kwargs) @@ -284,7 +289,7 @@ def test_with_within(self): def test_can_be_pickled_and_modified(self): g, v = self.graph() h = g.subgraph([v[1], v[2], v[3]]) - h2 = pickle.loads(pickle.dumps(h)) # nosec B301 + h2 = roundtrip(h) self.assertEqual(self.description(h2), self.description(h)) h.add_edge(h["v1"], h["v3"]) self.assertEqual(h.ne, 5)