anansi.algorithms.dijkstra source code

1 module anansi.algorithms.dijkstra;
2 
3 import anansi.container.priorityqueue;
4 import anansi.algorithms.bfs;
5 import anansi.algorithms.relax : relax;
6 import anansi.algorithms.vertex_queue : VertexQueue;
7 import anansi.traits;
8 import anansi.types;
9 import std.exception, std.math;
10 
11 
12 version(unittest) { import std.stdio; import anansi.adjacencylist; }
13 
14 /**
15  * A priority queue item for sorting vertices by the cost to get to them.
16  */
17 
18 package struct DijkstraQueue(GraphT, DistanceMapT) {
19     static assert(isGraph!GraphT);
20     static assert(isPropertyMap!(DistanceMapT, GraphT.VertexDescriptor, real));
21 
22     alias Vertex = GraphT.VertexDescriptor;
23 
24     package static struct Item {
25         public Vertex vertex;
26         public real cost;
27 
28         public int opCmp(ref const(Item) other) const { 
29             return cast(int) sgn(other.cost - this.cost);
30         }
31     }
32 
33     this(ref const(DistanceMapT) distances) {
34         _distanceMap = &distances;
35     }
36 
37     public void push(Vertex v) {
38         const auto distance = (*_distanceMap)[v];
39         _queue.push(Item(v, distance));
40     }
41 
42     public void pop() {
43         _queue.pop();
44     }
45 
46     public Vertex front() const  {
47         return _queue.front.vertex;
48     }
49 
50     public void updateVertex(Vertex v) {
51         const auto distance = (*_distanceMap)[v];
52         _queue.updateIf((ref Item x) => x.vertex == v, 
53                         (ref Item i) => i.cost = distance);
54     }
55 
56     @property public bool empty() const {
57         return _queue.empty;
58     }
59 
60     @property public size_t length() const {
61         return _queue.length;
62     } 
63 
64     PriorityQueue!Item _queue;
65     const (DistanceMapT*) _distanceMap; 
66 }
67 
68 unittest {
69     writeln("Dijkstra: queue items are ordered by increasing cost.");
70     alias G = AdjacencyList!();
71     alias V = G.VertexDescriptor;
72     alias Item = DijkstraQueue!(G, real[V]).Item;
73 
74     auto a = Item(0, 0.1);
75     auto b = Item(0, 0.2);
76 
77     assert (a > b);
78     assert (!(b > a));
79     assert (b < a);
80     assert (!(a < b));
81     assert (a != b);
82 }
83 
84 /**
85  * A default implementation of the Dijkstra visitor.
86  */
87 struct NullDijkstraVisitor(GraphT) {
88     alias Vertex = GraphT.VertexDescriptor;
89     alias Edge = GraphT.EdgeDescriptor;
90 
91     void initVertex(ref const(GraphT) g, Vertex v) {}
92     void discoverVertex(ref const(GraphT) g, Vertex v) {}
93     void examineVertex(ref const(GraphT) g, Vertex v) {}
94     void examineEdge(ref const(GraphT) g, Edge e) {}
95     void edgeRelaxed(ref const(GraphT) g, Edge e) {}
96     void edgeNotRelaxed(ref const(GraphT) g, Edge e) {}
97     void finishVertex(ref const(GraphT) g, Vertex e) {}
98 }
99 
100 /**
101  * A BFS visitor that transforms a normal breadth-first search algoritm
102  * into Dijkstra's shortest paths.
103  *
104  * Remark: Not every `sumFunction` works correctly.
105  */
106 package struct DijkstraBfsVisitor(GraphT,
107                                   QueueT,
108                                   DijkstraVisitorT,
109                                   DistanceMapT,
110                                   PredecessorMapT,
111                                   WeightMapT,
112                                   alias sumFunction = (a, b) => a+b) {
113     alias Vertex = GraphT.VertexDescriptor;
114     alias Edge = GraphT.EdgeDescriptor;
115 
116     static assert(isReadablePropertyMap!(WeightMapT, Edge, real));
117     static assert(isPropertyMap!(DistanceMapT, Vertex, real));
118     static assert(isPropertyMap!(PredecessorMapT, Vertex, Vertex));
119 
120     this(ref DijkstraVisitorT visitor,
121          ref DistanceMapT distanceMap,
122          ref const(WeightMapT) weightMap,
123          ref PredecessorMapT predecessorMap,
124          ref QueueT queue) {
125         _visitor = &visitor;
126         _distanceMap = &distanceMap;
127         _weightMap = &weightMap;
128         _predecessorMap = &predecessorMap;
129         _queue = &queue;
130     }
131 
132     /**
133      * Passes the call through to the supplied Dijkstra visitor.
134      */
135     void initVertex(ref const(GraphT) g, Vertex v) {
136         _visitor.initVertex(g, v);
137     }
138 
139     void discoverVertex(ref const(GraphT) g, Vertex v) {
140         _visitor.discoverVertex(g, v);
141     }
142 
143     void examineVertex(ref const(GraphT) g, Vertex v) {
144         _visitor.examineVertex(g, v);
145     }
146 
147     void examineEdge(ref const(GraphT) g, Edge e) {
148         // TODO: Document the behavior of this `static if`.
149         static if (__traits(compiles, (*_weightMap)[e])) {
150             auto weight = (*_weightMap)[e];
151         } else {
152             auto weight = (*_weightMap)(e);
153         }
154         enforce (weight >= 0.0);
155         _visitor.examineEdge(g, e);
156     }
157 
158     void treeEdge(ref const(GraphT) g, Edge e) {
159         bool decreased = relax(g, *_weightMap,
160                                   *_distanceMap,
161                                   *_predecessorMap, e);
162         if (decreased)
163           _visitor.edgeRelaxed(g, e);
164         else
165           _visitor.edgeNotRelaxed(g, e);
166     }
167 
168     void nonTreeEdge(ref const(GraphT) g, Edge e) {}
169 
170     void greyTarget(ref const(GraphT) g, Edge e) {
171         const Vertex v = g.target(e);
172         auto oldDistance = (*_distanceMap)[v];
173 
174         bool decreased = relax(g, *_weightMap,
175                                   *_distanceMap,
176                                   *_predecessorMap, e);
177         if (decreased) {
178             _queue.updateVertex(v);
179             _visitor.edgeRelaxed(g, e);
180         } else {
181             _visitor.edgeNotRelaxed(g, e);
182         }
183     }
184 
185     void blackTarget(ref const(GraphT) g, Edge e) {}
186 
187     void finishVertex(ref const(GraphT) g, Vertex v) {
188         _visitor.finishVertex(g, v);
189     }
190 
191     private DijkstraVisitorT* _visitor;
192     private DistanceMapT* _distanceMap;
193     private const(WeightMapT*) _weightMap;
194     private PredecessorMapT* _predecessorMap;
195     private QueueT* _queue;
196 }
197 
198 /**
199  *
200  */
201 template dijkstraShortestPaths(GraphT,
202                                VertexDescriptorT,
203                                VisitorT = NullDijkstraVisitor!GraphT,
204                                WeightMapT = real[VertexDescriptorT],
205                                PredecessorMapT = VertexDescriptorT[VertexDescriptorT],
206                                alias sumFunction = (a, b) => a+b) {
207     void dijkstraShortestPaths(ref const(GraphT) g,
208                                VertexDescriptorT src,
209                                ref const(WeightMapT) weights,
210                                ref PredecessorMapT predecessorMap,
211                                VisitorT visitor = VisitorT.init) {
212 
213         static if (is(VertexDescriptorT == size_t)) {
214             auto colourMap = new Colour[g.vertexCount];
215             auto distanceMap = new real[g.vertexCount];
216         }
217         else {
218             real[VertexDescriptorT] distanceMap;
219             Colour[VertexDescriptorT] colourMap;
220         }
221 
222         dijkstraShortestPaths(g, src,
223                               weights,
224                               predecessorMap,
225                               visitor,
226                               colourMap,
227                               distanceMap);
228     }
229 }
230 
231 /**
232  *
233  */
234 template dijkstraShortestPaths(GraphT,
235                                VisitorT,
236                                VertexDescriptorT,
237                                ColourMapT,
238                                PredecessorMapT,
239                                WeightMapT,
240                                DistanceMapT,
241                                alias sumFunction = (a, b) => a+b) {
242     void dijkstraShortestPaths(ref const(GraphT) g,
243                                VertexDescriptorT src,
244                                ref const(WeightMapT) weights,
245                                ref PredecessorMapT predecessorMap,
246                                VisitorT visitor,
247                                ref ColourMapT colourMap,
248                                ref DistanceMapT distanceMap) {
249 
250         foreach(v; g.vertices) {
251             visitor.initVertex(g, v);
252             distanceMap[v] = real.infinity;
253             predecessorMap[v] = v;
254             colourMap[v] = Colour.White;
255         }
256         distanceMap[src] = 0.0;
257 
258         dijkstraShortestPathsNoInit(g, src, visitor,
259                                     colourMap,
260                                     weights,
261                                     predecessorMap,
262                                     distanceMap);
263     }
264 }
265 
266 template dijkstraShortestPathsNoInit(GraphT,
267                                      VertexDescriptorT,
268                                      VisitorT = NullDijkstraVisitor!GraphT,
269                                      ColourMapT = Colour[VertexDescriptorT],
270                                      PredecessorMapT = VertexDescriptorT[VertexDescriptorT],
271                                      WeightMapT = real[VertexDescriptorT],
272                                      DistanceMapT = real[VertexDescriptorT],
273                                      alias sumFunction = (a, b) => a+b) {
274     static assert(isGraph!GraphT);
275 
276     alias EdgeDescriptorT = GraphT.EdgeDescriptor;
277 
278     static assert(isPropertyMap!(ColourMapT, VertexDescriptorT, Colour));
279     static assert(isPropertyMap!(PredecessorMapT, VertexDescriptorT, VertexDescriptorT));
280     static assert(isReadablePropertyMap!(WeightMapT, EdgeDescriptorT, real));
281     static assert(isPropertyMap!(DistanceMapT, VertexDescriptorT, real));
282 
283 
284     void dijkstraShortestPathsNoInit(ref const(GraphT) g,
285                                      VertexDescriptorT src,
286                                      VisitorT visitor,
287                                      ref ColourMapT colourMap,
288                                      ref const(WeightMapT) weights,
289                                      ref PredecessorMapT predecessorMap,
290                                      ref DistanceMapT distanceMap) {
291 
292         alias QueueT = VertexQueue!(GraphT, DistanceMapT);
293         auto queue = QueueT(distanceMap);
294 
295         alias Dijkstra = DijkstraBfsVisitor!(GraphT,
296                                              QueueT,
297                                              VisitorT,
298                                              DistanceMapT,
299                                              PredecessorMapT,
300                                              WeightMapT,
301                                              sumFunction);
302         breadthFirstVisit(
303             g, src, colourMap,
304             queue,
305             Dijkstra(visitor, distanceMap, weights, predecessorMap, queue));
306     }
307 }
308 
309 version (unittest) {
310     import std.array, std.algorithm, std.conv, std.stdio;
311     import anansi.adjacencylist, anansi.traits, anansi.container.set;
312     private alias G = AdjacencyList!(VecS, VecS, DirectedS, char, string);
313     private alias Vertex = G.VertexDescriptor;
314     private alias Edge = G.EdgeDescriptor;
315 }
316 
317 unittest {
318     writeln("Dijkstra: no edges shouldn't crash");
319     G g;
320     real[Edge] weights;
321     Vertex[Vertex] predecessors;
322     g.addVertex('a');
323     dijkstraShortestPaths(g, g.vertices.front, weights, predecessors);
324 }
325 
326 unittest {
327     writeln("Dijkstra: shorter paths with multiple hops are preferred over long path with one hop.");
328 
329     //  / ----- 4 ----- \
330     // a - 1 -> b - 1 -> c - 2 -> d
331     //  \        \ ----- 5 ----- /
332     //    3 -> e
333     //
334 
335     G g;
336     real[Edge] weights;
337     Vertex[Vertex] predecessor;
338 
339     auto a = g.addVertex('a');
340     auto b = g.addVertex('b');
341     weights[g.addEdge(a, b, "a -> b").edge] = 1.0;
342 
343     auto c = g.addVertex('c');
344     weights[g.addEdge(b, c, "b -> c").edge] = 1.0;
345     weights[g.addEdge(a, c, "a -> c").edge] = 4.0;
346 
347     auto d = g.addVertex('d');
348     weights[g.addEdge(c, d, "c -> d").edge] = 2.0;
349     weights[g.addEdge(b, d, "b -> d").edge] = 5.0;
350 
351     auto e = g.addVertex('e');
352     weights[g.addEdge(a, e, "a -> e").edge] = 3.0;
353 
354 
355     dijkstraShortestPaths(g, g.vertices.front, weights, predecessor);
356 
357     assert(predecessor[a] == a);
358     assert(predecessor[b] == a);
359     assert(predecessor[c] == b);
360     assert(predecessor[d] == c);
361     assert(predecessor[e] == a);
362 }