__init__.py 59.3 KB
Newer Older
Tiago Peixoto's avatar
Tiago Peixoto committed
1
2
3
4
5
#! /usr/bin/env python
# -*- coding: utf-8 -*-
#
# graph_tool -- a general graph manipulation python module
#
Tiago Peixoto's avatar
Tiago Peixoto committed
6
# Copyright (C) 2006-2014 Tiago de Paula Peixoto <tiago@skewed.de>
Tiago Peixoto's avatar
Tiago Peixoto committed
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <http://www.gnu.org/licenses/>.

"""
``graph_tool.search`` - Search algorithms
-----------------------------------------

This module includes several search algorithms, which are customizable to
arbitrary purposes. It is mostly a wrapper around the Visitor interface of the
`Boost Graph Library <http://www.boost.org/doc/libs/release/libs/graph/doc/visitor_concepts.html>`_,
and the respective search functions.


Summary
+++++++

.. autosummary::
   :nosignatures:

   bfs_search
   dfs_search
   dijkstra_search
   bellman_ford_search
   astar_search
   BFSVisitor
   DFSVisitor
   DijkstraVisitor
   BellmanFordVisitor
   AStarVisitor
   StopSearch

Examples
++++++++

In this module, most documentation examples will make use of the network
:download:`search_example.xml <search_example.xml>`, shown below.

>>> g = gt.load_graph("search_example.xml")
Tiago Peixoto's avatar
Tiago Peixoto committed
56
57
58
59
60
61
62
>>> name = g.vp["name"]
>>> weight = g.ep["weight"]
>>> pos = g.vp["pos"]
>>> gt.graph_draw(g, pos, vertex_text=name, vertex_font_size=12, vertex_shape="double_circle",
...               vertex_fill_color="#729fcf", vertex_pen_width=3,
...               edge_pen_width=weight, output="search_example.pdf")
<...>
63
64
65
66
67
68
69

.. testcode::
   :hide:

   gt.graph_draw(g, pos=pos, vertex_text=name, vertex_font_size=12, vertex_shape="double_circle",
                 vertex_fill_color="#729fcf", vertex_pen_width=3,
                 edge_pen_width=weight, output="search_example.png")
Tiago Peixoto's avatar
Tiago Peixoto committed
70

71
.. figure:: search_example.*
Tiago Peixoto's avatar
Tiago Peixoto committed
72
73
74
75
76
77
78
79
80
81
82
   :alt: search example
   :align: center

   This is the network used in the examples below. The width of the edges
   correspond to the values of the "weight" property map.


Contents
++++++++
"""

83
from __future__ import division, absolute_import, print_function
84
85
86
import sys
if sys.version_info < (3,):
    range = xrange
87

Tiago Peixoto's avatar
Tiago Peixoto committed
88
from .. dl_import import dl_import
89
dl_import("from . import libgraph_tool_search")
Tiago Peixoto's avatar
Tiago Peixoto committed
90

91
from .. import _prop, _python_type
Tiago Peixoto's avatar
Tiago Peixoto committed
92
93
94
95
96
97
98
99
import weakref

__all__ = ["bfs_search", "BFSVisitor", "dfs_search", "DFSVisitor",
           "dijkstra_search", "DijkstraVisitor", "bellman_ford_search",
           "BellmanFordVisitor", "astar_search", "AStarVisitor",
           "StopSearch"]


100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
class VisitorWrapper(object):
    def __init__(self, g, visitor, edge_members, vertex_members):
        self.visitor = visitor
        self.g = g
        self.edge_members = set(edge_members)
        self.vertex_members = set(vertex_members)

    def __getattr__(self, attr):
        try:
            orig_attr = self.visitor.__getattribute__(attr)
        except AttributeError:
            return object.__getattribute__(self, attr)
        if callable(orig_attr):
            def wrapped_visitor_member(*args, **kwargs):
                old_perms = dict(self.g._Graph__perms)
                perms ={"del_vertex": False, "del_edge": False, "add_edge": False}
                if attr in self.edge_members:
                    perms.update({"del_edge": True, "add_edge": True})
                elif attr in self.vertex_members:
                    perms.update({"add_vertex": False})
                self.g._Graph__perms.update(perms)
                try:
                    ret = orig_attr(*args, **kwargs)
                finally:
                    self.g._Graph__perms.update(old_perms)
                return ret
            return wrapped_visitor_member
        else:
            return orig_attr
Tiago Peixoto's avatar
Tiago Peixoto committed
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185

class BFSVisitor(object):
    r"""A visitor object that is invoked at the event-points inside the
    :func:`~graph_tool.search.bfs_search` algorithm. By default, it performs no
    action, and should be used as a base class in order to be useful."""

    def initialize_vertex(self, u):
        """This is invoked on every vertex of the graph before the start of the
        graph search. """
        return

    def discover_vertex(self, u):
        """This is invoked when a vertex is encountered for the first time."""
        return

    def examine_vertex(self, u):
        """This is invoked on a vertex as it is popped from the queue. This
        happens immediately before examine_edge() is invoked on each of the
        out-edges of vertex u."""
        return

    def examine_edge(self, e):
        """This is invoked on every out-edge of each vertex after it is
        discovered."""
        return

    def tree_edge(self, e):
        """This is invoked on each edge as it becomes a member of the edges that
        form the search tree."""
        return

    def non_tree_edge(self, e):
        """This is invoked on back or cross edges for directed graphs and cross
        edges for undirected graphs. """
        return

    def gray_target(self, e):
        """This is invoked on the subset of non-tree edges whose target vertex
        is colored gray at the time of examination. The color gray indicates
        that the vertex is currently in the queue."""
        return

    def black_target(self, e):
        """This is invoked on the subset of non-tree edges whose target vertex
        is colored black at the time of examination. The color black indicates
        that the vertex has been removed from the queue."""
        return

    def finish_vertex(self, u):
        """This invoked on a vertex after all of its out edges have been added
        to the search tree and all of the adjacent vertices have been discovered
        (but before the out-edges of the adjacent vertices have been examined).
        """
        return


def bfs_search(g, source, visitor=BFSVisitor()):
186
    r"""Breadth-first traversal of a directed or undirected graph.
Tiago Peixoto's avatar
Tiago Peixoto committed
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262

    Parameters
    ----------
    g : :class:`~graph_tool.Graph`
        Graph to be used.
    source : :class:`~graph_tool.Vertex`
        Source vertex.
    visitor : :class:`~graph_tool.search.BFSVisitor` (optional, default: ``BFSVisitor()``)
        A visitor object that is invoked at the event points inside the
        algorithm. This should be a subclass of
        :class:`~graph_tool.search.BFSVisitor`.

    See Also
    --------
    dfs_search: Depth-first search
    dijkstra_search: Dijkstra's search algorithm
    astar_search: :math:`A^*` heuristic search algorithm

    Notes
    -----

    A breadth-first traversal visits vertices that are closer to the source
    before visiting vertices that are further away. In this context "distance"
    is defined as the number of edges in the shortest path from the source
    vertex.

    The time complexity is :math:`O(V + E)`.

    The pseudo-code for the BFS algorithm is listed below, with the annotated
    event points, for which the given visitor object will be called with the
    appropriate method.

    ::

        BFS(G, source)
         for each vertex u in V[G]       initialize vertex u
           color[u] := WHITE
           d[u] := infinity
         end for
         color[source] := GRAY
         d[source] := 0
         ENQUEUE(Q, source)              discover vertex source
         while (Q != Ø)
           u := DEQUEUE(Q)               examine vertex u
           for each vertex v in Adj[u]   examine edge (u,v)
             if (color[v] = WHITE)       (u,v) is a tree edge
               color[v] := GRAY
               ENQUEUE(Q, v)             discover vertex v
             else                        (u,v) is a non-tree edge
               if (color[v] = GRAY)
                 ...                     (u,v) has a gray target
               else
                 ...                     (u,v) has a black target
           end for
           color[u] := BLACK             finish vertex u
         end while


    Examples
    --------

    We must define what should be done during the search by subclassing
    :class:`~graph_tool.search.BFSVisitor`, and specializing the appropriate
    methods. In the following we will keep track of the distance from the root,
    and the predecessor tree.

    .. testcode::

        class VisitorExample(gt.BFSVisitor):

            def __init__(self, name, pred, dist):
                self.name = name
                self.pred = pred
                self.dist = dist

            def discover_vertex(self, u):
263
                print("-->", self.name[u], "has been discovered!")
Tiago Peixoto's avatar
Tiago Peixoto committed
264
265

            def examine_vertex(self, u):
266
                print(self.name[u], "has been examined...")
Tiago Peixoto's avatar
Tiago Peixoto committed
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296

            def tree_edge(self, e):
                self.pred[e.target()] = int(e.source())
                self.dist[e.target()] = self.dist[e.source()] + 1

    With the above class defined, we can perform the BFS search as follows.

    >>> dist = g.new_vertex_property("int")
    >>> pred = g.new_vertex_property("int")
    >>> gt.bfs_search(g, g.vertex(0), VisitorExample(name, pred, dist))
    --> Bob has been discovered!
    Bob has been examined...
    --> Eve has been discovered!
    --> Chuck has been discovered!
    --> Carlos has been discovered!
    --> Isaac has been discovered!
    Eve has been examined...
    --> Imothep has been discovered!
    --> Carol has been discovered!
    Chuck has been examined...
    Carlos has been examined...
    --> Alice has been discovered!
    Isaac has been examined...
    Imothep has been examined...
    Carol has been examined...
    Alice has been examined...
    --> Oscar has been discovered!
    --> Dave has been discovered!
    Oscar has been examined...
    Dave has been examined...
297
    >>> print(dist.a)
Tiago Peixoto's avatar
Tiago Peixoto committed
298
    [0 2 2 1 1 3 1 1 3 2]
299
    >>> print(pred.a)
Tiago Peixoto's avatar
Tiago Peixoto committed
300
301
302
303
304
305
306
307
308
309
310
    [0 3 6 0 0 1 0 0 1 6]


    References
    ----------
    .. [bfs] Edward Moore, "The shortest path through a maze", International
             Symposium on the Theory of Switching, 1959
    .. [bfs-bgl] http://www.boost.org/doc/libs/release/libs/graph/doc/breadth_first_search.html
    .. [bfs-wikipedia] http://en.wikipedia.org/wiki/Breadth-first_search
    """

311
312
313
    visitor = VisitorWrapper(g, visitor,
                             ["initialize_vertex", "examine_vertex", "finish_vertex"],
                             ["initialize_vertex"])
Tiago Peixoto's avatar
Tiago Peixoto committed
314
315
316

    try:
        libgraph_tool_search.bfs_search(g._Graph__graph,
317
                                        weakref.ref(g),
Tiago Peixoto's avatar
Tiago Peixoto committed
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
                                        int(source), visitor)
    except StopSearch:
        pass


class DFSVisitor(object):
    r"""
    A visitor object that is invoked at the event-points inside the
    :func:`~graph_tool.search.dfs_search` algorithm. By default, it performs no
    action, and should be used as a base class in order to be useful.
    """

    def initialize_vertex(self, u):
        """
        This is invoked on every vertex of the graph before the start of the
        graph search.
        """
        return

    def start_vertex(self, u):
        """
        This is invoked on the source vertex once before the start of the
        search.
        """
        return

    def discover_vertex(self, u):
        """This is invoked when a vertex is encountered for the first time."""
        return

    def examine_edge(self, e):
        """
        This is invoked on every out-edge of each vertex after it is discovered.
        """
        return

    def tree_edge(self, e):
        """
        This is invoked on each edge as it becomes a member of the edges that
        form the search tree.
        """
        return

    def back_edge(self, e):
        """
        This is invoked on the back edges in the graph. For an undirected graph
        there is some ambiguity between tree edges and back edges since the edge
        (u,v) and (v,u) are the same edge, but both the
        :meth:`~graph_tool.search.DFSVisitor.tree_edge` and
        :meth:`~graph_tool.search..DFSVisitor.back_edge` functions will be
        invoked. One way to resolve this ambiguity is to record the tree edges,
        and then disregard the back-edges that are already marked as tree
        edges. An easy way to record tree edges is to record predecessors at the
        tree_edge event point.
        """
        return

    def forward_or_cross_edge(self, e):
        """
        This is invoked on forward or cross edges in the graph. In an undirected
        graph this method is never called.
        """
        return

    def finish_vertex(self, e):
        """
        This is invoked on vertex u after finish_vertex has been called for all
        the vertices in the DFS-tree rooted at vertex u. If vertex u is a leaf
        in the DFS-tree, then the
        :meth:`~graph_tool..search.DFSVisitor.finish_vertex` function is called
        on u after all the out-edges of u have been examined.
        """
        return


def dfs_search(g, source, visitor=DFSVisitor()):
394
    r"""Depth-first traversal of a directed or undirected graph.
Tiago Peixoto's avatar
Tiago Peixoto committed
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470

    Parameters
    ----------
    g : :class:`~graph_tool.Graph`
        Graph to be used.
    source : :class:`~graph_tool.Vertex`
        Source vertex.
    visitor : :class:`~graph_tool.search.DFSVisitor` (optional, default: ``DFSVisitor()``)
        A visitor object that is invoked at the event points inside the
        algorithm. This should be a subclass of
        :class:`~graph_tool.search.DFSVisitor`.

    See Also
    --------
    bfs_search: Breadth-first search
    dijkstra_search: Dijkstra's search algorithm
    astar_search: :math:`A^*` heuristic search algorithm

    Notes
    -----

    When possible, a depth-first traversal chooses a vertex adjacent to the
    current vertex to visit next. If all adjacent vertices have already been
    discovered, or there are no adjacent vertices, then the algorithm backtracks
    to the last vertex that had undiscovered neighbors. Once all reachable
    vertices have been visited, the algorithm selects from any remaining
    undiscovered vertices and continues the traversal. The algorithm finishes
    when all vertices have been visited.

    The time complexity is :math:`O(V + E)`.

    The pseudo-code for the DFS algorithm is listed below, with the annotated
    event points, for which the given visitor object will be called with the
    appropriate method.

    ::

        DFS(G)
          for each vertex u in V
            color[u] := WHITE                 initialize vertex u
          end for
          time := 0
          call DFS-VISIT(G, source)           start vertex s

        DFS-VISIT(G, u)
          color[u] := GRAY                    discover vertex u
          for each v in Adj[u]                examine edge (u,v)
            if (color[v] = WHITE)             (u,v) is a tree edge
              call DFS-VISIT(G, v)
            else if (color[v] = GRAY)         (u,v) is a back edge
              ...
            else if (color[v] = BLACK)        (u,v) is a cross or forward edge
              ...
          end for
          color[u] := BLACK                   finish vertex u


    Examples
    --------

    We must define what should be done during the search by subclassing
    :class:`~graph_tool.search.DFSVisitor`, and specializing the appropriate
    methods. In the following we will keep track of the discover time, and the
    predecessor tree.

    .. testcode::

        class VisitorExample(gt.DFSVisitor):

            def __init__(self, name, pred, time):
                self.name = name
                self.pred = pred
                self.time = time
                self.last_time = 0

            def discover_vertex(self, u):
471
                print("-->", self.name[u], "has been discovered!")
Tiago Peixoto's avatar
Tiago Peixoto committed
472
473
474
475
                self.time[u] = self.last_time
                self.last_time += 1

            def examine_edge(self, e):
476
477
                print("edge (%s, %s) has been examined..." % \
                    (self.name[e.source()], self.name[e.target()]))
Tiago Peixoto's avatar
Tiago Peixoto committed
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531

            def tree_edge(self, e):
                self.pred[e.target()] = int(e.source())


    With the above class defined, we can perform the DFS search as follows.

    >>> time = g.new_vertex_property("int")
    >>> pred = g.new_vertex_property("int")
    >>> gt.dfs_search(g, g.vertex(0), VisitorExample(name, pred, time))
    --> Bob has been discovered!
    edge (Bob, Eve) has been examined...
    --> Eve has been discovered!
    edge (Eve, Isaac) has been examined...
    --> Isaac has been discovered!
    edge (Isaac, Bob) has been examined...
    edge (Isaac, Chuck) has been examined...
    --> Chuck has been discovered!
    edge (Chuck, Eve) has been examined...
    edge (Chuck, Isaac) has been examined...
    edge (Chuck, Imothep) has been examined...
    --> Imothep has been discovered!
    edge (Imothep, Carol) has been examined...
    --> Carol has been discovered!
    edge (Carol, Eve) has been examined...
    edge (Carol, Imothep) has been examined...
    edge (Imothep, Carlos) has been examined...
    --> Carlos has been discovered!
    edge (Carlos, Eve) has been examined...
    edge (Carlos, Imothep) has been examined...
    edge (Carlos, Bob) has been examined...
    edge (Carlos, Alice) has been examined...
    --> Alice has been discovered!
    edge (Alice, Oscar) has been examined...
    --> Oscar has been discovered!
    edge (Oscar, Alice) has been examined...
    edge (Oscar, Dave) has been examined...
    --> Dave has been discovered!
    edge (Dave, Oscar) has been examined...
    edge (Dave, Alice) has been examined...
    edge (Alice, Dave) has been examined...
    edge (Alice, Carlos) has been examined...
    edge (Imothep, Chuck) has been examined...
    edge (Imothep, Eve) has been examined...
    edge (Chuck, Bob) has been examined...
    edge (Isaac, Eve) has been examined...
    edge (Eve, Imothep) has been examined...
    edge (Eve, Bob) has been examined...
    edge (Eve, Carol) has been examined...
    edge (Eve, Carlos) has been examined...
    edge (Eve, Chuck) has been examined...
    edge (Bob, Chuck) has been examined...
    edge (Bob, Carlos) has been examined...
    edge (Bob, Isaac) has been examined...
532
    >>> print(time.a)
Tiago Peixoto's avatar
Tiago Peixoto committed
533
    [0 7 5 6 3 9 1 2 8 4]
534
    >>> print(pred.a)
Tiago Peixoto's avatar
Tiago Peixoto committed
535
536
537
538
539
540
541
542
543
544
    [0 3 9 9 7 8 0 6 1 4]



    References
    ----------
    .. [dfs-bgl] http://www.boost.org/doc/libs/release/libs/graph/doc/depth_first_search.html
    .. [dfs-wikipedia] http://en.wikipedia.org/wiki/Depth-first_search
    """

545
546
547
    visitor = VisitorWrapper(g, visitor,
                             ["initialize_vertex", "discover_vertex", "finish_vertex",
                              "start_vertex"], ["initialize_vertex"])
Tiago Peixoto's avatar
Tiago Peixoto committed
548
549
550

    try:
        libgraph_tool_search.dfs_search(g._Graph__graph,
551
                                        weakref.ref(g),
Tiago Peixoto's avatar
Tiago Peixoto committed
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
                                        int(source), visitor)
    except StopSearch:
        pass


class DijkstraVisitor(object):
    r"""A visitor object that is invoked at the event-points inside the
    :func:`~graph_tool.search.dijkstra_search` algorithm. By default, it
    performs no action, and should be used as a base class in order to be
    useful.
    """

    def initialize_vertex(self, u):
        """
        This is invoked on every vertex of the graph before the start of the
        graph search.
        """
        return

    def examine_vertex(self, u):
        """
        This is invoked on a vertex as it is popped from the queue. This happens
        immediately before :meth:`~graph_tool.DijsktraVisitor.examine_edge` is
        invoked on each of the out-edges of vertex u.
        """
        return

    def examine_edge(self, e):
        """
        This is invoked on every out-edge of each vertex after it is discovered.
        """
        return

    def discover_vertex(self, u):
        """This is invoked when a vertex is encountered for the first time."""
        return

    def edge_relaxed(self, e):
        """
        Upon examination, if the following condition holds then the edge is
        relaxed (its distance is reduced), and this method is invoked.

        ::

            (u, v) = tuple(e)
            assert(compare(combine(d[u], weight[e]), d[v]))

        """
        return

    def edge_not_relaxed(self, e):
        """
        Upon examination, if the edge is not relaxed (see
        :meth:`~graph_tool.search.DijsktraVisitor.edge_relaxed`) then this
        method is invoked.
        """
        return

    def finish_vertex(self, u):
        """
        This invoked on a vertex after all of its out edges have been added to
        the search tree and all of the adjacent vertices have been discovered
        (but before their out-edges have been examined).
        """
        return


def dijkstra_search(g, source, weight, visitor=DijkstraVisitor(), dist_map=None,
                    pred_map=None, combine=lambda a, b: a + b,
                    compare=lambda a, b: a < b, zero=0, infinity=float('inf')):
622
    r"""Dijsktra traversal of a directed or undirected graph, with non-negative weights.
Tiago Peixoto's avatar
Tiago Peixoto committed
623
624
625
626
627
628
629
630
631

    Parameters
    ----------
    g : :class:`~graph_tool.Graph`
        Graph to be used.
    source : :class:`~graph_tool.Vertex`
        Source vertex.
    weight : :class:`~graph_tool.PropertyMap`
        Edge property map with weight values.
Tiago Peixoto's avatar
Tiago Peixoto committed
632
    visitor : :class:`~graph_tool.search.DijkstraVisitor` (optional, default: ``DijkstraVisitor()``)
Tiago Peixoto's avatar
Tiago Peixoto committed
633
634
635
        A visitor object that is invoked at the event points inside the
        algorithm. This should be a subclass of
        :class:`~graph_tool.search.DijkstraVisitor`.
Tiago Peixoto's avatar
Tiago Peixoto committed
636
    dist_map : :class:`~graph_tool.PropertyMap` (optional, default: ``None``)
Tiago Peixoto's avatar
Tiago Peixoto committed
637
638
        A vertex property map where the distances from the source will be
        stored.
Tiago Peixoto's avatar
Tiago Peixoto committed
639
    pred_map : :class:`~graph_tool.PropertyMap` (optional, default: ``None``)
Tiago Peixoto's avatar
Tiago Peixoto committed
640
641
        A vertex property map where the predecessor map will be
        stored (must have value type "int").
Tiago Peixoto's avatar
Tiago Peixoto committed
642
    combine : binary function (optional, default: ``lambda a, b: a + b``)
Tiago Peixoto's avatar
Tiago Peixoto committed
643
644
        This function is used to combine distances to compute the distance of a
        path.
Tiago Peixoto's avatar
Tiago Peixoto committed
645
    compare : binary function (optional, default: ``lambda a, b: a < b``)
Tiago Peixoto's avatar
Tiago Peixoto committed
646
647
        This function is use to compare distances to determine which vertex is
        closer to the source vertex.
648
649
650
651
652
653
    zero : int or float (optional, default: ``0``)
         Value assumed to correspond to a distance of zero by the combine and
         compare functions.
    infinity : int or float (optional, default: ``float('inf')``)
         Value assumed to correspond to a distance of infinity by the combine and
         compare functions.
Tiago Peixoto's avatar
Tiago Peixoto committed
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731

    Returns
    -------
    dist_map : :class:`~graph_tool.PropertyMap`
        A vertex property map with the computed distances from the source.
    pred_map : :class:`~graph_tool.PropertyMap`
        A vertex property map with the predecessor tree.

    See Also
    --------
    bfs_search: Breadth-first search
    dfs_search: Depth-first search
    astar_search: :math:`A^*` heuristic search algorithm

    Notes
    -----

    Dijkstra's algorithm finds all the shortest paths from the source vertex to
    every other vertex by iteratively "growing" the set of vertices S to which
    it knows the shortest path. At each step of the algorithm, the next vertex
    added to S is determined by a priority queue. The queue contains the
    vertices in V - S prioritized by their distance label, which is the length
    of the shortest path seen so far for each vertex. The vertex u at the top of
    the priority queue is then added to S, and each of its out-edges is relaxed:
    if the distance to u plus the weight of the out-edge (u,v) is less than the
    distance label for v then the estimated distance for vertex v is
    reduced. The algorithm then loops back, processing the next vertex at the
    top of the priority queue. The algorithm finishes when the priority queue is
    empty.

    The time complexity is :math:`O(V \log V)`.

    The pseudo-code for Dijkstra's algorithm is listed below, with the annotated
    event points, for which the given visitor object will be called with the
    appropriate method.

    ::

        DIJKSTRA(G, source, weight)
          for each vertex u in V                      initialize vertex u
            d[u] := infinity
            p[u] := u
          end for
          d[source] := 0
          INSERT(Q, source)                           discover vertex s
          while (Q != Ø)
            u := EXTRACT-MIN(Q)                       examine vertex u
            for each vertex v in Adj[u]               examine edge (u,v)
              if (weight[(u,v)] + d[u] < d[v])        edge (u,v) relaxed
                d[v] := weight[(u,v)] + d[u]
                p[v] := u
                DECREASE-KEY(Q, v)
              else                                    edge (u,v) not relaxed
                ...
              if (d[v] was originally infinity)
                INSERT(Q, v)                          discover vertex v
            end for                                   finish vertex u
          end while
          return d

    Examples
    --------

    We must define what should be done during the search by subclassing
    :class:`~graph_tool.search.DijkstraVisitor`, and specializing the
    appropriate methods. In the following we will keep track of the discover
    time, and the predecessor tree.

    .. testcode::

        class VisitorExample(gt.DijkstraVisitor):

            def __init__(self, name, time):
                self.name = name
                self.time = time
                self.last_time = 0

            def discover_vertex(self, u):
732
                print("-->", self.name[u], "has been discovered!")
Tiago Peixoto's avatar
Tiago Peixoto committed
733
734
735
736
                self.time[u] = self.last_time
                self.last_time += 1

            def examine_edge(self, e):
737
738
                print("edge (%s, %s) has been examined..." % \
                    (self.name[e.source()], self.name[e.target()]))
Tiago Peixoto's avatar
Tiago Peixoto committed
739
740

            def edge_relaxed(self, e):
741
742
                print("edge (%s, %s) has been relaxed..." % \
                    (self.name[e.source()], self.name[e.target()]))
Tiago Peixoto's avatar
Tiago Peixoto committed
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801


    With the above class defined, we can perform the Dijkstra search as follows.

    >>> time = g.new_vertex_property("int")
    >>> dist, pred = gt.dijkstra_search(g, g.vertex(0), weight, VisitorExample(name, time))
    --> Bob has been discovered!
    edge (Bob, Eve) has been examined...
    edge (Bob, Eve) has been relaxed...
    --> Eve has been discovered!
    edge (Bob, Chuck) has been examined...
    edge (Bob, Chuck) has been relaxed...
    --> Chuck has been discovered!
    edge (Bob, Carlos) has been examined...
    edge (Bob, Carlos) has been relaxed...
    --> Carlos has been discovered!
    edge (Bob, Isaac) has been examined...
    edge (Bob, Isaac) has been relaxed...
    --> Isaac has been discovered!
    edge (Eve, Isaac) has been examined...
    edge (Eve, Imothep) has been examined...
    edge (Eve, Imothep) has been relaxed...
    --> Imothep has been discovered!
    edge (Eve, Bob) has been examined...
    edge (Eve, Carol) has been examined...
    edge (Eve, Carol) has been relaxed...
    --> Carol has been discovered!
    edge (Eve, Carlos) has been examined...
    edge (Eve, Chuck) has been examined...
    edge (Isaac, Bob) has been examined...
    edge (Isaac, Chuck) has been examined...
    edge (Isaac, Eve) has been examined...
    edge (Chuck, Eve) has been examined...
    edge (Chuck, Isaac) has been examined...
    edge (Chuck, Imothep) has been examined...
    edge (Chuck, Bob) has been examined...
    edge (Carlos, Eve) has been examined...
    edge (Carlos, Imothep) has been examined...
    edge (Carlos, Bob) has been examined...
    edge (Carlos, Alice) has been examined...
    edge (Carlos, Alice) has been relaxed...
    --> Alice has been discovered!
    edge (Imothep, Carol) has been examined...
    edge (Imothep, Carlos) has been examined...
    edge (Imothep, Chuck) has been examined...
    edge (Imothep, Eve) has been examined...
    edge (Alice, Oscar) has been examined...
    edge (Alice, Oscar) has been relaxed...
    --> Oscar has been discovered!
    edge (Alice, Dave) has been examined...
    edge (Alice, Dave) has been relaxed...
    --> Dave has been discovered!
    edge (Alice, Carlos) has been examined...
    edge (Carol, Eve) has been examined...
    edge (Carol, Imothep) has been examined...
    edge (Oscar, Alice) has been examined...
    edge (Oscar, Dave) has been examined...
    edge (Dave, Oscar) has been examined...
    edge (Dave, Alice) has been examined...
802
    >>> print(time.a)
Tiago Peixoto's avatar
Tiago Peixoto committed
803
    [0 7 6 3 2 9 1 4 8 5]
804
    >>> print(pred.a)
Tiago Peixoto's avatar
Tiago Peixoto committed
805
    [0 3 6 0 0 1 0 0 1 6]
806
    >>> print(dist.a)
Tiago Peixoto's avatar
Tiago Peixoto committed
807
808
809
810
811
812
813
814
815
816
817
    [  0.           8.91915887   9.27141329   4.29277116   4.02118246
      12.23513866   3.23790211   3.45487436  11.04391549   7.74858396]

    References
    ----------
    .. [dijkstra] E. Dijkstra, "A note on two problems in connexion with
        graphs", Numerische Mathematik, 1:269-271, 1959.
    .. [dijkstra-bgl] http://www.boost.org/doc/libs/release/libs/graph/doc/dijkstra_shortest_paths_no_color_map.html
    .. [dijkstra-wikipedia] http://en.wikipedia.org/wiki/Dijkstra's_algorithm
    """

818
819
820
821
    visitor = VisitorWrapper(g, visitor,
                             ["initialize_vertex", "examine_vertex", "finish_vertex"],
                             ["initialize_vertex"])

Tiago Peixoto's avatar
Tiago Peixoto committed
822
823
824
825
826
827
828
829
830
831
    if visitor is None:
        visitor = DijkstraVisitor()
    if dist_map is None:
        dist_map = g.new_vertex_property(weight.value_type())
    if pred_map is None:
        pred_map = g.new_vertex_property("int")
    if pred_map.value_type() != "int32_t":
        raise ValueError("pred_map must be of value type 'int32_t', not '%s'." % \
                             pred_map.value_type())

832
    try:
833
834
        if dist_map.value_type() != "python::object":
            zero = _python_type(dist_map.value_type())(zero)
835
836
837
838
839
    except OverflowError:
        zero = (weight.a.max() + 1) * g.num_vertices()
        zero = _python_type(dist_map.value_type())(zero)

    try:
840
841
        if dist_map.value_type() != "python::object":
            infinity = _python_type(dist_map.value_type())(infinity)
842
843
844
845
    except OverflowError:
        infinity = (weight.a.max() + 1) * g.num_vertices()
        infinity = _python_type(dist_map.value_type())(infinity)

Tiago Peixoto's avatar
Tiago Peixoto committed
846
847
    try:
        libgraph_tool_search.dijkstra_search(g._Graph__graph,
848
                                             weakref.ref(g),
Tiago Peixoto's avatar
Tiago Peixoto committed
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
                                             int(source),
                                             _prop("v", g, dist_map),
                                             _prop("v", g, pred_map),
                                             _prop("e", g, weight), visitor,
                                             compare, combine, zero, infinity)
    except StopSearch:
        pass

    return dist_map, pred_map


class BellmanFordVisitor(object):
    r"""A visitor object that is invoked at the event-points inside the
    :func:`~graph_tool.search.bellman_ford_search` algorithm. By default, it
    performs no action, and should be used as a base class in order to be
    useful.
    """

    def examine_edge(self, e):
        """
        This is invoked on every edge in the graph ``|V|`` times.
        """
        return

    def edge_relaxed(self, e):
        """
        This is invoked when the distance label for the target vertex is
        decreased. The edge (u,v) that participated in the last relaxation for
        vertex v is an edge in the shortest paths tree.
        """
        return

    def edge_not_relaxed(self, e):
        """
        This is invoked if the distance label for the target vertex is not
        decreased.
        """
        return

    def edge_minimized(self, e):
        """
        This is invoked during the second stage of the algorithm, during the
        test of whether each edge was minimized. If the edge is minimized then
        this function is invoked.
        """
        return

    def edge_not_minimized(self, e):
        """
        This is invoked during the second stage of the algorithm, during the
        test of whether each edge was minimized. If the edge was not minimized,
        this function is invoked. This happens when there is a negative cycle in
        the graph.
        """
        return


def bellman_ford_search(g, source, weight, visitor=BellmanFordVisitor(),
                        dist_map=None, pred_map=None,
                        combine=lambda a, b: a + b,
                        compare=lambda a, b: a < b, zero=0,
                        infinity=float('inf')):
911
    r"""Bellman-Ford traversal of a directed or undirected graph, with negative weights.
Tiago Peixoto's avatar
Tiago Peixoto committed
912
913
914
915
916
917
918
919
920

    Parameters
    ----------
    g : :class:`~graph_tool.Graph`
        Graph to be used.
    source : :class:`~graph_tool.Vertex`
        Source vertex.
    weight : :class:`~graph_tool.PropertyMap`
        Edge property map with weight values.
Tiago Peixoto's avatar
Tiago Peixoto committed
921
    visitor : :class:`~graph_tool.search.DijkstraVisitor` (optional, default: ``DijkstraVisitor()``)
Tiago Peixoto's avatar
Tiago Peixoto committed
922
923
924
        A visitor object that is invoked at the event points inside the
        algorithm. This should be a subclass of
        :class:`~graph_tool.search.DijkstraVisitor`.
Tiago Peixoto's avatar
Tiago Peixoto committed
925
    dist_map : :class:`~graph_tool.PropertyMap` (optional, default: ``None``)
Tiago Peixoto's avatar
Tiago Peixoto committed
926
927
        A vertex property map where the distances from the source will be
        stored.
Tiago Peixoto's avatar
Tiago Peixoto committed
928
    pred_map : :class:`~graph_tool.PropertyMap` (optional, default: ``None``)
Tiago Peixoto's avatar
Tiago Peixoto committed
929
930
        A vertex property map where the predecessor map will be
        stored (must have value type "int").
Tiago Peixoto's avatar
Tiago Peixoto committed
931
    combine : binary function (optional, default: ``lambda a, b: a + b``)
Tiago Peixoto's avatar
Tiago Peixoto committed
932
933
        This function is used to combine distances to compute the distance of a
        path.
Tiago Peixoto's avatar
Tiago Peixoto committed
934
    compare : binary function (optional, default: ``lambda a, b: a < b``)
Tiago Peixoto's avatar
Tiago Peixoto committed
935
936
        This function is use to compare distances to determine which vertex is
        closer to the source vertex.
937
938
939
940
941
942
    zero : int or float (optional, default: ``0``)
         Value assumed to correspond to a distance of zero by the combine and
         compare functions.
    infinity : int or float (optional, default: ``float('inf')``)
         Value assumed to correspond to a distance of infinity by the combine and
         compare functions.
Tiago Peixoto's avatar
Tiago Peixoto committed
943
944
945
946
947


    Returns
    -------
    minimized : bool
Tiago Peixoto's avatar
Tiago Peixoto committed
948
        True if all edges were successfully minimized, or False if there is a
Tiago Peixoto's avatar
Tiago Peixoto committed
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
1001
1002
1003
1004
1005
1006
1007
1008
1009
1010
1011
1012
1013
1014
1015
1016
1017
1018
1019
1020
1021
1022
1023
1024
1025
1026
1027
1028
1029
1030
1031
1032
1033
1034
1035
1036
1037
1038
        negative loop in the graph.
    dist_map : :class:`~graph_tool.PropertyMap`
        A vertex property map with the computed distances from the source.
    pred_map : :class:`~graph_tool.PropertyMap`
        A vertex property map with the predecessor tree.

    See Also
    --------
    bfs_search: Breadth-first search
    dfs_search: Depth-first search
    dijsktra_search: Dijkstra search
    astar_search: :math:`A^*` heuristic search

    Notes
    -----

    The Bellman-Ford algorithm [bellman-ford]_ solves the single-source shortest
    paths problem for a graph with both positive and negative edge weights. If
    you only need to solve the shortest paths problem for positive edge weights,
    :func:`~graph_tool.search.dijkstra_search` provides a more efficient
    alternative. If all the edge weights are all equal, then
    :func:`~graph_tool.search.bfs_search` provides an even more efficient
    alternative.

    The Bellman-Ford algorithm proceeds by looping through all of the edges in
    the graph, applying the relaxation operation to each edge. In the following
    pseudo-code, ``v`` is a vertex adjacent to ``u``, ``w`` maps edges to their
    weight, and ``d`` is a distance map that records the length of the shortest
    path to each vertex seen so far. ``p`` is a predecessor map which records
    the parent of each vertex, which will ultimately be the parent in the
    shortest paths tree

    ::

        RELAX(u, v, w, d, p)
          if (w(u,v) + d[u] < d[v])
            d[v] := w(u,v) + d[u]       relax edge (u,v)
            p[v] := u
          else
            ...                         edge (u,v) is not relaxed


    The algorithm repeats this loop ``|V|`` times after which it is guaranteed
    that the distances to each vertex have been reduced to the minimum possible
    unless there is a negative cycle in the graph. If there is a negative cycle,
    then there will be edges in the graph that were not properly minimized. That
    is, there will be edges ``(u,v)`` such that ``w(u,v) + d[u] < d[v]``. The
    algorithm loops over the edges in the graph one final time to check if all
    the edges were minimized, returning true if they were and returning false
    otherwise.

    ::


        BELLMAN-FORD(G)
          for each vertex u in V
            d[u] := infinity
            p[u] := u
          end for
          for i := 1 to |V|-1
            for each edge (u,v) in E          examine edge (u,v)
              RELAX(u, v, w, d, p)
            end for
          end for
          for each edge (u,v) in E
            if (w(u,v) + d[u] < d[v])
              return (false, , )              edge (u,v) was not minimized
            else
              ...                             edge (u,v) was minimized
          end for
          return (true, p, d)

    The time complexity is :math:`O(V E)`.

    Examples
    --------

    We must define what should be done during the search by subclassing
    :class:`~graph_tool.search.BellmanFordVisitor`, and specializing the
    appropriate methods. In the following we will keep track of the edge
    minimizations.

    .. testcode::

        class VisitorExample(gt.BellmanFordVisitor):

            def __init__(self, name):
                self.name = name

            def edge_minimized(self, e):
1039
1040
                print("edge (%s, %s) has been minimized..." % \
                      (self.name[e.source()], self.name[e.target()]))
Tiago Peixoto's avatar
Tiago Peixoto committed
1041
1042

            def edge_not_minimized(self, e):
1043
1044
                print("edge (%s, %s) has not been minimized..." % \
                      (self.name[e.source()], self.name[e.target()]))
Tiago Peixoto's avatar
Tiago Peixoto committed
1045
1046
1047
1048
1049
1050
1051
1052
1053
1054
1055
1056
1057
1058
1059
1060
1061
1062
1063
1064
1065
1066
1067
1068
1069


    With the above class defined, we can perform the Bellman-Ford search as
    follows.

    >>> nweight = g.copy_property(weight)
    >>> nweight.a[6] *= -1   # include negative weight in edge (Carlos, Alice)
    >>> minimized, dist, pred = gt.bellman_ford_search(g, g.vertex(0), nweight, VisitorExample(name))
    edge (Bob, Eve) has been minimized...
    edge (Bob, Chuck) has been minimized...
    edge (Bob, Carlos) has been minimized...
    edge (Bob, Isaac) has been minimized...
    edge (Alice, Oscar) has been minimized...
    edge (Alice, Dave) has been minimized...
    edge (Alice, Carlos) has been minimized...
    edge (Carol, Eve) has been minimized...
    edge (Carol, Imothep) has been minimized...
    edge (Carlos, Eve) has been minimized...
    edge (Carlos, Imothep) has been minimized...
    edge (Chuck, Eve) has been minimized...
    edge (Chuck, Isaac) has been minimized...
    edge (Chuck, Imothep) has been minimized...
    edge (Dave, Oscar) has been minimized...
    edge (Eve, Isaac) has been minimized...
    edge (Eve, Imothep) has been minimized...
1070
    >>> print(minimized)
Tiago Peixoto's avatar
Tiago Peixoto committed
1071
    True
1072
    >>> print(pred.a)
Tiago Peixoto's avatar
Tiago Peixoto committed
1073
    [3 3 9 1 6 1 3 6 1 3]
1074
    >>> print(dist.a)
Tiago Peixoto's avatar
Tiago Peixoto committed
1075
1076
1077
1078
1079
1080
1081
1082
1083
1084
1085
    [-28.42555934 -37.34471821 -25.20438243 -41.97110592 -35.20316571
     -34.02873843 -36.58860946 -33.55645565 -35.2199616  -36.0029274 ]

    References
    ----------
    .. [bellman-ford] R. Bellman, "On a routing problem", Quarterly of Applied
       Mathematics, 16(1):87-90, 1958.
    .. [bellman-ford-bgl] http://www.boost.org/doc/libs/release/libs/graph/doc/bellman_ford_shortest.html
    .. [bellman-ford-wikipedia] http://en.wikipedia.org/wiki/Bellman-Ford_algorithm
    """

1086
1087
    visitor = VisitorWrapper(g, visitor, [], [])

Tiago Peixoto's avatar
Tiago Peixoto committed
1088
1089
1090
1091
1092
1093
1094
    if dist_map is None:
        dist_map = g.new_vertex_property(weight.value_type())
    if pred_map is None:
        pred_map = g.new_vertex_property("int")
    if pred_map.value_type() != "int32_t":
        raise ValueError("pred_map must be of value type 'int32_t', not '%s'." % \
                             pred_map.value_type())
1095
1096

    try:
1097
1098
        if dist_map.value_type() != "python::object":
            zero = _python_type(dist_map.value_type())(zero)
1099
1100
1101
1102
1103
    except OverflowError:
        zero = (weight.a.max() + 1) * g.num_vertices()
        zero = _python_type(dist_map.value_type())(zero)

    try:
1104
1105
        if dist_map.value_type() != "python::object":
            infinity = _python_type(dist_map.value_type())(infinity)
1106
1107
1108
1109
    except OverflowError:
        infinity = (weight.a.max() + 1) * g.num_vertices()
        infinity = _python_type(dist_map.value_type())(infinity)

Tiago Peixoto's avatar
Tiago Peixoto committed
1110
1111
1112
1113
    minimized = False
    try:
        minimized = \
            libgraph_tool_search.bellman_ford_search(g._Graph__graph,
1114
                                                     weakref.ref(g),
Tiago Peixoto's avatar
Tiago Peixoto committed
1115
1116
1117
1118
1119
1120
1121
1122
1123
1124
1125
1126
1127
1128
1129
1130
1131
1132
1133
1134
1135
1136
1137
1138
1139
1140
1141
1142
1143
1144
1145
1146
1147
1148
1149
1150
1151
1152
1153
1154
1155
1156
1157
1158
1159
1160
1161
1162
1163
1164
1165
1166
1167
1168
1169
1170
1171
1172
1173
1174
1175
1176
1177
1178
1179
1180
1181
1182
1183
1184
1185
1186
1187
1188
1189
1190
1191
1192
1193
1194
1195
1196
1197
1198
1199
1200
1201
1202
                                                     int(source),
                                                     _prop("v", g, dist_map),
                                                     _prop("v", g, pred_map),
                                                     _prop("e", g, weight),
                                                     visitor, compare, combine,
                                                     zero, infinity)
    except StopSearch:
        pass

    return minimized, dist_map, pred_map


class AStarVisitor(object):
    r"""A visitor object that is invoked at the event-points inside the
    :func:`~graph_tool.search.astar_search` algorithm. By default, it
    performs no action, and should be used as a base class in order to be
    useful.
    """

    def initialize_vertex(self, u):
        """
        This is invoked on every vertex of the graph before the start of the
        graph search.
        """
        return

    def examine_vertex(self, u):
        """ This is invoked on a vertex as it is popped from the queue (i.e. it
        has the lowest cost on the ``OPEN`` list). This happens immediately
        before examine_edge() is invoked on each of the out-edges of vertex u.
        """
        return

    def examine_edge(self, e):
        """
        This is invoked on every out-edge of each vertex after it is examined.
        """
        return

    def discover_vertex(self, u):
        """
        This is invoked when a vertex is first discovered and is added to the
        ``OPEN`` list.
        """
        return

    def edge_relaxed(self, e):
        """
        Upon examination, if the following condition holds then the edge is
        relaxed (its distance is reduced), and this method is invoked.

        ::

            (u, v) = tuple(e)
            assert(compare(combine(d[u], weight[e]), d[v]))

        """
        return

    def edge_not_relaxed(self, e):
        """
        Upon examination, if the edge is not relaxed (see
        :meth:`~graph_tool.search.AStarVisitor.edge_relaxed`) then this
        method is invoked.
        """
        return

    def black_target(self, e):
        """ This is invoked when a vertex that is on the ``CLOSED`` list is
        "rediscovered" via a more efficient path, and is re-added to the
        ``OPEN`` list.
        """
        return

    def finish_vertex(self, u):
        """
        This is invoked on a vertex when it is added to the CLOSED list,
        which happens after all of its out edges have been examined.
        """
        return


def astar_search(g, source, weight, visitor=AStarVisitor(),
                 heuristic=lambda v: 1, dist_map=None, pred_map=None,
                 cost_map=None, combine=lambda a, b: a + b,
                 compare=lambda a, b: a < b, zero=0,
                 infinity=float('inf'), implicit=False):
    r"""
1203
    Heuristic :math:`A^*` search on a weighted, directed or undirected graph for the case where all edge weights are non-negative.
Tiago Peixoto's avatar
Tiago Peixoto committed
1204
1205
1206
1207
1208
1209
1210
1211
1212

    Parameters
    ----------
    g : :class:`~graph_tool.Graph`
        Graph to be used.
    source : :class:`~graph_tool.Vertex`
        Source vertex.
    weight : :class:`~graph_tool.PropertyMap`
        Edge property map with weight values.
Tiago Peixoto's avatar
Tiago Peixoto committed
1213
    visitor : :class:`~graph_tool.search.AStarVisitor` (optional, default: ``AStarVisitor()``)
Tiago Peixoto's avatar
Tiago Peixoto committed
1214
1215
1216
1217
1218
1219
1220
        A visitor object that is invoked at the event points inside the
        algorithm. This should be a subclass of
        :class:`~graph_tool.search.AStarVisitor`.
    heuristic : unary function (optional, default: ``lambda v: 1``)
        The heuristic function that guides the search. It should take a single
        argument which is a :class:`~graph_tool.Vertex`, and output an estimated
        distance from the source vertex.
Tiago Peixoto's avatar
Tiago Peixoto committed
1221
    dist_map : :class:`~graph_tool.PropertyMap` (optional, default: ``None``)
Tiago Peixoto's avatar
Tiago Peixoto committed
1222
1223
        A vertex property map where the distances from the source will be
        stored.
Tiago Peixoto's avatar
Tiago Peixoto committed
1224
    pred_map : :class:`~graph_tool.PropertyMap` (optional, default: ``None``)
Tiago Peixoto's avatar
Tiago Peixoto committed
1225
1226
        A vertex property map where the predecessor map will be
        stored (must have value type "int").
Tiago Peixoto's avatar
Tiago Peixoto committed
1227
    cost_map : :class:`~graph_tool.PropertyMap` (optional, default: ``None``)
Tiago Peixoto's avatar
Tiago Peixoto committed
1228
1229
1230
1231
1232
1233
1234
1235
1236
        A vertex property map where the vertex costs will be stored. It must
        have the same value type as ``dist_map``. This parameter is only used if
        ``implicit`` is True.
    combine : binary function (optional, default: ``lambda a, b: a + b``)
        This function is used to combine distances to compute the distance of a
        path.
    compare : binary function (optional, default: ``lambda a, b: a < b``)
        This function is use to compare distances to determine which vertex is
        closer to the source vertex.
Tiago Peixoto's avatar
Tiago Peixoto committed
1237
    implicit : bool (optional, default: ``False``)
Tiago Peixoto's avatar
Tiago Peixoto committed
1238
1239
        If true, the underlying graph will be assumed to be implicit
        (i.e. constructed during the search).
1240
1241
1242
1243
1244
1245
    zero : int or float (optional, default: ``0``)
         Value assumed to correspond to a distance of zero by the combine and
         compare functions.
    infinity : int or float (optional, default: ``float('inf')``)
         Value assumed to correspond to a distance of infinity by the combine and
         compare functions.
Tiago Peixoto's avatar
Tiago Peixoto committed
1246
1247
1248
1249
1250
1251
1252
1253
1254
1255
1256
1257
1258
1259
1260
1261
1262
1263
1264
1265
1266
1267
1268
1269
1270
1271
1272
1273
1274
1275
1276
1277
1278
1279
1280
1281
1282
1283
1284
1285
1286
1287
1288
1289
1290
1291
1292
1293
1294
1295
1296
1297
1298
1299
1300
1301
1302
1303
1304
1305
1306
1307
1308
1309
1310
1311
1312
1313
1314
1315
1316
1317
1318
1319
1320
1321
1322
1323
1324
1325
1326
1327
1328
1329
1330
1331
1332
1333
1334
1335
1336
1337
1338
1339
1340
1341
1342
1343
1344
1345
1346
1347
1348
1349
1350
1351
1352
1353
1354
1355
1356
1357
1358
1359
1360
1361
1362
1363
1364
1365
1366
1367
1368
1369
1370
1371
1372
1373
1374
1375
1376
1377
1378
1379
1380
1381
1382
1383
1384
1385
1386
1387
1388
1389
1390
1391


    Returns
    -------
    dist_map : :class:`~graph_tool.PropertyMap`
        A vertex property map with the computed distances from the source.

    See Also
    --------
    bfs_search: Breadth-first search
    dfs_search: Depth-first search
    dijkstra_search: Dijkstra's search algorithm

    Notes
    -----

    The :math:`A^*` algorithm is a heuristic graph search algorithm: an
    :math:`A^*` search is "guided" by a heuristic function. A heuristic function
    :math:`h(v)` is one which estimates the cost from a non-goal state (v) in
    the graph to some goal state, t. Intuitively, :math:`A^*` follows paths
    (through the graph) to the goal that are estimated by the heuristic function
    to be the best paths. Unlike best-first search, :math:`A^*` takes into
    account the known cost from the start of the search to v; the paths
    :math:`A^*` takes are guided by a function :math:`f(v) = g(v) + h(v)`, where
    :math:`h(v)` is the heuristic function, and :math:`g(v)` (sometimes denoted
    :math:`c(s, v)`) is the known cost from the start to v. Clearly, the
    efficiency of :math:`A^*` is highly dependent on the heuristic function with
    which it is used.

    The time complexity is :math:`O((E + V) \log V)`.

    The pseudo-code for the :math:`A^*` algorithm is listed below, with the
    annotated event points, for which the given visitor object will be called
    with the appropriate method.

    ::

        A*(G, source, h)
          for each vertex u in V                 initialize vertex u
            d[u] := f[u] := infinity
            color[u] := WHITE
          end for
          color[s] := GRAY
          d[s] := 0
          f[s] := h(source)
          INSERT(Q, source)                      discover vertex source
          while (Q != Ø)
            u := EXTRACT-MIN(Q)                  examine vertex u
            for each vertex v in Adj[u]          examine edge (u,v)
              if (w(u,v) + d[u] < d[v])
                d[v] := w(u,v) + d[u]            edge (u,v) relaxed
                f[v] := d[v] + h(v)
                if (color[v] = WHITE)
                  color[v] := GRAY
                  INSERT(Q, v)                   discover vertex v
                else if (color[v] = BLACK)
                  color[v] := GRAY
                  INSERT(Q, v)                   reopen vertex v
                end if
              else
                ...                              edge (u,v) not relaxed
            end for
            color[u] := BLACK                    finish vertex u
          end while


    Examples
    --------

    We will use an irregular two-dimensional lattice as an example, where the
    heuristic function will be based on the euclidean distance to the target.

    The heuristic function will be defined as:

    .. testcode::

        def h(v, target, pos):
            return sqrt(sum((pos[v].a - pos[target].a) ** 2))

    where ``pos`` is the vertex position in the plane.

    We must define what should be done during the search by subclassing
    :class:`~graph_tool.search.AStarVisitor`, and specializing the appropriate
    methods. In the following we will keep track of the discovered vertices, and
    which edges were examined, as well as the predecessor tree. We will also
    abort the search when a given target vertex is found, by raising the
    :class:`~graph_tool.search.StopSearch` exception.

    .. testcode::

        class VisitorExample(gt.AStarVisitor):

            def __init__(self, touched_v, touched_e, target):
                self.touched_v = touched_v
                self.touched_e = touched_e
                self.target = target

            def discover_vertex(self, u):
                self.touched_v[u] = True

            def examine_edge(self, e):
                self.touched_e[e] = True

            def edge_relaxed(self, e):
                if e.target() == self.target:
                    raise gt.StopSearch()


    With the above class defined, we can perform the :math:`A^*` search as
    follows.

    .. testsetup::

        from numpy.random import seed, random
        import matplotlib.cm
        seed(42)

    >>> points = random((500, 2)) * 4
    >>> points[0] = [-0.01, 0.01]
    >>> points[1] = [4.01, 4.01]
    >>> g, pos = gt.triangulation(points, type="delaunay")
    >>> weight = g.new_edge_property("double") # Edge weights corresponding to
    ...                                        # Euclidean distances
    >>> for e in g.edges():
    ...    weight[e] = sqrt(sum((pos[e.source()].a -
    ...                          pos[e.target()].a) ** 2))
    >>> touch_v = g.new_vertex_property("bool")
    >>> touch_e = g.new_edge_property("bool")
    >>> target = g.vertex(1)
    >>> dist, pred = gt.astar_search(g, g.vertex(0), weight,
    ...                              VisitorExample(touch_v, touch_e, target),
    ...                              heuristic=lambda v: h(v, target, pos))

    We can now observe the best path found, and how many vertices and edges were
    visited in the process.

    >>> ecolor = g.new_edge_property("string")
    >>> ewidth = g.new_edge_property("double")
    >>> ewidth.a = 1
    >>> for e in g.edges():
    ...    ecolor[e] = "blue" if touch_e[e] else "black"
    >>> v = target
    >>> while v != g.vertex(0):
    ...     p = g.vertex(pred[v])
    ...     for e in v.out_edges():
    ...         if e.target() == p:
1392
    ...             ecolor[e] = "#a40000"
Tiago Peixoto's avatar
Tiago Peixoto committed
1393
1394
    ...             ewidth[e] = 3
    ...     v = p
1395
1396
1397
    >>> gt.graph_draw(g, pos=pos, output_size=(300, 300), vertex_fill_color=touch_v,
    ...               vcmap=matplotlib.cm.binary, edge_color=ecolor,
    ...               edge_pen_width=ewidth, output="astar-delaunay.pdf")
Tiago Peixoto's avatar
Tiago Peixoto committed
1398
1399
    <...>

Tiago Peixoto's avatar
Tiago Peixoto committed
1400
1401
1402
1403
1404
1405
1406
    .. testcode::
       :hide:

       gt.graph_draw(g, pos=pos, output_size=(300, 300), vertex_fill_color=touch_v,
                     vcmap=matplotlib.cm.binary, edge_color=ecolor,
                     edge_pen_width=ewidth, output="astar-delaunay.png")

1407
    .. figure:: astar-delaunay.*
Tiago Peixoto's avatar
Tiago Peixoto committed
1408
1409
1410
1411
1412
1413
1414
1415
1416
1417
1418
1419
1420
1421
1422
1423
1424
1425
1426
1427
1428
1429
1430
1431
1432
1433
1434
1435
1436
1437
1438
1439
1440
1441
1442
1443
1444
1445
1446
1447
1448
       :align: center

       The shortest path is shown in red. The visited edges are shown in blue,
       and the visited vertices in black.


    The :math:`A^*` algorithm is very useful for searching *implicit* graphs,
    i.e. graphs which are not entirely stored in memory and are generated
    "on-the-fly" during the search. In the following example we will carry out a
    search in a hamming hypercube of 10 bits witch has random weights on its
    edges in the range :math:`[0,1]`. The vertices of the hypercube will be
    created during the search.

    The heuristic function will use the Hamming distance between vertices:

    .. testcode::

        def h(v, target, state):
            return sum(abs(state[v].a - target)) / 2


    In the following visitor we will keep growing the graph on-the-fly, and
    abort the search when a given target vertex is found, by raising the
    :class:`~graph_tool.search.StopSearch` exception.

    .. testcode::

        from numpy.random import random

        class HammingVisitor(gt.AStarVisitor):

            def __init__(self, g, target, state, weight, dist, cost):
                self.g = g
                self.state = state
                self.target = target
                self.weight = weight
                self.dist = dist
                self.cost = cost
                self.visited = {}

            def examine_vertex(self, u):
1449
                for i in range(len(self.state[u])):
Tiago Peixoto's avatar
Tiago Peixoto committed
1450
1451
1452
1453
1454
1455
1456
1457
1458
1459
1460
1461
1462
1463
1464
1465
1466
1467
1468
1469
1470
1471
1472
1473
1474
                    nstate = list(self.state[u])
                    nstate[i] ^= 1
                    if tuple(nstate) in self.visited:
                        v = self.visited[tuple(nstate)]
                    else:
                        v = self.g.add_vertex()
                        self.visited[tuple(nstate)] = v
                        self.state[v] = nstate
                        self.dist[v] = self.cost[v] = float('inf')
                    for e in u.out_edges():
                        if e.target() == v:
                            break
                    else:
                        e = self.g.add_edge(u, v)
                        self.weight[e] = random()
                self.visited[tuple(self.state[u])] = u

            def edge_relaxed(self, e):
                if self.state[e.target()] == self.target:
                    self.visited[tuple(self.target)] = e.target()
                    raise gt.StopSearch()

    With the above class defined, we can perform the :math:`A^*` search as
    follows.

1475
1476
1477
1478
1479
1480
    .. testsetup::

        from numpy.random import seed, random
        seed(42)
        gt.seed_rng(42)

Tiago Peixoto's avatar
Tiago Peixoto committed
1481
1482
1483
1484
1485
1486
1487
1488
1489
1490
1491
1492
1493
1494
1495
1496
1497
1498
1499
    >>> g = gt.Graph(directed=False)
    >>> state = g.new_vertex_property("vector<bool>")
    >>> v = g.add_vertex()
    >>> state[v] = [0] * 10
    >>> target = [1] * 10
    >>> weight = g.new_edge_property("double")
    >>> dist = g.new_vertex_property("double")
    >>> cost = g.new_vertex_property("double")
    >>> visitor = HammingVisitor(g, target, state, weight, dist, cost)
    >>> dist, pred = gt.astar_search(g, g.vertex(0), weight, visitor, dist_map=dist,
    ...                              cost_map=cost, heuristic=lambda v: h(v, array(target), state),
    ...                              implicit=True)

    We can now observe the best path found, and how many vertices and edges were
    visited in the process.

    >>> ecolor = g.new_edge_property("string")
    >>> vcolor = g.new_vertex_property("string")
    >>> ewidth = g.new_edge_property("double")
1500
    >>> ewidth.a = 1
Tiago Peixoto's avatar
Tiago Peixoto committed
1501
1502
1503
1504
1505
1506
1507
1508
1509
1510
    >>> for e in g.edges():
    ...    ecolor[e] = "black"
    >>> for v in g.vertices():
    ...    vcolor[v] = "white"
    >>> v = visitor.visited[tuple(target)]
    >>> while v != g.vertex(0):
    ...     vcolor[v] = "black"
    ...     p = g.vertex(pred[v])
    ...     for e in v.out_edges():
    ...         if e.target() == p:
1511
1512
    ...             ecolor[e] = "#a40000"
    ...             ewidth[e] = 3
Tiago Peixoto's avatar
Tiago Peixoto committed
1513
1514
    ...     v = p
    >>> vcolor[v] = "black"
1515
1516
1517
1518
1519
1520
1521
1522
1523
1524
    >>> pos = gt.graph_draw(g, output_size=(300, 300), vertex_fill_color=vcolor, edge_color=ecolor,
    ...                     edge_pen_width=ewidth, output="astar-implicit.pdf")

    .. testcode::
       :hide:

       gt.graph_draw(g, pos=pos, output_size=(300, 300), vertex_fill_color=vcolor,
                     edge_color=ecolor, edge_pen_width=ewidth,
                     output="astar-implicit.png")

Tiago Peixoto's avatar
Tiago Peixoto committed
1525

1526
    .. figure:: astar-implicit.*
Tiago Peixoto's avatar
Tiago Peixoto committed
1527
1528
1529
1530
1531
1532
1533
1534
1535
1536
1537
       :align: center

       The shortest path is shown in red, and the vertices which belong to it
       are in black. Note that the number of vertices visited is much smaller
       than the total number :math:`2^{10} = 1024`.

    References
    ----------

    .. [astar] Hart, P. E.; Nilsson, N. J.; Raphael, B. "A Formal Basis for the
       Heuristic Determination of Minimum Cost Paths". IEEE Transactions on
1538
       Systems Science and Cybernetics SSC4 4 (2): 100-107, 1968.
Tiago Peixoto's avatar
Tiago Peixoto committed
1539
1540
1541
1542
1543
       :doi:`10.1109/TSSC.1968.300136`
    .. [astar-bgl] http://www.boost.org/doc/libs/release/libs/graph/doc/astar_search.html
    .. [astar-wikipedia] http://en.wikipedia.org/wiki/A*_search_algorithm
    """

1544
1545
1546
1547
    visitor = VisitorWrapper(g, visitor,
                             ["initialize_vertex", "examine_vertex", "finish_vertex"],
                             ["initialize_vertex"])

Tiago Peixoto's avatar
Tiago Peixoto committed
1548
1549
1550
1551
1552
1553
1554
1555
1556
1557
1558
1559
    if dist_map is None:
        dist_map = g.new_vertex_property(weight.value_type())
    if pred_map is None:
        pred_map = g.new_vertex_property("int")
    if pred_map.value_type() != "int32_t":
        raise ValueError("pred_map must be of value type 'int32_t', not '%s'." % \
                             pred_map.value_type())

    dist_type = dist_map.python_value_type()
    if dist_type is not object:
        h = lambda v: dist_type(heuristic(v))
    else:
1560
        h = heuristic
Tiago Peixoto's avatar
Tiago Peixoto committed
1561

1562
    try:
1563
1564
        if dist_map.value_type() != "python::object":
            zero = _python_type(dist_map.value_type())(zero)
1565
1566
1567
1568
1569
    except OverflowError:
        zero = (weight.a.max() + 1) * g.num_vertices()
        zero = _python_type(dist_map.value_type())(zero)

    try:
1570
1571
        if dist_map.value_type() != "python::object":
            infinity = _python_type(dist_map.value_type())(infinity)
1572
1573
1574
1575
1576
    except OverflowError:
        infinity = (weight.a.max() + 1) * g.num_vertices()
        infinity = _python_type(dist_map.value_type())(infinity)


Tiago Peixoto's avatar
Tiago Peixoto committed
1577
1578
1579
1580
1581
1582
    try:
        if not implicit:
            g._Graph__perms.update({"del_vertex": False, "del_edge": False,
                                    "add_edge": False})

            libgraph_tool_search.astar_search(g._Graph__graph,
1583
                                              weakref.ref(g),
Tiago Peixoto's avatar
Tiago Peixoto committed
1584
1585
1586
1587
1588
1589
1590
1591
1592
1593
1594
1595
1596
                                              int(source), _prop("v", g, dist_map),
                                              _prop("v", g, pred_map),
                                              _prop("e", g, weight), visitor,
                                              compare, combine, zero, infinity,
                                              h)
        else:
            if cost_map is None:
                cost_map = g.new_vertex_property(dist_map.value_type())
            elif cost_map.value_type() != dist_map.value_type():
                raise ValueError("The cost_map value type must be the same as" +
                                 " dist_map.")
            g._Graph__perms.update({"del_vertex": False})
            libgraph_tool_search.astar_search_implicit\
1597
                (g._Graph__graph, weakref.ref(g), int(source),
Tiago Peixoto's avatar
Tiago Peixoto committed
1598
1599
1600
1601
1602
1603
1604
1605
1606
1607
1608
1609
1610
1611
                 _prop("v", g, dist_map), _prop("v", g, pred_map),
                 _prop("v", g, cost_map), _prop("e", g, weight), visitor,
                 compare, combine, zero, infinity, h)
    except StopSearch:
        g._Graph__perms.update({"del_vertex": True, "del_edge": True,
                                "add_edge": True})
    finally:
        g._Graph__perms.update({"del_vertex": True, "del_edge": True,
                                "add_edge": True})

    return dist_map, pred_map


class StopSearch(Exception):
Tiago Peixoto's avatar
Tiago Peixoto committed
1612
    """If this exception is raised from inside any search visitor object, the search is aborted."""
Tiago Peixoto's avatar
Tiago Peixoto committed
1613
    pass