1 | /* $Id$ $Revision$ */ |
2 | /* vim:set shiftwidth=4 ts=8: */ |
3 | |
4 | /************************************************************************* |
5 | * Copyright (c) 2011 AT&T Intellectual Property |
6 | * All rights reserved. This program and the accompanying materials |
7 | * are made available under the terms of the Eclipse Public License v1.0 |
8 | * which accompanies this distribution, and is available at |
9 | * http://www.eclipse.org/legal/epl-v10.html |
10 | * |
11 | * Contributors: See CVS logs. Details at http://www.graphviz.org/ |
12 | *************************************************************************/ |
13 | |
14 | |
15 | /****************************************** |
16 | |
17 | Dijkstra algorithm |
18 | Computes single-source distances for |
19 | weighted graphs |
20 | |
21 | ******************************************/ |
22 | |
23 | |
24 | #include "bfs.h" |
25 | #include "dijkstra.h" |
26 | #include <limits.h> |
27 | #include <stdlib.h> |
28 | /* #include <math.h> */ |
29 | |
30 | #define MAX_DIST (double)INT_MAX |
31 | |
32 | typedef DistType Word; |
33 | |
34 | #define LOOP while(TRUE) |
35 | |
36 | /* This heap class is suited to the Dijkstra alg. |
37 | data[i]=vertexNum <==> index[vertexNum]=i |
38 | */ |
39 | |
40 | #define left(i) (2*(i)) |
41 | #define right(i) (2*(i)+1) |
42 | #define parent(i) ((i)/2) |
43 | #define insideHeap(h,i) ((i)<h->heapSize) |
44 | #define greaterPriority(h,i,j,dist) (dist[h->data[i]]<dist[h->data[j]]) |
45 | #define assign(h,i,j,index) {h->data[i]=h->data[j]; index[h->data[i]]=i;} |
46 | #define exchange(h,i,j,index) {int temp; \ |
47 | temp=h->data[i]; \ |
48 | h->data[i]=h->data[j]; \ |
49 | h->data[j]=temp; \ |
50 | index[h->data[i]]=i; \ |
51 | index[h->data[j]]=j; \ |
52 | } |
53 | |
54 | typedef struct { |
55 | int *data; |
56 | int heapSize; |
57 | } heap; |
58 | |
59 | static void heapify(heap * h, int i, int index[], Word dist[]) |
60 | { |
61 | int l, r, largest; |
62 | while (1) { |
63 | l = left(i); |
64 | r = right(i); |
65 | if (insideHeap(h, l) && greaterPriority(h, l, i, dist)) |
66 | largest = l; |
67 | else |
68 | largest = i; |
69 | if (insideHeap(h, r) && greaterPriority(h, r, largest, dist)) |
70 | largest = r; |
71 | |
72 | if (largest == i) |
73 | break; |
74 | |
75 | exchange(h, largest, i, index); |
76 | i = largest; |
77 | } |
78 | } |
79 | |
80 | #ifdef OBSOLETE |
81 | /* originally, the code called mkHeap to allocate space, then |
82 | * initHeap to realloc the space. This is silly. |
83 | * Now we just call initHeap. |
84 | */ |
85 | static void mkHeap(heap * h, int size) |
86 | { |
87 | h->data = N_GNEW(size, int); |
88 | h->heapSize = 0; |
89 | } |
90 | #endif |
91 | |
92 | static void freeHeap(heap * h) |
93 | { |
94 | if (h->data) free(h->data); |
95 | } |
96 | |
97 | static void |
98 | initHeap(heap * h, int startVertex, int index[], Word dist[], int n) |
99 | { |
100 | int i, count; |
101 | int j; /* We cannot use an unsigned value in this loop */ |
102 | /* h->data=(int*) realloc(h->data, (n-1)*sizeof(int)); */ |
103 | if (n == 1) h->data = NULL; |
104 | else h->data = N_GNEW(n - 1, int); |
105 | h->heapSize = n - 1; |
106 | |
107 | for (count = 0, i = 0; i < n; i++) |
108 | if (i != startVertex) { |
109 | h->data[count] = i; |
110 | index[i] = count; |
111 | count++; |
112 | } |
113 | |
114 | for (j = (n - 1) / 2; j >= 0; j--) |
115 | heapify(h, j, index, dist); |
116 | } |
117 | |
118 | static boolean (heap * h, int *max, int index[], Word dist[]) |
119 | { |
120 | if (h->heapSize == 0) |
121 | return FALSE; |
122 | |
123 | *max = h->data[0]; |
124 | h->data[0] = h->data[h->heapSize - 1]; |
125 | index[h->data[0]] = 0; |
126 | h->heapSize--; |
127 | heapify(h, 0, index, dist); |
128 | |
129 | return TRUE; |
130 | } |
131 | |
132 | static void |
133 | increaseKey(heap * h, int increasedVertex, Word newDist, int index[], |
134 | Word dist[]) |
135 | { |
136 | int placeInHeap; |
137 | int i; |
138 | |
139 | if (dist[increasedVertex] <= newDist) |
140 | return; |
141 | |
142 | placeInHeap = index[increasedVertex]; |
143 | |
144 | dist[increasedVertex] = newDist; |
145 | |
146 | i = placeInHeap; |
147 | while (i > 0 && dist[h->data[parent(i)]] > newDist) { /* can write here: greaterPriority(i,parent(i),dist) */ |
148 | assign(h, i, parent(i), index); |
149 | i = parent(i); |
150 | } |
151 | h->data[i] = increasedVertex; |
152 | index[increasedVertex] = i; |
153 | } |
154 | |
155 | void dijkstra(int vertex, vtx_data * graph, int n, DistType * dist) |
156 | { |
157 | int i; |
158 | heap H; |
159 | int closestVertex, neighbor; |
160 | DistType closestDist, prevClosestDist = INT_MAX; |
161 | static int *index; |
162 | |
163 | #ifdef OBSOLETE |
164 | mkHeap(&H, n); |
165 | #endif |
166 | index = (int *) realloc(index, n * sizeof(int)); |
167 | |
168 | /* initial distances with edge weights: */ |
169 | for (i = 0; i < n; i++) |
170 | dist[i] = (DistType) MAX_DIST; |
171 | dist[vertex] = 0; |
172 | for (i = 1; i < graph[vertex].nedges; i++) |
173 | dist[graph[vertex].edges[i]] = (DistType) graph[vertex].ewgts[i]; |
174 | |
175 | initHeap(&H, vertex, index, dist, n); |
176 | |
177 | while (extractMax(&H, &closestVertex, index, dist)) { |
178 | closestDist = dist[closestVertex]; |
179 | if (closestDist == MAX_DIST) |
180 | break; |
181 | for (i = 1; i < graph[closestVertex].nedges; i++) { |
182 | neighbor = graph[closestVertex].edges[i]; |
183 | increaseKey(&H, neighbor, |
184 | closestDist + |
185 | (DistType) graph[closestVertex].ewgts[i], index, |
186 | dist); |
187 | } |
188 | prevClosestDist = closestDist; |
189 | } |
190 | |
191 | /* For dealing with disconnected graphs: */ |
192 | for (i = 0; i < n; i++) |
193 | if (dist[i] == MAX_DIST) /* 'i' is not connected to 'vertex' */ |
194 | dist[i] = prevClosestDist + 10; |
195 | freeHeap(&H); |
196 | } |
197 | |
198 | /* Dijkstra bounded to nodes in *unweighted* radius */ |
199 | int |
200 | dijkstra_bounded(int vertex, vtx_data * graph, int n, DistType * dist, |
201 | int bound, int *visited_nodes) |
202 | /* make dijkstra, but consider only nodes whose *unweighted* distance from 'vertex' */ |
203 | /* is at most 'bound' */ |
204 | /* MON-EFFICIENT implementation, see below. */ |
205 | { |
206 | int num_visited_nodes; |
207 | int i; |
208 | static boolean *node_in_neighborhood = NULL; |
209 | static int size = 0; |
210 | static int *index; |
211 | Queue Q; |
212 | heap H; |
213 | int closestVertex, neighbor; |
214 | DistType closestDist; |
215 | int num_found = 0; |
216 | |
217 | /* first, perform BFS to find the nodes in the region */ |
218 | mkQueue(&Q, n); |
219 | /* remember that dist should be init. with -1's */ |
220 | for (i = 0; i < n; i++) { |
221 | dist[i] = -1; /* far, TOO COSTLY (O(n))! */ |
222 | } |
223 | num_visited_nodes = |
224 | bfs_bounded(vertex, graph, n, dist, &Q, bound, visited_nodes); |
225 | if (size < n) { |
226 | node_in_neighborhood = |
227 | (boolean *) realloc(node_in_neighborhood, n * sizeof(boolean)); |
228 | for (i = size; i < n; i++) { |
229 | node_in_neighborhood[i] = FALSE; |
230 | } |
231 | size = n; |
232 | } |
233 | for (i = 0; i < num_visited_nodes; i++) { |
234 | node_in_neighborhood[visited_nodes[i]] = TRUE; |
235 | } |
236 | |
237 | |
238 | #ifdef OBSOLETE |
239 | mkHeap(&H, n); |
240 | #endif |
241 | index = (int *) realloc(index, n * sizeof(int)); |
242 | |
243 | /* initial distances with edge weights: */ |
244 | for (i = 0; i < n; i++) /* far, TOO COSTLY (O(n))! */ |
245 | dist[i] = (DistType) MAX_DIST; |
246 | dist[vertex] = 0; |
247 | for (i = 1; i < graph[vertex].nedges; i++) |
248 | dist[graph[vertex].edges[i]] = (DistType) graph[vertex].ewgts[i]; |
249 | |
250 | /* again, TOO COSTLY (O(n)) to put all nodes in heap! */ |
251 | initHeap(&H, vertex, index, dist, n); |
252 | |
253 | while (num_found < num_visited_nodes |
254 | && extractMax(&H, &closestVertex, index, dist)) { |
255 | if (node_in_neighborhood[closestVertex]) { |
256 | num_found++; |
257 | } |
258 | closestDist = dist[closestVertex]; |
259 | if (closestDist == MAX_DIST) |
260 | break; |
261 | for (i = 1; i < graph[closestVertex].nedges; i++) { |
262 | neighbor = graph[closestVertex].edges[i]; |
263 | increaseKey(&H, neighbor, |
264 | closestDist + |
265 | (DistType) graph[closestVertex].ewgts[i], index, |
266 | dist); |
267 | } |
268 | } |
269 | |
270 | /* restore initial false-status of 'node_in_neighborhood' */ |
271 | for (i = 0; i < num_visited_nodes; i++) { |
272 | node_in_neighborhood[visited_nodes[i]] = FALSE; |
273 | } |
274 | freeHeap(&H); |
275 | freeQueue(&Q); |
276 | return num_visited_nodes; |
277 | } |
278 | |
279 | static void heapify_f(heap * h, int i, int index[], float dist[]) |
280 | { |
281 | int l, r, largest; |
282 | while (1) { |
283 | l = left(i); |
284 | r = right(i); |
285 | if (insideHeap(h, l) && greaterPriority(h, l, i, dist)) |
286 | largest = l; |
287 | else |
288 | largest = i; |
289 | if (insideHeap(h, r) && greaterPriority(h, r, largest, dist)) |
290 | largest = r; |
291 | |
292 | if (largest == i) |
293 | break; |
294 | |
295 | exchange(h, largest, i, index); |
296 | i = largest; |
297 | } |
298 | } |
299 | |
300 | static void |
301 | initHeap_f(heap * h, int startVertex, int index[], float dist[], int n) |
302 | { |
303 | int i, count; |
304 | int j; /* We cannot use an unsigned value in this loop */ |
305 | h->data = N_GNEW(n - 1, int); |
306 | h->heapSize = n - 1; |
307 | |
308 | for (count = 0, i = 0; i < n; i++) |
309 | if (i != startVertex) { |
310 | h->data[count] = i; |
311 | index[i] = count; |
312 | count++; |
313 | } |
314 | |
315 | for (j = (n - 1) / 2; j >= 0; j--) |
316 | heapify_f(h, j, index, dist); |
317 | } |
318 | |
319 | static boolean (heap * h, int *max, int index[], float dist[]) |
320 | { |
321 | if (h->heapSize == 0) |
322 | return FALSE; |
323 | |
324 | *max = h->data[0]; |
325 | h->data[0] = h->data[h->heapSize - 1]; |
326 | index[h->data[0]] = 0; |
327 | h->heapSize--; |
328 | heapify_f(h, 0, index, dist); |
329 | |
330 | return TRUE; |
331 | } |
332 | |
333 | static void |
334 | increaseKey_f(heap * h, int increasedVertex, float newDist, int index[], |
335 | float dist[]) |
336 | { |
337 | int placeInHeap; |
338 | int i; |
339 | |
340 | if (dist[increasedVertex] <= newDist) |
341 | return; |
342 | |
343 | placeInHeap = index[increasedVertex]; |
344 | |
345 | dist[increasedVertex] = newDist; |
346 | |
347 | i = placeInHeap; |
348 | while (i > 0 && dist[h->data[parent(i)]] > newDist) { /* can write here: greaterPriority(i,parent(i),dist) */ |
349 | assign(h, i, parent(i), index); |
350 | i = parent(i); |
351 | } |
352 | h->data[i] = increasedVertex; |
353 | index[increasedVertex] = i; |
354 | } |
355 | |
356 | /* dijkstra_f: |
357 | * Weighted shortest paths from vertex. |
358 | * Assume graph is connected. |
359 | */ |
360 | void dijkstra_f(int vertex, vtx_data * graph, int n, float *dist) |
361 | { |
362 | int i; |
363 | heap H; |
364 | int closestVertex = 0, neighbor; |
365 | float closestDist; |
366 | int *index; |
367 | |
368 | #ifdef OBSOLETE |
369 | mkHeap(&H, n); |
370 | #endif |
371 | index = N_GNEW(n, int); |
372 | |
373 | /* initial distances with edge weights: */ |
374 | for (i = 0; i < n; i++) |
375 | dist[i] = MAXFLOAT; |
376 | dist[vertex] = 0; |
377 | for (i = 1; i < graph[vertex].nedges; i++) |
378 | dist[graph[vertex].edges[i]] = graph[vertex].ewgts[i]; |
379 | |
380 | initHeap_f(&H, vertex, index, dist, n); |
381 | |
382 | while (extractMax_f(&H, &closestVertex, index, dist)) { |
383 | closestDist = dist[closestVertex]; |
384 | if (closestDist == MAXFLOAT) |
385 | break; |
386 | for (i = 1; i < graph[closestVertex].nedges; i++) { |
387 | neighbor = graph[closestVertex].edges[i]; |
388 | increaseKey_f(&H, neighbor, |
389 | closestDist + graph[closestVertex].ewgts[i], |
390 | index, dist); |
391 | } |
392 | } |
393 | |
394 | freeHeap(&H); |
395 | free(index); |
396 | } |
397 | |