two sided dijkstra

Project.toml

@@ -7,13 +7,15 @@ ArnoldiMethod = "ec485272-7323-5ecc-a04f-4719b315124d"
 Compat = "34da2185-b29b-5c13-b0c7-acf172513d20"
 DataStructures = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8"
 Distributed = "8ba89e20-285c-5b6f-9357-94700520ee1b"
+Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
 Inflate = "d25df0c9-e2be-5dd7-82c8-3ad0b3e990b9"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 SharedArrays = "1a1011a3-84de-559e-8e89-a11a2f7dc383"
 SimpleTraits = "699a6c99-e7fa-54fc-8d76-47d257e15c1d"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
+StatsPlots = "f3b207a7-027a-5e70-b257-86293d7955fd"
 
 [compat]
 Aqua = "0.6"

src/Graphs.jl

@@ -248,6 +248,7 @@ export
     transitivereduction,
     yen_k_shortest_paths,
     desopo_pape_shortest_paths,
+    bidijkstra_shortest_path,
 
     # centrality
     betweenness_centrality,

src/shortestpaths/dijkstra.jl

@@ -169,3 +169,115 @@ function dijkstra_shortest_paths(
         g, [src;], distmx; allpaths=allpaths, trackvertices=trackvertices, maxdist=maxdist
     )
 end
+
+function relax(u::U, 
+               v::U, 
+               distmx::AbstractMatrix{T}, 
+               dists::Vector{T}, 
+               parents::Vector{U}, 
+               visited::Vector{Bool}, 
+               Q::PriorityQueue{U,T}
+) where {T<:Real} where {U<:Integer}
+     alt = dists[u] + distmx[u, v]
+
+     if !visited[v]
+         visited[v] = true
+         dists[v] = alt
+         parents[v] = u
+         Q[v] = alt
+     elseif alt < dists[v]
+         dists[v] = alt
+         parents[v] = u
+         Q[v] = alt
+     end
+end
+
+"""
+    bidijkstra_shortest_paths(g, src, dst, distmx=weights(g));
+
+Perform [Bidirectional Dijkstra's algorithm](https://www.homepages.ucl.ac.uk/~ucahmto/math/2020/05/30/bidirectional-dijkstra.html)
+on a graph, computing the shortest path between `src` and `dst`.
+
+# Examples
+```jldoctest
+julia> using Graphs
+
+julia> bidijkstra_shortest_path(cycle_graph(5), 1, 4);
+3-element Vector{Int64}:
+ 1
+ 5
+ 4
+
+julia> bidijkstra_shortest_path(path_graph(5), 1, 4);
+4-element Vector{Int64}:
+ 1
+ 2
+ 3
+ 4
+```
+"""
+function bidijkstra_shortest_path(
+    g::AbstractGraph,
+    src::U,
+    dst::U,
+    distmx::AbstractMatrix{T}=weights(g)
+) where {T<:Real} where {U<:Integer}
+    if src == dst
+        return Int64[]
+    end
+    # keep weight of the best seen path and the midpoint vertex
+    mu, mid_v = typemax(T), -1
+
+    nvg = nv(g)
+    dists_f, dists_b= fill(typemax(T), nvg), fill(typemax(T), nvg)
+
+    parents_f, parents_b= zeros(U, nvg), zeros(U, nvg)
+
+    visited_f, visited_b = zeros(Bool, nvg),zeros(Bool, nvg)
+
+    preds_f, preds_b = fill(Vector{U}(), nvg), fill(Vector{U}(), nvg)
+
+    Qf, Qb = PriorityQueue{U,T}(), PriorityQueue{U,T}()
+
+    dists_f[src] = zero(T)
+    visited_f[src] = true
+    Qf[src] = zero(T)
+
+    dists_b[dst] = zero(T)
+    visited_b[dst] = true
+    Qb[dst] = zero(T)
+
+    while !isempty(Qf) && !isempty(Qb)
+        uf, ub = dequeue!(Qf), dequeue!(Qb)
+
+        for v in outneighbors(g, uf)
+            relax(uf, v, distmx, dists_f, parents_f, visited_f, Qf)
+            if visited_b[v] && (dists_f[uf]+distmx[uf,v]+dists_b[v]) < mu
+                # we have found an edge between the forward and backward exploration
+                mu = dists_f[uf]+distmx[uf,v]+dists_b[v]
+                mid_v = v
+            end
+        end
+
+        for v in inneighbors(g, ub)
+            relax(ub, v, distmx, dists_b, parents_b, visited_b, Qb)
+            if visited_f[v] && (dists_f[v]+distmx[v,ub]+dists_b[ub]) < mu
+                # we have found an edge between the forward and backward exploration
+                mu = dists_f[v]+distmx[v,ub]+dists_b[ub]
+                mid_v = v
+            end
+        end
+        if dists_f[uf]+dists_b[ub] >= mu
+            break
+        end
+    end
+    if mid_v == -1
+        # no path exists between source and destination
+        return Int64[]
+    end
+    ds_f = DijkstraState{T,U}(parents_f, dists_f, preds_f, zeros(nvg), Vector{U}())
+    ds_b = DijkstraState{T,U}(parents_b, dists_b, preds_b, zeros(nvg), Vector{U}())
+    path = vcat(enumerate_paths(ds_f, mid_v), reverse(enumerate_paths(ds_b, mid_v)[1:end-1]))
+    return path
+end
+

test/shortestpaths/dijkstra.jl

@@ -111,4 +111,20 @@
             ds = @inferred(dijkstra_shortest_paths(g, 3, m;maxdist=3.0))
             @test ds.dists == [2, 1, 0, Inf, Inf, 3]
         end
+
+    # bidijkstra_shortest_path
+    g4 = path_digraph(5)
+    d1 = float([0 1 2 3 4; 5 0 6 7 8; 9 10 0 11 12; 13 14 15 0 16; 17 18 19 20 0])
+    d2 = sparse(float([0 1 2 3 4; 5 0 6 7 8; 9 10 0 11 12; 13 14 15 0 16; 17 18 19 20 0]))
+
+    for g in testdigraphs(g4)
+        x = @inferred(dijkstra_shortest_paths(g, 2, d1))
+        p = enumerate_paths(x, 4)
+        y = @inferred(bidijkstra_shortest_path(g, 2, 4, d1))
+        z = @inferred(bidijkstra_shortest_path(g, 2, 4, d2))
+
+        @test p == y == z
+    end
+
+
 end