<table><thead><tr><th>N/1000</th><th>insertion</th><th>quicksort</th><th>quicksort (log stack)</th><th>quicksort (rand)</th><th>mergesort</th><th>mergesort (reuse)</th><th>List.sort</th><th>heapsort (arrays)</th><th>heapsort (trees)</th></tr></thead><tbody><tr><td>2<sup>4</sup></td><td>14.94</td><td>-</td><td>8.94</td><td>0.02</td><td>0.02</td><td><strong>0.01</strong></td><td><strong>0.01</strong></td><td>0.03</td><td>0.03</td></tr><tr><td>2<sup>5</sup></td><td>-</td><td>-</td><td>-</td><td>0.06</td><td>0.05</td><td><strong>0.02</strong></td><td><strong>0.02</strong></td><td>0.06</td><td>0.07</td></tr><tr><td>2<sup>6</sup></td><td>-</td><td>-</td><td>-</td><td>0.13</td><td>0.11</td><td>0.05</td><td><strong>0.04</strong></td><td>0.12</td><td>0.19</td></tr><tr><td>2<sup>7</sup></td><td>-</td><td>-</td><td>-</td><td>0.26</td><td>0.27</td><td><strong>0.11</strong></td><td><strong>0.11</strong></td><td>0.25</td><td>0.49</td></tr><tr><td>2<sup>8</sup></td><td>-</td><td>-</td><td>-</td><td>0.77</td><td>0.62</td><td><strong>0.25</strong></td><td>0.27</td><td>0.54</td><td>1.18</td></tr><tr><td>2<sup>9</sup></td><td>-</td><td>-</td><td>-</td><td>1.66</td><td>1.38</td><td><strong>0.58</strong></td><td>0.59</td><td>1.16</td><td>2.76</td></tr><tr><td>2<sup>10</sup></td><td>-</td><td>-</td><td>-</td><td>-</td><td>3.07</td><td>1.28</td><td><strong>1.26</strong></td><td>2.43</td><td>6.36</td></tr></tbody></table>

<figure><img src="image_125.png" /><figcaption>Figure 12.9: Comparison of sorting algorithms on already sorted lists.</figcaption></figure>

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
<table><caption>Figure 12.10: Sorting algorithms on arrays of integers.</caption><tr><td>insertion</td><td>insertion sort of program 90</td></tr><tr><td>quicksort</td><td>quicksort of program 92</td></tr><tr><td>quicksort (insertion)</td><td>quicksort using an insertion sort on small segments (see exercise 12.7)</td></tr><tr><td>quicksort (log stack)</td><td>quicksort with the first call on the smallest of the two segments (see exercise 12.6)</td></tr><tr><td>quicksort (rand)</td><td>quicksort with a pivot chosen randomly (see exercise 12.8)</td></tr><tr><td>mergesort</td><td>merge sort of program 94</td></tr><tr><td>merge sort (swap)</td><td>merge sort with the optimization consisting in avoiding the use of Array.blit (see exercise 12.12)</td></tr><tr><td>Array.stable_sort</td><td>the function Array.stable_sort of the OCaml standard library (version 3.11.2)</td></tr><tr><td>heapsort</td><td>heapsort of program 97</td></tr><tr><td>Array.sort</td><td>heapsort of the OCaml standard library</td></tr></table>

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
<table><thead><tr><th>N/1000</th><th>insertion</th><th>quicksort</th><th>quicksort (insertion)</th><th>quicksort (log stack)</th><th>quicksort (rand)</th><th>mergesort</th><th>mergesort (swap)</th><th>Array.stable_sort</th><th>heapsort</th><th>Array.sort</th></tr></thead><tbody><tr><td>2<sup>4</sup></td><td>1.65</td><td>0.01</td><td>0.01</td><td>0.01</td><td>0.01</td><td>0.01</td><td>0.01</td><td>0.01</td><td>0.01</td><td>0.01</td></tr><tr><td>2<sup>5</sup></td><td>-</td><td>0.02</td><td>0.02</td><td>0.02</td><td>0.02</td><td>0.02</td><td>0.02</td><td><strong>0.01</strong></td><td>0.02</td><td>0.02</td></tr><tr><td>2<sup>6</sup></td><td>-</td><td>0.04</td><td>0.04</td><td><strong>0.03</strong></td><td>0.04</td><td>0.04</td><td><strong>0.03</strong></td><td><strong>0.03</strong></td><td>0.04</td><td>0.05</td></tr><tr><td>2<sup>7</sup></td><td>-</td><td>0.08</td><td>0.08</td><td>0.08</td><td>0.08</td><td>0.09</td><td>0.07</td><td><strong>0.06</strong></td><td>0.10</td><td>0.10</td></tr><tr><td>2<sup>8</sup></td><td>-</td><td>0.16</td><td>0.15</td><td>0.16</td><td>0.18</td><td>0.19</td><td>0.14</td><td><strong>0.13</strong></td><td>0.21</td><td>0.21</td></tr><tr><td>2<sup>9</sup></td><td>-</td><td>0.33</td><td>0.33</td><td>0.33</td><td>0.36</td><td>0.41</td><td>0.29</td><td><strong>0.26</strong></td><td>0.46</td><td>0.46</td></tr><tr><td>2<sup>10</sup></td><td>-</td><td>0.74</td><td>0.70</td><td>0.71</td><td>0.72</td><td>0.84</td><td>0.61</td><td><strong>0.57</strong></td><td>1.01</td><td>1.02</td></tr><tr><td>2<sup>11</sup></td><td>-</td><td>1.43</td><td>1.40</td><td>1.41</td><td>1.64</td><td>1.76</td><td>1.28</td><td><strong>1.10</strong></td><td>2.36</td><td>2.35</td></tr><tr><td>2<sup>12</sup></td><td>-</td><td>3.05</td><td>3.01</td><td>3.05</td><td>3.41</td><td>3.74</td><td>2.69</td><td><strong>2.34</strong></td><td>5.37</td><td>5.40</td></tr></tbody></table>

<figure><img src="image_126.png" /><figcaption>Figure 12.11: Comparison of sorting on random arrays.
— openly licensed via CC BY SA 4.0 —</figcaption></figure>
<table><thead><tr><th>N/1000</th><th>insertion</th><th>quicksort</th><th>quicksort (insertion)</th><th>quicksort (log stack)</th><th>quicksort (rand)</th><th>mergesort</th><th>mergesort (swap)</th><th>Array.stable_sort</th><th>heapsort</th><th>Array.sort</th></tr></thead><tbody><tr><td>2<sup>4</sup></td><td><strong>0.00</strong></td><td>1.98</td><td>2.06</td><td>2.07</td><td>0.01</td><td>0.01</td><td>0.01</td><td><strong>0.00</strong></td><td>0.01</td><td>0.01</td></tr><tr><td>2<sup>5</sup></td><td><strong>0.00</strong></td><td>-</td><td>-</td><td>-</td><td>0.02</td><td>0.02</td><td><strong>0.01</strong></td><td><strong>0.01</strong></td><td>0.02</td><td>0.02</td></tr><tr><td>2<sup>6</sup></td><td><strong>0.00</strong></td><td>-</td><td>-</td><td>-</td><td>0.04</td><td>0.04</td><td><strong>0.02</strong></td><td><strong>0.02</strong></td><td>0.04</td><td>0.04</td></tr><tr><td>2<sup>7</sup></td><td><strong>0.00</strong></td><td>-</td><td>-</td><td>-</td><td>0.07</td><td>0.07</td><td>0.05</td><td><strong>0.04</strong></td><td>0.10</td><td>0.08</td></tr><tr><td>2<sup>8</sup></td><td><strong>0.00</strong></td><td>-</td><td>-</td><td>-</td><td>0.16</td><td>0.15</td><td>0.10</td><td><strong>0.08</strong></td><td>0.21</td><td>0.16</td></tr><tr><td>2<sup>9</sup></td><td><strong>0.01</strong></td><td>-</td><td>-</td><td>-</td><td>0.32</td><td>0.33</td><td>0.21</td><td><strong>0.18</strong></td><td>0.45</td><td>0.37</td></tr><tr><td>2<sup>10</sup></td><td><strong>0.02</strong></td><td>-</td><td>-</td><td>-</td><td>0.71</td><td>0.67</td><td>0.43</td><td><strong>0.37</strong></td><td>0.94</td><td>0.73</td></tr><tr><td>2<sup>11</sup></td><td><strong>0.05</strong></td><td>-</td><td>-</td><td>-</td><td>1.47</td><td>1.38</td><td>0.90</td><td><strong>0.76</strong></td><td>1.98</td><td>1.53</td></tr><tr><td>2<sup>12</sup></td><td><strong>0.11</strong></td><td>-</td><td>-</td><td>-</td><td>3.11</td><td>2.98</td><td>1.94</td><td><strong>1.65</strong></td><td>4.28</td><td>3.34</td></tr></tbody></table>

<figure><img src="image_127.png" /><figcaption>Figure 12.12: Comparison of sorting on already sorted arrays.
— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —</figcaption></figure>
— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —
13

Graph Algorithms

In this chapter, we present graph algorithms. These algorithms are independent
of the underlying graph data structure. They can be used, in particular, with all
the graph structures of chapter 7. Below, we will take as given a module *G* whose
signature includes at least the types, functions, and modules of program 98.
Graphs are of type *graph*, and their vertices are of type *vertex*. The function
*iter saf e* iter saf e iter saf e iter saf e iter saf e iter saf e iter saf e
iter saf e iterates over the vertices of a given graph, applying a function
passed as argument to each of them. The function iteruites over the
successors of a given vertex. The module H provides a hash table structure with
vertices as keys.

**13.1 Breadth-First Search**

The first traversal that we present is called *breadth-first search* (BFS). It consists in exploring the graph by moving outwards “in concentric circles,” from a particular vertex *s*, called the root. We first traverse the vertices that are exactly one arc away from the root *s*, then those situated two arcs away, and so forth. Consider, for example, the graph on the left in figure 13.1. Beginning from vertex 5, we first visit vertices 4 and 1, which are directly linked to vertex 5.

— openly licensed via CC BY SA 4.0 —
<caption>**Program 98** — Minimal signature of a graph structure.</caption>

```h
    type graph
    type vertex
    val iter vertex: (vertex -> unit) -> graph -> unit
    val iter succ: (vertex -> unit) -> graph -> vertex -> unit
module H: HashTable with type key = vertex
```

We then visit vertices 0, 6 and 2, followed by vertex 7, and finally vertex 3. The graph, with distances to the root indicated as exponents, is shown on the right.

<figure><img src="image_128.png" /><figcaption>Figure 13.1: Breadth-first search of a graph.</figcaption></figure>

We write breadth-first search as a function `iter_bfs` that takes as arguments a graph `g`, a root `s`, and a function `f` that will be applied to the elements of `g` as these are discovered. We begin by creating a hash table that will contain the vertices already discovered by the traversal.

```h
let iter_bfs f g s =
    let visited = H.create () in
```

The search itself is based on the use of a *queue* into which vertices are inserted as they are discovered. The idea is that, at any point in time, the queue will contain the vertices situated at distance `d` from the root, followed by the vertices situated at distance `d + 1`:

```
  vertices at distance d
  vertices at distance d + 1
```

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0/) —
This is the concretization of our notion of “concentric circles,” more precisely of the two consecutive concentric circles that are being considered. This property is crucial for the correctness of breadth-first search<fnref n="1" />. Here, we use the module `Queue` to implement the queue. Initially, it contains only the root. It is useful to write a function `add` that adds a vertex simultaneously in the hash table and the queue. We begin by using it on the root `s`.

```python
let queue = Queue.create() in
let add w = H.add visited w (); Queue.add w queue in
add s;
```

Another important property is that each vertex contained in the queue is also contained in the table `visited`, which the function `add` guarantees in an obvious manner. We use a loop to examine vertices in the order in which they are discovered, until all vertices reachable from the root have been examined. As long as the queue is not empty, we extract the first element `v` of the queue and apply the visitor function `f` to it.

```python
while not (Queue.is_empty queue) do
    let v = Queue.pop queue in
    f v;
```

For each successor `w` of `v` not yet found, that is, not contained in the table `visited`, we use `add` to add it to both the queue and the table.

```python
iter (fun w -> if not (H.mem visited w) then add w) g v
done
```

This completes the `while` loop and the function `iter_bfs`. The complete code is given in program 99 below.

This code applies equally well to undirected and directed graphs. Note that some vertices may not be discovered by this algorithm. These are the vertices `v` for which there is no path between the root and `v`. In other words, breadth-first search determines the set of vertices that can be reached from the root and even determines their minimal distance from the root in number of arcs (see exercise 13.1).

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —

<footnote n="1">¹For a detailed proof of the correctness of breadth-first search, consult [7, chap. 23].</footnote>
<caption>Program 99 — Breadth-first search.</caption>

```pascal
let iter_bfs f g s =
    let visited = H.create () in
    let queue = Queue.create () in
    let add w = H.add visited w (); Queue.add w queue in
    add s;
    while not (Queue.is_empty queue) do
        let v = Queue.pop queue in
        f v;
        iter (fun w -> if not (H.mem visited w) then add w) g v
done
```

The complexity of breadth-first search is easy to determine. Each vertex is placed in the queue at most once and is therefore examined at most once. Each arc is considered at most once, when its starting vertex is examined. For a graph with $N$ vertices and $E$ arcs, the complexity is therefore $O(N + E)$, which is optimal. The complexity in space is $O(N)$, since the queue, like the hash table, may contain (almost) all vertices in the worst case.

There is another way to implement breadth-first search, without using a queue. As noted above, the queue has a very particular structure, with vertices at distance $d$ followed by vertices at distance $d+1$. We see, therefore, that the queue structure is not really necessary. All we need are two “buckets,” one containing the vertices at distance $d$, and the other, those at distance $d+1$. We may concretize these, for instance, using lists. When we are done with bucket $d$, we swap the two buckets. This does not affect the complexity.

## 13.2 Depth-First Search

The second graph traversal that we present is *depth-first search* (DFS). It uses a backtracking algorithm: As long as we can follow an arc, we do so, and when it is no longer possible, we backtrack. As with breadth-first search, we mark vertices

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
as they are discovered, to avoid getting stuck in a cycle. Consider the example of the graph on the left in figure 13.2. We begin a depth-first search from vertex 2. Two arcs come out of this vertex, 2 → 4 and 2 → 5. We arbitrarily choose to consider the arc 2 → 5 first and, accordingly, proceed to vertex 5. No arc comes out of 5, so we cannot continue. We therefore backtrack to vertex 2, and consider the second outgoing arc, 2 → 4. From 4, we can only follow the arc 4 → 3. Similarly, from 3, we can only follow the arc 3 → 1. Two arcs emerge from vertex 1: 1 → 0 and 1 → 4. We choose to follow the arc 1 → 4 first. It leads to a vertex that has already been visited. We therefore backtrack. Returning to 1, we consider the other arc, 1 → 0, which leads to 0. From there, the only outgoing arc leads to 3, which has also already been visited. We therefore return to 0, then to 1, then to 3, then to 4, and finally to 2. The traversal is complete. If we redraw the graph with the order in which the vertices were discovered as exponents, we get the graph on the right in figure 13.2.

<figure><img src="image_129.png" /><figcaption>Figure 13.2: Example of depth-first search.</figcaption></figure>

As with breadth-first search, we will make use of a hash table containing the vertices already discovered.

```python
let iter_dfs f g s =
    let visited = H.create () in
```

Depth-first search is written as a recursive local function `visit` taking a vertex `v` as argument. Its code is a few lines long:

```python
let rec visit v =
    if not (H.mem visited v) then begin
        H.add visited v ();
    f v;
```

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
Program 100 — Depth-first search.

```pascal
let iter_dfs f g s =
  let visited = H.create () in
  let rec visit v =
    if not (H.mem visited v) then begin
      H.add visited v ();
      f v;
      iter visit g v
    end
  in
  visit s
```

```pascal
  iterapi visit g v
end
```

If vertex v has already been discovered, there is nothing more to be done. Otherwise, we mark it and visit it with the function f. Then, we consider each successor w, recursively triggering a depth-first search on each one. One detail is crucial: We added v to the table *visited before* considering its successors. This is what keeps us from getting stuck in a cycle. To conclude the function *iter_dfs*, we have only to call *visit* on the root s:

```pascal
in
  visit s
```

The complete code is given in program 100.

The complexity is $O(N + E)$, by the same logic as in case of breadth-first search. Space complexity is subtler since it is the stack of recursive calls that contains the vertices being visited, and that plays the role of the queue in breadth-first search. Since in the worst case all the vertices can be present in the stack or hash table, it follows that the space complexity is $O(N)$.

Like breadth-first search, depth-first search allows us to determine the set of vertices that can be reached from the root s. Figure 13.3 illustrates another

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —
example, where depth-first search is triggered from vertex 1. Vertices 2 and 5 are not discovered. In particular, depth-first search, like breadth-first search, is a way to determine the existence of a path between a particular vertex, namely, the root, and other vertices of the graph. If that is the only objective (if, for example, the minimal distance does not interest us), depth-first search is generally more efficient. Indeed, its memory use (that is, the call stack) is often inferior to that of breadth-first search. The typical example is that of a tree, in which memory use is limited to the height of the tree for depth-first search, but could be as large as the entire tree in case of breadth-first search. Depth-first search has many other applications that lie beyond the scope of this book. See, for example, *Introduction to Algorithms* [7].

<figure><img src="image_130.png" /><figcaption>Figure 13.3: Another depth-first search.</figcaption></figure>

## 13.3 Shortest Path

We consider here graphs whose arcs are labeled with weights (representing distances, for example). We seek the shortest path between two vertices, where the length is no longer the number of arcs but rather the sum of weights along the path. Thus, if we consider the weighted graph of figure 13.4, the shortest path from vertex 2 to vertex 0 is of length 5. This is the path $2 \to 4 \to 3 \to 1 \to 0$, which is shorter than the path $2 \to 1 \to 0$, of length 6, even if the latter contains fewer arcs. Similarly, the shortest path from vertex 2 to vertex 5, via vertex 4, is of length 2.

As with breadth- and depth-first search, we take as given a vertex *s*, from which we seek to determine the shortest paths to the other vertices. We also require the weight of each arc. One possibility is to modify the graph structure

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —
<figure><img src="image_131.png" /><figcaption>Figure 13.4: Example of a weighted graph.</figcaption></figure>

to include this. Here, however, we take as given a weight function `weight` of the
following type, where weights are of type `float`:

```hkl
val weight: graph -> vertex -> vertex -> float
```

This function will only be called with arguments `g`, `x`, and `y` if there is an arc between `x` and `y` in the graph `g`.

We present two algorithms to solve the problem. The first, Dijkstra's algorithm, assumes that the weights are always non-negative. The second, the Bellman-Ford algorithm, involves no such hypothesis regarding the weights.

## Dijkstra's Algorithm

Dijkstra's algorithm is a generalization of breadth-first search. Here, too, we proceed in "concentric circles." The difference is that the radii of these circles represent distance in terms of total weight rather than number of arcs. Thus, in the example above, starting from the root 2, we first reach the vertices at distance 1 (namely, 4), then those at distance 2 (namely, 3 and 5), then those at distance 3 (namely, 1), and finally those at distance 5 (namely, 0). A possible complication is that we may reach a vertex at a certain distance (for example, vertex 5 with the arc $2 \to 5$), and then later find a shorter path by taking different arcs (for example, $2 \to 4 \to 5$). A queue no longer suffices as it did in breadth-first search. Instead, we will use a *priority queue* (see chapter 6). The priority queue will contain the vertices already found, ordered by distance from the root. When a better path to a vertex is found, the vertex is reinserted in

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —
the queue with a higher priority, that is, a shorter distance<fnref n="2" />.

We now describe the code of Dijkstra's algorithm. To implement the priority queue, we assume a functor `PriorityQueue`, parametrized by a type `t` equipped with a total order. (Such functors were presented in chapter 6.) Here, the type `t` will be that of pairs `(v, d)`, where `v` is a vertex and `d` is its distance from the root. The order on these pairs is that obtained by comparing the distances. We introduce the following module:

```python
module VertexDistance = struct
    type t = vertex * float
    let compare (., d1) (., d2) = Stdlib comparisons d1 d2
end
```

It only remains to apply the functor `PriorityQueue` to the module `VertexDistance`, that is:

```python
module P = PriorityQueue VertexDistance
```

We assume here that the module `P` provides a function `create` to construct a new priority queue, a function `add` to add an element, and a function `extract_min` that removes and returns the smallest element for the order defined by the function `compare`.

We now arrive at the code of the algorithm itself. We write it as a function `dijkstra` that takes as arguments a graph `g`, a root `s`, and a function `f`. The function `f` will be applied to elements of `g` as their shortest distance to `s` is determined, the vertex being passed to `f` as first argument, and its distance as second argument. We begin by creating a hash table that will contain the vertices already reached by the traversal.

```python
let rec dijkstra f g s =
    let visited = H.create () in
```

We also create a hash table, `distance`, and the priority queue, `queue`. The hash table will contain the distance currently known for each vertex.

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —

<footnote n="2">²Another solution would be to use a priority queue structure where it is possible to *modify* the priority of an element already in the queue. Such structures exist but are complex to implement and, although asymptotically better, do not necessarily bring any practical gains. The solution presented here is a good compromise.</footnote>
```pascal
let distance = H.create () in
let queue = P.create () in
```

It is convenient to introduce a function `add` that will add a vertex to both these data structures at the same time. We begin by applying it to the root `s`.

```pascal
let add v d = H.replace distance v d; P.add queue (v, d) in
add s 0.;
```

We will see below why `H.replace` is preferred over `H.add`. We then continue with a loop as long as the queue is not empty.

```pascal
while not (P.is_empty queue) do
```

We extract the first element of the queue, `u`, together with the distance to the root, `du`. If the vertex `u` is already in `visited`, then we have already determined the distance of `s` to `u`, and there is nothing more to be done.

```pascal
let (u, du) = P.extract_min queue in
if not (H.mem visited u) then begin
```

This situation can indeed occur when a first path is found and another, shorter, path is found subsequently. The latter then passes before the former in the priority queue. When the longer path finally comes out of the queue, it should be ignored. If the vertex is not in `visited`, however, that means that we have just found a shortest path from `s` to `u`, which we indicate by adding `u` to the table `visited` and by applying the visitor function `f` to `u` and `du`.

```pascal
H.add visited u ();
f u du;
```

Next, we will examine all arcs emerging from `u`. For each arc `u → v`, the distance to `v` obtained by following the corresponding arc is the sum of the distance `du` and the arc weight, that is, `weight u v`. Several scenarios are possible for the vertex `v`. Either it is the first time that it is reached, or it was already in `distance`. In the latter case, we may have decreased the distance to `v` by passing through `u`. If so, it suffices to use the function `add` to update the tables. We handle each successor `v` in the function `visit` below:

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0/) —
```pascal
let visit v =
    let d = du +. weight g u v in
    if not (H.mem distance v) || d < H.find distance v then
        add v d
```

It only remains to apply this function to all arcs emerging from u, which concludes the body of the `while` loop.

```pascal
in
iter visit g u
end
done
```

Once we exit the `while` loop, all vertices that may be reached from s have been visited. The complete code is given in program 101 above.

Figure 13.5 shows the result of Dijkstra's algorithm on the graph given in the previous example, with 2 as root. The graph is drawn with the final distances obtained for each vertex as exponents. In this example, all vertices have been reached by the traversal. As for the breadth- and depth-first searches, this is not always the case: only vertices for which there is a path from the root will be reached.

## Complexity

Let us evaluate the cost of Dijkstra's algorithm in the worst case. The priority queue can contain up to *E* elements, because the algorithm visits each arc at most once, and each arc can lead to the insertion of an element in the queue. Assuming that the priority queue operations `add` and `extract_min` have logarithmic cost (as is the case for the priority queues described in chapter 6), each operation on the queue will have cost *O*(log *E*), that is, *O*(log *N*), since *E* ≤ *N*². The total cost is therefore *O*(*E* log *N*).

## Bellman-Ford Algorithm

We present here a second solution to the shortest path problem, this time without assuming that the weights are non-negative. Omitting this assumption

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —
Program 101 — Shortest path (Dijkstra's algorithm).

```python
module VertexDistance = struct
    type t = vertex * float
    let compare (_, d1) (_, d2) = Stdlib compare d1 d2
end

module P = PriorityQueue VertexDistance

let rec dijkstra f g s =
    let visited = H.create () in
    let distance = H.create () in
    let queue = P.create () in
    let add v d = H.replace distance v d; P.add queue (v, d) in
    add s 0;
    while not (P.is_empty queue) do
        let (u, du) = P.extract_min queue in
        if not (H.mem visited u) then begin
            H.add visited u ();
            f u du;
            let visit v =
                let d = du +. weight g u v in
                if not (H.mem distance v) || d < H.find distance v then
                    add v d
                in
            iter visit g u
        end
    done
```

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
<figure><img src="image_132.png" /><figcaption>Figure 13.5: Result of Dijkstra's algorithm.</figcaption></figure>

introduces a subtlety. If the graph contains a cycle with strictly negative total weight, and if this cycle can be reached from the root, then none of the vertices that can be reached from the cycle will have a shortest path because we can take this cycle as many times as we wish, reducing the total weight of the path each time. The Bellman-Ford algorithm computes the shortest paths from the root and detects any negative cycle that may be reached from the root.

This algorithm consists in repeating the following operation $N$ times: Consider each arc $u \to v$ of the graph and determine if taking it will reduce the distance currently known between the source and $v$. The idea is that a shortest path between $s$ and $v$ cannot contain more than $N-1$ arcs, without which this path would contain a cycle, which we could remove from the path if it is of non-negative weight. This is why $N$ iterations suffice.

Let us now write the code of the Bellman-Ford algorithm. We begin by defining a function `iter_edge` that traverses all the arcs of a graph $g$.

```pascal
let iter_edge f g =
    iterVertex (fun u ->
                iterSucc (fun v -> f u v (weight g u v)) g u) g
```

For each arc $u \to v$, the function `iter_edge` applies the function $f$ to the vertices $u$ and $v$, and the distance of the arc $u \to v$ given by the function `weight`.

The function `bellman_ford` takes a graph $g$ and a root vertex $s$ as arguments. It begins by creating a hash table $h$ containing the known distances for each vertex.

```pascal
let bellman_ford g s =
    let h = H.create () in
```

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
We begin by initializing this table, giving $s$ distance 0 and every other vertex an infinite distance.

```pascal
iterVertex (fun v -> H.add h v maxFloat) g;
H.add h s 0.;
```

We use here the constant `maxFloat` to represent an infinite distance. Then we repeat the following operation $N-1$ times: We examine each arc, using a `for` loop and the function `iterEdge`.

```pascal
for i = 1 to nbVertex g - 1 do
  iterEdge (fun u v w ->
```

For each arc $u \rightarrow v$, we compute the distance $d$ from the root to $v$ that is obtained by taking this arc. If it is smaller than the best distance currently known, we record it.

```pascal
let d = H.find h u +. w in
if d < H.find h v then H.replace h v d
```

Once this `for` loop is complete, we redo the operation an $N$-th and final time. This time, the shortest paths can only be improved by the presence of a negative cycle. We will signal this case by raising an exception.

```pascal
done;
iterEdge (fun u v w ->
  if H.find h u +. w < H.find h v then
    raise NegativeCycle
```

This completes the Bellman-Ford algorithm. We can return the table $h$ that gives the distance of each vertex to the source. This is what is done in program 102.

We illustrate the Bellman-Ford algorithm with an example. Consider the graph of figure 13.6, for which we seek to determine the distances from vertex 5. The table to the right of the graph gives the distance found for each vertex after each step. Initially, all the vertices are at an infinite distance, except the root 5. At each step, we assume the vertices are traversed in the order of their numbering. We observe, for example, that the distance to vertex 2 initially

— openly licensed via CC BY SA 4.0 —
Program 102 — The Bellman-Ford algorithm.

```pascal
let iter_edge f g =
  iterVertex (fun u ->
    iterSucc (fun v -> f u v (weight g u v)) g u) g

exception NegativeCycle

let bellman_ford g s =
  let h = H.create () in
  iterVertex (fun v -> H.add h v max_float) g;
  H.add h s 0.;
  for i = 1 to nbVertex g - 1 do
    iterEdge (fun u v w ->
      let d = H.find h u +. w in
      if d < H.find h v then H.replace h v d
    ) g
  done;
  iterEdge (fun u v w ->
    if H.find h u +. w < H.find h v then
      raise NegativeCycle
  ) g;
  h
```

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
<figure><img src="image_133.png" /></figure>

<table><caption>Figure 13.6: Illustration of the Bellman-Ford algorithm.</caption><thead><tr><th>vertex</th><th colspan="5">distance after <i>i</i> steps</th></tr></thead><tbody><tr><td>0</td><td>&infin;</td><td>&infin;</td><td>2</td><td>-3</td><td>-3</td></tr><tr><td>1</td><td>&infin;</td><td>1</td><td>-3</td><td>-3</td><td>-3</td></tr><tr><td>2</td><td>&infin;</td><td>3</td><td>2</td><td>-2</td><td>-2</td></tr><tr><td>3</td><td>&infin;</td><td>&infin;</td><td>-2</td><td>-2</td><td>-2</td></tr><tr><td>4</td><td>&infin;</td><td>-4</td><td>-4</td><td>-4</td><td>-4</td></tr><tr><td>5</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td></tr></tbody></table>

takes the value 3 (the arc $5 \to 2$), then the value 2 (passing through vertex 1), and finally the value -2 (passing through the vertices 4 and 1). After that, it no longer varies. More generally, we observe here that none of the distances is improved during the last two steps of the algorithm. In particular, there is no negative cycle.

## Complexity

The cost of the Bellman-Ford algorithm in the worst case is easy to evaluate. We perform exactly *N* passes, and each pass examines each of the *E* arcs once and only once. The complexity is therefore $O(EN)$, that is, $O(N^3)$ in the worst case. In practice, however, we can significantly improve performance by observing that many arcs are examined unnecessarily when the distance from the root to the starting vertex has not been modified. Another simple optimization consists in stopping once a pass has not improved any of the distances. The complexity in the worst case remains, nevertheless, the same.

## 13.4 Minimum spanning tree

We consider here an undirected graph whose arcs are labeled with weights. We also assume that the graph is connected, that is, that each pair of vertices is linked by a path. A *spanning tree* is a subset of the arcs, that does not contain a cycle, and that links all vertices together. It is easy to convince ourselves that

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —
such a tree exists. For example, a depth-first search beginning at an arbitrary vertex will determine such a tree: it is formed by all the arcs that led to the discovery of a new vertex.

We seek then to solve the problem of finding a spanning tree whose total weight is minimal. We present here a solution known as *Kruskal's algorithm*. This algorithm consists in traversing the arcs in increasing order of weight. For each arc $u-v$, we determine if a cycle is created by adding the arc to the set of arcs already selected to form the spanning tree. If this is so, we ignore the arc. Otherwise, we select it. Once all the arcs have been examined, the set of selected arcs form a minimum spanning tree (see exercises 13.24 and 13.25).

To implement this algorithm, we need to find an efficient means to determine if selecting a new arc introduces a cycle among the already selected arcs. Luckily, there is a simple solution to this problem, in the form of the structure of disjoint sets studied in chapter 8. Indeed, it suffices to construct an initial structure in which each vertex constitutes a class by itself. For any selected arc $u-v$, we merge the classes of $u$ and $v$. This materializes the fact that there now exists a path between every pair of vertices in the new class. To test if an arc introduces a cycle, it suffices therefore to test whether its extremities belong to the same class.

Let us now write the code for Kruskal's algorithm, as a function `spanning_tree` that takes a graph as argument and returns a minimum spanning tree in the form of a list of arcs. We take as given a module `UF` for the structure of disjoint sets (see chapter 8), with the following signature:

```python
module UF: sig
    val create : vertex list -> t
    val find : t -> vertex -> vertex
    val union : t -> vertex -> vertex -> unit
end
```

We assume there is a function `vertices` that returns the list of vertices of a graph. Such a function is easy to obtain from `iterVertex`. Similarly, we assume there is a function `edges` that returns the list of all arcs of a graph, each arc being a triple $(x, y, w)$, where $w$ is the weight of the arc $x-y$. Here, too, this function is easy to obtain from `iterEdge` and `weight`.

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0/) —
The code of Kruskal's algorithm begins by creating a *union-find* structure based on the list of vertices of the graph:

```h
let spanning_tree g =
    let uf = UF.create (vertices g) in
```

Next, we sort the list of arcs of the graph in increasing order of weight, using the library function `List.sort`. The comparison used here is that of weights, with the usual order on floating-point numbers, that is, `Stdlib comparisons`.

```h
let compare (_,_,w1) (_,_,w2) = Stdlib comparisons w1 w2 in
let edges = List.sort compare (edges g) in
```

The list of arcs selected to form the spanning tree is stored in a local reference `st`.

```h
let st = ref [] in
```

We then write a function `cover` that examines an arc `(u,v,w)` and determines if it should be selected. For this, we compare the classes of vertices `u` and `v` using the function `UF.find`. If they are different, we add the arc to the list `st` and merge the two classes, using the function `UF.union`.

```h
let cover ((u, v, w) as e) =
    if UF.find uf u <> UF.find uf v then begin
        UF.union uf u v;
        st := e :: !st
    end
```

Finally, it only remains to examine each arc of the list `edges` using this function, and then return the resulting list.

```h
in
List.iter cover edges;
!st
```

The code of Kruskal's algorithm is given in program 103.

We illustrate Kruskal's algorithm with an example. Consider the graph on the top left of figure 13.7. Arcs are considered in increasing order of weight, as presented in the table on the right in the figure. For each arc, we indicate

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
Program 103 — Kruskal's algorithm.

```pascal
let spanning_tree g =
  let uf = UF.create (vertices g) in
  let compare (_,_,w1) (_,_,w2) = Stdlib.find compare w1 w2 in
  let edges = List.sort compare (edges g) in
  let st = ref [] in
  let cover ((u, v, w) as e) =
    if UF.find uf u <> UF.find uf v then begin
      UF.union uf u v;
      st := e :: !st
    end
  in
  List.iter cover edges;
!st
```

if it is selected (column “added”) and the resulting partition in the *union-find* structure. The spanning tree obtained in the end is illustrated at the bottom left of the figure. Its total weight is -2. It is important to note that arcs of the same weight may be considered in any order. Thus, we have examined here the arc 0 - 1 before the arc 1 - 4, but we could equally have done it the other way around. The spanning tree would have been different, but its total weight would have been the same.

A number of variations and optimizations are possible. Exercise 13.26 thus proposes interrupting the arc traversal once we have selected *N* - 1 arcs. Similarly, exercise 13.27 proposes using a priority queue rather than a sorting algorithm. Be that as it may, the key idea remains the use of the disjoint-set data structure.

Complexity

The cost of Kruskal's algorithm is split between that of sorting the set *E* of arcs and that of examining each arc. The first cost is *O(E log E)* in the worst case, if

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0/) —
<figure><img src="image_134.png" /><figcaption>Figure 13.7: Illustration of Kruskal's algorithm.</figcaption></figure>

we use a sorting algorithm of optimal complexity (see chapter 12). The second cost is $O(E)$ if we consider that the operations on the *union-find* structure are constant time (see chapter 8). The complexity of Kruskal's algorithm is therefore $O(E \log E)$.

## 13.5 Exercises

### Breadth-first Search

13.1 Modify the function `iter_bfs` to compute the distance, in number of arcs, between the root and each vertex discovered by the traversal. A simple approach consists in associating this distance with each vertex in the table `visited`. Return this table as the result of the function `iter_bfs`.

13.2 Modify the function `iter_bfs` so as to return a shortest path between the root and the vertex, for each vertex discovered by the traversal. For each vertex discovered, store the vertex that led to it, for instance, in a hash table. We can then recover the path from the vertex discovered to the root “by going backwards.”

13.3 Write a function that performs a breadth-first traversal of the nodes of a tree by imitating the breadth-first search of a graph.

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
13.4 Implement a cursor for binary trees (see section 9.3 *Cursors*) corresponding to breadth-first search.

13.5 Use the previous exercise to write a program that colors a graph with *k* colors without using the same color for two vertices linked by an arc. Proceed by *backtracking*, using the persistent nature of the cursor to backtrack. If no color assignment is possible, raise the exception *Not_found*.

Depth-first Search

13.6 Redo exercise 13.2 for depth-first search.

13.7 In many applications of depth-first search, we wish to traverse *all* the vertices of the graph and not merely those that can be reached from a vertex *s*. Using the function `iter vertex` that traverses all vertices of the graph, modify the function `iter DFS` so that every vertex of the graph is visited exactly once.

13.8 For an undirected graph, a *connected component* is a maximal set of vertices pairwise linked by a path. Explain how the variant of depth-first search proposed in the previous exercise determines the connected components of an undirected graph. Write the corresponding code.

13.9 Rewrite the function `iter DFS` using a `while` loop rather than a recursive function. Hint: Use a stack containing the vertices from which the depth-first search is to be performed. The code should resemble that of a breadth-first search, with the stack replacing the queue. There is, however, a difference in the treatment of the vertices already visited. Show that the vertices are not necessarily visited in the same order as in the recursive version.

13.10 Depth-first search can be modified to detect the presence of a cycle in the graph. When the function `iter DFS` comes across an already visited vertex, we do not know a *priori* if we have found a cycle. It could be that the vertex has already been reached by a parallel path. Instead of two states (reached/not reached), the marking of the vertices should therefore be modified to use three: reached/being visited/visited. Write a function `has_cycle` that determines the presence of a cycle in a graph. Modify the function so that it returns the list of vertices of the cycle found, if any.

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
**13.11** Let $G$ be a directed graph that does not contain a cycle. We call such graphs *directed acyclic graphs* (DAG). A *topological sort* of $G$ is a traversal of its vertices compatible with the arcs, that is, where a vertex $x$ is visited before a vertex $y$ if there is an arc $x \to y$. Modify program 100 so that it performs a topological sort, in the form of a function with the following type:

```h
topological_sort: (vertex -> unit) -> graph -> unit
```

Use a stack to which a vertex $v$ is added once the call to `visit v` terminates. Once the depth-first search is completed, apply the function passed as argument to all the elements of the stack.

**13.12** For a directed graph, a *strongly connected component* is any maximal set of vertices pairwise linked by a path. Kosaraju-Sharir's algorithm computes the set of strongly connected components of a graph $G$, in time $O(N+E)$, using two depth-first searches. It proceeds as follows:

* We begin by constructing the transposition $G^R$ of the graph $G$, that is, the graph having the same vertices as $G$ but with the arcs reversed (see exercise 7.1).

* Next, we traverse the vertices of $G^R$ according to a topological sort (see previous exercise). For each vertex $v$, we perform a depth-first search in $G$, starting from $v$, if $v$ is not already part of a strongly connected component. All the vertices that are reached by this search are then placed in a new component.

Write a function `scc: graph -> vertex list list` that uses this algorithm to compute the strongly connected components of a graph.

**13.13** We can use depth-first search to construct a perfect maze in an $n \times m$ grid. A perfect maze is one in which there is one and only one path between any two cells. We begin by considering the graph in which all cells of the grid are linked to their neighbors:

<figure><img src="image_135.png" /></figure>

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —
We then perform a depth-first search beginning from an arbitrary vertex (for example, the top left one). When we traverse the successors of a vertex, we do so in random order. Once the traversal is completed, the maze is obtained by declaring that we can pass from one vertex to another if the corresponding arc was taken during the depth-first search. Exercise 8.7 proposes another way of constructing a perfect maze.

**13.14** We consider here the problem of playing peg solitaire. This is a game of patience consisting of 32 pegs on a board with 33 holes. Initially, the pegs occupy all the holes, except the one at the center of the board.

<figure><img src="image_136.png" /></figure>

A peg can be removed from the board by having one of the four adjacent pegs jump over it, provided that the hole where the peg is to land is empty. Here are two possible moves starting from the initial position.

<figure><img src="image_137.png" /></figure>

The aim of the game is to arrive at a situation in which there is only one peg remaining, which occupies the central position.

<figure><img src="image_138.png" /></figure>

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
We can easily solve this problem using depth-first search on a graph in which
each vertex represents a board configuration.

* A board configuration can be represented by an integer whose $(i+7j)$-th bit indicates the presence of a peg at position $(i, j)$ for $0 \le i, j < 7$. Assuming a 64-bit machine, we let `type state = int`. Define two constants, `initial` and `final`, representing respectively the initial and final configurations of the game.

* A move can be represented by two configurations, the first one indicating the two pegs involved (the one that is moved and the one jumped over), and the second indicating the final position of the peg that was moved. For example, the first move above is represented as follows:

<figure><img src="image_139.png" /></figure>

Let `type move = state * state`. Construct a list containing the 88 possible moves (independent of what the configuration may be).

* Write a function `possible_move` to test if a move is possible in a given configuration and a function `move` that performs the move, if one is possible. Derive from it a function `iter S U C C` that traverses all possible moves in a given configuration.

* Modify the function `iter_dfs` of program 100 to interrupt the traversal once the configuration `final` is reached and return the path found, as in exercise 13.6 above.

* Finally, write a function `print` that displays a configuration on the screen, and use it to display the solution found, in the form of 32 successive configurations.

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
Note: The problem of peg solitaire is one of the problems discussed in the book by Cousineau and Mauny *The Functional Approach to Programming* [9, page 216].

**13.15** In the preceding exercise, depth-first search was a good solution because all the paths were of the same length. When this is not the case, and if we seek a path of minimal length, it is necessary to turn to breadth-first search instead. However, this may require a lot of memory. In this case, *Iterative deepening search* (IDS) offers a solution that has low memory requirements, but which is still able to find a path of minimal length. The idea is to execute successive depth-first searches, limited to maximal depths that are larger and larger.

* Explain why we can no longer memoize the vertices already visited during the depth-first searches. Consider the following graph and assume that we traverse it beginning at vertex 0.

<figure><img src="image_140.png" /></figure>

* Implement iterative deepening search as a function `ids: (vertex -> bool) -> graph -> vertex` that stops once it reaches a vertex that satisfies the function passed as argument.

* In case of the traversal of a complete binary tree starting at the root, show that iterative deepening search is no costlier than breadth-first and depth-first searches.

## Shortest Paths

**13.16** Generalize the function `dijkstra` for distances of an arbitrary type, for example in the form of a functor.

**13.17** Redo exercise 13.2 for Dijkstra's algorithm. Note: When a path is improved, the corresponding table must be updated.

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —
**13.18** Redo the steps of the Bellman-Ford algorithm, illustrated in figure 13.6, assuming that the weight of the arc $0 \to 4$ is -3.

**13.19** Warshall's algorithm (exercise 7.6) can be easily adapted to compute all shortest paths in a graph. This is the Floyd-Warshall algorithm. Let $V = \{0, \dots, N-1\}$. Initially, the distance $d_{i,j}$ between vertices $i$ and $j$ equals the weight of the corresponding arc, if there is one, and is otherwise $\infty$. The algorithm then computes the minimal distance $d_{i,j}$ from $i$ to $j$ in $O(N^3)$ time as follows:

$$ \begin{array}{l} \text{for } k \text{ from } 0 \text{ to } N-1 \\ \quad \text{for } i \text{ from } 0 \text{ to } N-1 \\ \quad \quad \text{for } j \text{ from } 0 \text{ to } N-1 \\ \quad \quad \quad d_{i,j} \leftarrow \min(d_{i,j}, d_{i,k} + d_{k,j}) \end{array} $$

Write a function that implements this algorithm. How can we treat the general case, when vertices are not necessarily the integers $0, \dots, N-1$?

**13.20** Although Dijkstra's algorithm allows us to find a path of minimal cost, it may visit a number of vertices in the graph unnecessarily. Suppose we use Dijkstra's algorithm to find the shortest path (in terms of distance) between two cities on a road network. The algorithm will explore the network in larger and larger concentric circles centered around the starting city, until it reaches the destination city. If we take the example of the French road network, and if we seek the shortest path between Paris and Nice, it is likely that almost all the cities of France will have been visited before Nice is reached.

The A\* algorithm allows us to solve this problem through a finer control of the order in which vertices are visited. This control takes the form of a heuristic that estimates the distance remaining between a vertex and the destination. The A\* algorithm is identical to Dijkstra's algorithm, except that the queue no longer uses distance to the source as priority, but rather the sum of this distance and the value given by the heuristic. If the heuristic is *admissible*, that is, if it never overestimates the remaining distance, then the A\* algorithm will find the shortest path. An admissible heuristic for the road-network example is the distance as a crow flies. If the heuristic always returns zero, we recover Dijkstra's algorithm.

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
Let us consider the example of an infinite grid in which each point is linked to its eight direct neighbors.

<figure><img src="image_141.png" /></figure>

The cost here is the Euclidean distance, that is, four neighbors are at distance 1, and the other four at distance $\sqrt{2}$. We can therefore use as heuristic the Euclidean distance between the point in question $(x, y)$ and the destination point $(x_d, y_d)$, that is, $\sqrt{(x-x_d)^2 + (y-y_d)^2}$.

<figure><img src="image_142.png" /></figure>

The figure above shows the result of the A\* algorithm finding a path from point (0,0), labelled 1, to point (5,10), labelled 103, the black nodes being excluded from the grid. As we see, only 103 nodes were visited. Dijkstra's

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —
algorithm would also have found the shortest path, but the total number of visited nodes would have been much larger (over 800).

Assume we are given a function **heuristic**: graph -> vertex -> vertex -> float, and implement the $A^*$ algorithm by modifying the code of program 101. Apply this to the preceding example and verify that we do indeed arrive at the result of the figure above.

**13.21** "French solitaire" differs from that of exercise 13.14, called "English solitaire," by virtue of four additional pegs present in positions (1, 1), (5, 1), (1, 5), and (5, 5).

<figure><img src="image_143.png" /></figure>

French solitaire is, surprisingly, much harder to solve than English solitaire. Furthermore, there exist combinations of starting/final configurations that do not admit of a solution, for example, the configuration in which the central peg is removed at the start and the last peg must be placed at the center<fnref n="3" />. We will accordingly loosen the conditions on the positions of the initial and final pegs.

One approach to solving French solitaire consists in using the $A^*$ algorithm presented in the preceding exercise. Modify the code of exercise 13.14 as follows:

* Add the 8 new moves linked to the 4 new pegs to the list of all possible moves.

* Write a heuristic function that evaluates how "compactly" pegs are laid out on the board. Given the representation of the grid as an integer, we may compute a rough approximation by counting the number of bit

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —

<footnote n="3">³It is only recently that we have rediscovered how French solitaire was played in the 18th century: Initially, *all* cells contain a peg, and the player removes a peg of her choice. This peg can then be placed back on the board at any moment during the game. With these rules, there is then a solution in which the central peg is taken out at the start and the final peg is placed in the center. Here, however, we maintain the same rules as in English solitaire.</footnote>
inversions, that is, the number of times that we pass from a 1 bit to a
0 bit, and vice versa.

* Finally, replace depth-first search with the A\* algorithm. Here is an example of initial and final configurations for which there exists a solution.

<figure><img src="image_144.png" /></figure>

As in exercise 13.14, present the solution by displaying all successive configura-
tions of the board.

**13.22** It is possible to combine the ideas of iterative deepening search (exercise 13.15) and the A\* algorithm (exercise 13.20), to obtain an algorithm known as IDA\* (iterative deepening A\*). The idea is as follows: At each iteration, we execute a depth-first search of all the states where the heuristic does not exceed a given threshold, without seeking to memoize already visited states. Once the search is over, we take as the new threshold the smallest value of the heuristic obtained during the previous search. In other words, each iteration pushes the threshold of the heuristic a little further until a path is found or no vertex takes a value greater than the threshold. To initiate the process, it suffices to take as initial threshold the heuristic for the root vertex. As with the A\* algorithm, an admissible heuristic, that is, one that never overestimates the remaining distance, guarantees that IDA\* will find a shortest path. The IDA\* algorithm was discovered by Korf [17] in 1985.

Assume we are given a function **heuristic**: graph -> vertex -> vertex -> float,
and modify the code of exercise 13.15 to implement the IDA\* algorithm. The
following exercise proposes an application of this algorithm.

**13.23** A good application of the IDA\* algorithm proposed in the preceding exercise is the *game of fifteen*. The game consists in a 4 × 4 grid containing 15 tiles, numbered from 1 to 15. The sixteenth position is empty, and we may place in it one of the adjacent tiles. Initially, the tiles are scrambled, for example as follows:

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —
<table><tr><td>14</td><td>1</td><td>9</td><td>6</td></tr><tr><td>4</td><td>8</td><td>12</td><td>5</td></tr><tr><td>7</td><td>2</td><td>3</td><td></td></tr><tr><td>10</td><td>11</td><td>13</td><td>15</td></tr></table>

The aim is to place the tiles in increasing order<fnref n="4" />.

<table><tr><td></td><td>1</td><td>2</td><td>3</td></tr><tr><td>4</td><td>5</td><td>6</td><td>7</td></tr><tr><td>8</td><td>9</td><td>10</td><td>11</td></tr><tr><td>12</td><td>13</td><td>14</td><td>15</td></tr></table>

To solve the game, we will use the following heuristic: For each tile, we compute the distance that separates it from its final position, in terms of number of rows and columns. This is called the Manhattan distance. In the preceding example, tile 7 is initially at distance 4 from its final position because it has to be moved up by one row, and right by three columns. For a given configuration, we then take the sum of the Manhattan distances of the fifteen tiles, which is a lower bound on the number of moves that remain to be done. This is an admissible heuristic, and we will thus be able to find a smallest solution. Solve the game using the IDA\* algorithm. The example of the game of fifteen is presented in the article by Korf that introduced the IDA\* algorithm.

## Minimum Spanning Tree

**13.24** Show that the set of arcs returned by Kruskal's algorithm does indeed form a tree.

**13.25** Show that the result of Kruskal's algorithm is indeed a minimum spanning tree.

**13.26** Improve the code of program 103 so that it stops as soon as $N - 1$ arcs have been added to the spanning tree.

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0/) —

<footnote n="4">⁴Not all initial configurations admit a solution.</footnote>
**13.27** Reimplement program 103 using a priority queue rather than a sorting algorithm. How does it differ from the previous exercise?

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —
— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —
Bibliography

[1] Stephen Adams. Functional pearls: Efficient sets – a balancing act. *Journal of Functional Programming*, 3(4):553–561, October 1993. Expanded version available as Technical Report CSTR 92-10, University of Southampton.

[2] G. M. Adel'son-Vel'skiï and E. M. Landis. *An algorithm for the organization of information*. *Soviet Mathematics-Doklady*, 3(5):1259–1263, September 1962.

[3] Alfred V. Aho, John E. Hopcroft, and Jeffrey Ullman. *Data Structures and Algorithms*. Addison-Wesley Longman Publishing Co., Inc., Boston, MA, USA, 1983.

[4] Henry G. Baker. Shallow binding makes functional arrays fast. *SIGPLAN Not.*, 26(8):145–147, 1991.

[5] Hans-Juergen Boehm, Russell R. Atkinson, and Michael F. Plass. Ropes: An alternative to strings. *Software - Practice and Experience*, 25(12):1315–1330, 1995.

— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by-sa/4.0) —
[6] John Conway. Doomsday rule. http://en.wikipedia.org/wiki/Doomsday_algorithm.

[7] Thomas H. Cormen, Charles E. Leiserson, Ronald L. Rivest, and Clifford Stein. *Introduction to Algorithms, Second Edition*. The MIT Press, September 2001.

[8] Guy Cousineau and Michel Mauny. *Approche fonctionnelle de la programmation*. Ediscience International, 1995.

[9] Guy Cousineau and Michel Mauny. *The Functional Approach to Programming*. Cambridge University Press, 1998. Traduction anglaise de [8].

[10] Jeremy Gibbons and Oege de Moor, editors. *The Fun of Programming*. Cornerstones in Computing. Palgrave, 2003.

[11] Jr. Henry G. Baker. Shallow binding in Lisp 1.5. *Commun. ACM*, 21(7):565–569, 1978.

[12] Gérard Huet. The Zipper. *Journal of Functional Programming*, 7(5):549–554, September 1997.

[13] D. E. Knuth. *The Art of Computer Programming. Volume 1: Fundamental Algorithms*. Addison-Wesley, 1968.

[14] D. E. Knuth. *The Art of Computer Programming. Volume 2: Seminumerical Algorithms*. Addison-Wesley, 1969.

[15] D. E. Knuth. *The Art of Computer Programming. Volume 3: Sorting and Searching*. Addison-Wesley, 1973.

[16] Donald E. Knuth. *The Art of Computer Programming, volume 4A: Combinatorial Algorithms, Part 1*. Addison-Wesley Professional, 1st edition, 2011.

[17] Richard E. Korf. Iterative-deepening-A*: An optimal admissible tree search. In *IJCAI*, pages 1034–1036, 1985.
[18] Donald R. Morrison. *PATRICIA—Practical Algorithm To Retrieve Information Coded in Alphanumeric*. J. ACM, 15(4):514–534, 1968.

[19] Chris Okasaki. *Purely Functional Data Structures*. Cambridge University Press, 1998.

[20] Chris Okasaki and Andrew Gill. *Fast Mergeable Integer Maps*. In *Workshop on ML*, pages 77–86, September 1998.

[21] Robert Sedgewick and Kevin Wayne. *Algorithms, 4th Edition*. Addison-Wesley, 2011.

[22] Daniel Dominic Sleator and Robert Endre Tarjan. *Self Adjusting Heaps*. SIAM J. Comput., 15(1):52–69, February 1986.

[23] Robert Endre Tarjan. *Efficiency of a Good But Not Linear Set Union Algorithm*. J. ACM, 22(2):215–225, 1975.

[24] Henry S. Warren. *Hacker's Delight*. Addison-Wesley Professional, July 2002.
— openly licensed via [CC BY SA 4.0](https://creativecommons.org/licenses/by/4.0/) —
Index

O (big O notation), 159

&& , 13

() , 33

:= , 20

<> , 153

= , 153

== , 153

|| , 13

!, 20

_, 12, 52

accumulator, 79

algorithm

Bellman-Ford, 423

Dijkstra, 420

Euclid, 343

extended Euclidean, 344

Floyd-Warshall, 310, 438

Kosaraju-Sharir, 434

Kruskal, 428

aliasing, 46

and , 39

appel

terminal, 376

append (List), 136

arc, 295

arithmetic, 343

Array (OCaml library)

init , 48

iteri , 177, 299

length , 45

make , 42

make_matrix , 46

set , 42

sort , 49

array , 40, 44, 165

persistent , 197

resizeable , 165, 277

asr , 13

— openly licensed via CC BY SA 4.0 —
assert, 65

Assert_failure, 65

A*, 438, 441

AVL, 188

AVL tree, 216

backtracking, 120

backtracking, 118

Baker, Henry G., 201

balancing, 188, 217

begin, 40

Bellman-Ford algorithm, 423

Bézout, 344

BFS, 413

big *O*, notation, 159

big-endian, 248

binary search, 134

binomial coefficient, 364

bit

low order, set, 178

lowest-order set, 252

operations on bits, 13

bit vector, 171

block, 42, 150

bool (type), 13

call

function, 8

stack, 66

tail, 66, 73

cardioid, 27

char (type), 14

closure, 155

transitive closure, 310

command line, 10

compare (Stdlib), 223

compilation, 83

separate, 93

compilation unit, 85

complexity

amortized, 236, 274

complexity, 157

amortized, 160, 170, 289

component

connected, 433

strongly connected, 434

constructor, 104

Cousineau, Guy, 437

Curry, Haskell, 54

cursor, 335, 433

cycle, 296

negative, 425

cycles

detection of, 433

DAG, 117, 434

dequeue, 290

dereference, 20

DFS, 416

dictionary, 225, 237, 247, 260

Dijkstra

algorithm, 420

Dutch national flag, 404

disjoint sets, 313

divide and conquer, 378

dynamic programming, 357

edge, 295

multiple edges, 311

encapsulation, 95, 99

endianness, 248

eprintf (Printf), 16

equality

physical, 153

structural, 153

Eratosthenes, 40
Euclid, 343

evaluation, order, 47

exception, 58, 59, 62

execution time, 157

exit, 74

exn, 149

exponentiation by squaring, 345

false, 13

Fibonacci, 197, 205, 225, 344, 355, 357

file, 61

filter (List), 137

float (Random), 20

float (type), 21

floating-point number, 20

floating-point numbers, 16

Floyd-Warshall algorithm, 310, 438

fold right (List), 79

for, 16, 25

function

anonymous, 47, 58

call, 8

declaration, 33

first class, 47

recursive, 33

function, 70

functor, 95, 96

game of fifteen, 441

garbage collector, 155

GC, 155

gcd, 343

Gill, Andrew, 260

graph, 295, 413

adjacency lists, 301

adjacency matrix, 297

dense, 300

strongly connected components, 434

topological sort, 434

Graphics, 27, 32

Hanoi, tower of, 133

hash table, 227

hash-consing, 260, 366

heap, 274, 276, 285

heapsort, 290

heuristic, 438

Huet, Gérard, 338

IDA\*, 441

IDS, 437, 441

ignore, 27

immutable, 66

in_channel, 60

include, 142

infinity, 22

init (Array), 48

input/output, 58

int (Random), 134

int (type), 12

interface, 83, 92

iter (List), 70

iterator, 74

iteri (Array), 177, 299

Knuth shuffle, 134, 405

Korf, Richard E., 441

Kosaraju-Sharir, algorithm, 434

Kruskal's algorithm, 428

label, 296

Lamé, 344

land, 13

leap (year), 8

length (Array), 45

let, 8, 9, 16, 23

library
Array, 42, 45, 46, 48, 49, 177, 299
List, 70, 78, 79, 136, 137, 325, 377, 430
Format, 16
Queue, 272
Random, 16, 20, 134
Seq, 338
Stdlib, 11, 223
Unix, 157
Weak, 369

linking, 91

Lisp, 367

List (OCaml library)

append, 136

filter, 137

fold-right, 79

iter, 70

map, 78

mem, 136

rev_append, 325, 377

sort, 430

list, 66

adjacency, 301

linked, 267

little-endian, 248

lnot, 13

local variable, 16

Logo, programming language, 95

loop

body, 19

for, 19, 25

index, 19

while, 44

lor, 13

lsl, 13

lsr, 13

lxor, 13

make (Array), 42

make_matrix (Array), 46

Mandelbrot set, 33

Manhattan distance, 442

map (List), 78

marshalling, 61

matching

pattern, 107

matrix, 45, 137, 138, 350

adjacency, 297

Mauny, Michel, 437

max_int, 12

maze, 321, 434

McCarthy, John, 73

mem (List), 136

memoization, 357

Michie, Donald, 359

min_int, 12

modular arithmetic, 346

module, 83

Monte-Carlo, 16

Morrison, Donald R., 260

mutable, 56

n-tuple, 47

nan, 22

neg_infinity, 22

not, 13

ntz, number of trailing zeros, 179

Okasaki, Chris, 260

opam, 3

open, 27

open, 27

order

evaluation, 57

higher, 47, 48

out_channel, 60
pair, 47, 51
partial application, 54
Pascal's triangle, 365
Pascal, triangle, 365
path
  in a graph, 296
  shortest, 419
Patricia tree, 248
pattern, 49
  wildcard, 52
pattern matching, 10, 49, 66, 69, 70, 102, 106, 108
  or patterns, 110
persistence, 118
π, 16
polymorphism, 74, 79
pop, population count, 179, 204
print_char, 16
print_float, 16
print_int, 16
print_string, 15
printf (OCaml library)
  eprintf, 16
private, 100

queens, problem of N, 118
Queue (OCaml library), 272
queue, 267
  double-ended, 290
  persistent, 272
  priority, 267, 274

Random (OCaml library), 16
  float, 20
  int, 134
record, 54
ref, 16, 57
reference, 16, 19, 57

rev_append (List), 325, 377
rope, 180
rosace, 132

search
  binary, in an array, 134
  breadth-first, 413
  depth-first, 416
  iterative deepening, 437
  shortest path, 419

Seq (OCaml library), 338
sequence
  expression, 16, 20

set (Array), 42
set, Mandelbrot, 33
sharing, 112, 117

shuffle, 134, 405
side effects, 16, 18
Sierpinski, Wacław, 142
sieve of Eratosthenes, 40

smart constructor, 117, 217, 219, 255
solitaire, 435

French, 440
sort

heap, 290, 393
insertion, 376
merge, 386
quick, 378
topological, 434

sort (Array), 49
sort (List), 430
sorting, 375

optimal complexity, 400
stack

call, 71
overflow, 73

Stdlib (OCaml library), 11
  compare, 223

— openly licensed via CC BY SA 4.0 —
string, 180

string (type), 14

Swift, Jonathan, 248

Sys, 11

tail recursion, 376

Tarjan, Robert Endre, 318

text editor, 194

The Art of Computer Programming, 206,
344, 354, 400

times (Unix), 157

tree, 112

2-3-4, 227

AVL, 216

binary, 112, 118

binary search, 208

Patricia, 248

prefix, 238, 239

red-black, 227

self-balancing, 285

spanning, 428

trie, 238

true, 13

truncate, 22

TryOCaml, 6

type

abstract, 95, 96

algebraic, 102, 104

enumerated, 104

private, 95, 101

variable, 111

union-find, 313

unit, 33

Unix (OCaml library)

times, 157

variable

global, 23

local, 18, 23

scope, 18

vector

bit, 171

vertex, 295

von Koch, snowflake, 132

Warren, Henry S., 180

Weak (OCaml library), 369

weak reference, 369

while, 40, 44

zipper, 323