Some enhancements for v0.4

CITATION.bib

@@ -2,7 +2,7 @@ @misc{ConstrainedShortestPaths.jl
 	author  = {Léo Baty and contributors},
 	title   = {ConstrainedShortestPaths.jl},
 	url     = {https://github.com/BatyLeo/ConstrainedShortestPaths.jl},
-	version = {v0.3.0},
+	version = {v0.4.0},
 	year    = {2024},
 	month   = {8}
 }

Project.toml

@@ -1,7 +1,7 @@
 name = "ConstrainedShortestPaths"
 uuid = "b3798467-87dc-4d99-943d-35a1bd39e395"
 authors = ["Léo Baty and contributors"]
-version = "0.3.0"
+version = "0.4.0"
 
 [deps]
 DataStructures = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8"
@@ -22,6 +22,7 @@ julia = "1.10"
 
 [extras]
 Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595"
+Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
 GLPK = "60bf3e95-4087-53dc-ae20-288a0d20c6a6"
 Graphs = "86223c79-3864-5bf0-83f7-82e725a168b6"
 JET = "c3a54625-cd67-489e-a8e7-0a5a0ff4e31b"
@@ -35,4 +36,4 @@ Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 UnicodePlots = "b8865327-cd53-5732-bb35-84acbb429228"
 
 [targets]
-test = ["Aqua", "GLPK", "Graphs", "JET", "JuMP", "JuliaFormatter", "Plots", "Random", "SparseArrays", "StableRNGs", "Test", "UnicodePlots"]
+test = ["Aqua", "Documenter", "GLPK", "Graphs", "JET", "JuMP", "JuliaFormatter", "Plots", "Random", "SparseArrays", "StableRNGs", "Test", "UnicodePlots"]

README.md

@@ -8,7 +8,10 @@
 
 ## Overview
 
-This package implements algorithms for solving **Generalized Constrained Shortest Paths** problems. It implements a simplified version of the framework from [Parmentier 2017](https://arxiv.org/abs/1504.07880), restricted to acyclic directed graphs.
+This package implements algorithms for solving (resource) **Constrained Shortest Paths** problems. 
+It implements a generalized A star algorithm with label dominance and optional bounding.
+It is currently restricted to acyclic directed graphs.
+Reference: [https://arxiv.org/abs/1504.07880](https://arxiv.org/abs/1504.07880).
 
 Let $D=(V, A)$ an **acyclic directed graph**, $o, d\in V$ **origin** and **destination** vertices, $c$ a **cost function**, and $\mathcal{P} \subset \mathcal{P}_{od}$ a subset of $o-d$ paths in $G$. This package can compute the corresponding **constrained shortest path**:
 
@@ -23,8 +26,8 @@ See the [documentation](https://batyleo.github.io/ConstrainedShortestPaths.jl) f
 
 ## Installation
 
-To install this package, open a julia REPL and run the following command in pkg mode:
+To install this package, open a julia REPL and run the following command:
 
 ```bash
-add ConstrainedShortestPaths
+]add ConstrainedShortestPaths
 ```

docs/src/literate/custom.jl

@@ -115,6 +115,15 @@ fb = [BackwardExpansionFunction(d[i, j], w[i, j]) for (i, j) in zip(If, Jf)]
 FF = sparse(If, Jf, ff);
 FB = sparse(If, Jf, fb);
 
-instance = CSPInstance(graph, 1, nb_vertices, resource, resource, Cost(), FF, FB)
+instance = CSPInstance(;
+    graph,
+    origin_vertex=1,
+    destination_vertex=nb_vertices,
+    origin_forward_resource=resource,
+    destination_backward_resource=resource,
+    cost_function=Cost(),
+    forward_functions=FF,
+    backward_functions=FB,
+)
 (; p_star, c_star) = generalized_constrained_shortest_path(instance; W=W)
 @info "Result" c_star p_star

src/ConstrainedShortestPaths.jl

@@ -3,12 +3,22 @@ module ConstrainedShortestPaths
 using DataStructures: PriorityQueue, enqueue!, dequeue!, isempty
 using DocStringExtensions: TYPEDEF, TYPEDFIELDS, TYPEDSIGNATURES
 using Graphs:
-    AbstractGraph, is_directed, is_cyclic, nv, weights, src, dst, edges, outneighbors
+    AbstractGraph,
+    is_directed,
+    is_cyclic,
+    nv,
+    weights,
+    src,
+    dst,
+    edges,
+    outneighbors,
+    induced_subgraph
 using PiecewiseLinearFunctions: PiecewiseLinearFunction
 using SparseArrays: sparse
 using Statistics: mean
 
 include("utils/utils.jl")
+include("interface.jl")
 include("algorithms.jl")
 include("examples/basic_shortest_path.jl")
 include("examples/resource_shortest_path.jl")

src/algorithms.jl

@@ -1,75 +1,22 @@
-"""
-$TYPEDEF
-
-# Fields
-$TYPEDFIELDS
-"""
-struct CSPInstance{T,G<:AbstractGraph{T},FR,BR,C,FF<:AbstractMatrix,BF<:AbstractMatrix}
-    "acyclic digraph in which to compute the shortest path"
-    graph::G
-    "origin vertex of path"
-    origin_vertex::T
-    "destination vertex of path"
-    destination_vertex::T
-    "forward resource at the origin vertex"
-    origin_forward_resource::FR
-    "backward resource at the destination vertex"
-    destination_backward_resource::BR
-    "cost function"
-    cost_function::C
-    "forward functions along edges"
-    forward_functions::FF
-    "backward functions along edges"
-    backward_functions::BF
-end
-
-"""
-$TYPEDSIGNATURES
-
-Constructor for [`CSPInstance`](@ref).
-"""
-function CSPInstance(;
-    graph,
-    origin_vertex,
-    destination_vertex,
-    origin_forward_resource,
-    destination_backward_resource,
-    cost_function,
-    forward_functions,
-    backward_functions,
-)
-    @assert is_directed(graph) "`graph` must be a directed graph"
-    @assert !is_cyclic(graph) "`graph` must be acyclic"
-    return CSPInstance(
-        graph,
-        origin_vertex,
-        destination_vertex,
-        origin_forward_resource,
-        destination_backward_resource,
-        cost_function,
-        forward_functions,
-        backward_functions,
-    )
-end
-
 """
 $TYPEDSIGNATURES
 
 Compute backward bounds of instance (see [Computing bounds](@ref)).
 """
-function compute_bounds(instance::CSPInstance{T,G}; kwargs...) where {T,G<:AbstractGraph{T}}
-    (; graph, origin_vertex, destination_vertex) = instance
+function compute_bounds(instance::CSPInstance; kwargs...)
+    (; graph, destination_vertex, topological_ordering, is_useful) = instance
 
-    vertices_order = topological_order(graph, origin_vertex, destination_vertex)
+    vertices_order = topological_ordering
+    @assert vertices_order[1] == destination_vertex
 
     bounds = Dict{Int,typeof(instance.destination_backward_resource)}()
-    # bounds = Vector{typeof(instance.destination_backward_resource)}(undef, nb_vertices)
+    # bounds = Vector{typeof(instance.destination_backward_resource)}(undef, nv(graph))
     bounds[destination_vertex] = instance.destination_backward_resource
 
     for vertex in vertices_order[2:end]
         vector = [
             instance.backward_functions[vertex, neighbor](bounds[neighbor]; kwargs...) for
-            neighbor in outneighbors(graph, vertex) if haskey(bounds, neighbor)
+            neighbor in outneighbors(graph, vertex) if is_useful[neighbor]
         ]
         bounds[vertex] = minimum(vector)
     end
@@ -83,10 +30,8 @@ $TYPEDSIGNATURES
 Perform generalized A star algorithm on instnace using bounds
 (see [Generalized `A^\\star`](@ref)).
 """
-function generalized_a_star(
-    instance::CSPInstance{T,G}, bounds::AbstractDict; kwargs...
-) where {T,G<:AbstractGraph{T}}
-    (; graph, origin_vertex, destination_vertex) = instance
+function generalized_a_star(instance::CSPInstance, bounds; kwargs...)
+    (; graph, origin_vertex, destination_vertex, is_useful) = instance
     nb_vertices = nv(graph)
 
     empty_path = [origin_vertex]
@@ -107,7 +52,7 @@ function generalized_a_star(
         p = dequeue!(L)
         v = p[end]
         for w in outneighbors(graph, v)
-            if !haskey(bounds, w)
+            if !is_useful[w]
                 continue
             end
             q = copy(p)
@@ -137,12 +82,65 @@ end
 """
 $TYPEDSIGNATURES
 
+Label dominance dynamic programming without bounding.
+"""
+function generalized_a_star(instance::ForwardCSPInstance; kwargs...)
+    (; graph, origin_vertex, destination_vertex, is_useful) = instance
+    nb_vertices = nv(graph)
+
+    empty_path = [origin_vertex]
+
+    forward_resources = Dict(empty_path => instance.origin_forward_resource)
+    L = PriorityQueue{Vector{Int},Float64}(
+        empty_path => instance.cost_function(forward_resources[empty_path])
+    )
+
+    forward_type = typeof(forward_resources[empty_path])
+    M = [forward_type[] for _ in 1:nb_vertices]
+    push!(M[origin_vertex], forward_resources[empty_path])
+    c_star = Inf
+    p_star = [origin_vertex]
+
+    while !isempty(L)
+        p = dequeue!(L)
+        v = p[end]
+        for w in outneighbors(graph, v)
+            if !is_useful[w]
+                continue
+            end
+            q = copy(p)
+            push!(q, w)
+            rp = forward_resources[p]
+            rq, is_feasible = instance.forward_functions[v, w](rp; kwargs...)
+            if !is_feasible
+                continue
+            end
+            forward_resources[q] = rq
+            c = instance.cost_function(rq) # compute partial cost
+            if w == destination_vertex # if destination is reached
+                if c < c_star
+                    c_star = c
+                    p_star = copy(q)
+                end
+            elseif !is_dominated(rq, M[w]) # else add path to queue if not dominated
+                remove_dominated!(M[w], rq)
+                push!(M[w], rq)
+                enqueue!(L, q => c)
+            end
+        end
+    end
+    return (; p_star, c_star)
+end
+
+"""
+$TYPEDSIGNATURES
+
 Compute all paths below threshold.
 """
 function generalized_a_star_with_threshold(
-    instance::CSPInstance{T,G}, bounds::AbstractDict, threshold::Float64; kwargs...
-) where {T,G<:AbstractGraph}
-    (; graph, origin_vertex, destination_vertex) = instance
+    instance::CSPInstance, bounds, threshold::Float64; kwargs...
+)
+    (; graph, origin_vertex, destination_vertex, is_useful) = instance
 
     empty_path = [origin_vertex]
 
@@ -159,6 +157,9 @@ function generalized_a_star_with_threshold(
         p = dequeue!(L)
         v = p[end]
         for w in outneighbors(graph, v)
+            if !is_useful[w]
+                continue
+            end
             q = copy(p)
             push!(q, w)
             rp = forward_resources[p]
@@ -187,13 +188,15 @@ $TYPEDSIGNATURES
 
 Compute the shortest path of `instance`.
 """
-function generalized_constrained_shortest_path(
-    instance::CSPInstance{T,G}; kwargs...
-) where {T,G<:AbstractGraph{T}}
+function generalized_constrained_shortest_path(instance::CSPInstance; kwargs...)
     bounds = compute_bounds(instance; kwargs...)
     return generalized_a_star(instance, bounds; kwargs...)
 end
 
+function generalized_constrained_shortest_path(instance::ForwardCSPInstance; kwargs...)
+    return generalized_a_star(instance; kwargs...)
+end
+
 """
 $TYPEDSIGNATURES
 

src/examples/basic_shortest_path.jl

@@ -16,10 +16,6 @@ function (f::BSPBackwardExtensionFunction)(q::BSPResource)
     return f.c + q
 end
 
-function BSP_cost(qf::BSPResource, qb::BSPResource)
-    return qf + qb
-end
-
 function remove_dominated!(Mw::Vector{BSPResource}, rq::BSPResource)
     empty!(Mw)
     return Mw = push!(Mw, rq)
@@ -42,7 +38,11 @@ Compute shortest path between vertices `s` and `t` of graph `graph`.
 - `c_star::Float64`: length of path `p_star`.
 """
 function basic_shortest_path(
-    graph::AbstractGraph{T}, s::T, t::T, distmx::AbstractMatrix=weights(graph)
+    graph::AbstractGraph{T},
+    s::T,
+    t::T,
+    distmx::AbstractMatrix=weights(graph);
+    bounding=true,
 ) where {T}
     # origin forward resource and backward forward resource set to 0
     resource = 0.0
@@ -51,19 +51,31 @@ function basic_shortest_path(
     If = [src(e) for e in edges(graph)]
     Jf = [dst(e) for e in edges(graph)]
     ff = [BSPForwardExtensionFunction(distmx[i, j]) for (i, j) in zip(If, Jf)]
-    fb = [BSPBackwardExtensionFunction(distmx[i, j]) for (i, j) in zip(If, Jf)]
     FF = sparse(If, Jf, ff)
-    FB = sparse(If, Jf, fb)
 
-    instance = CSPInstance(;
-        graph,
-        origin_vertex=s,
-        destination_vertex=t,
-        origin_forward_resource=resource,
-        destination_backward_resource=resource,
-        cost_function=BSP_cost,
-        forward_functions=FF,
-        backward_functions=FB,
-    )
+    instance = if bounding
+        fb = [BSPBackwardExtensionFunction(distmx[i, j]) for (i, j) in zip(If, Jf)]
+        FB = sparse(If, Jf, fb)
+
+        CSPInstance(;
+            graph,
+            origin_vertex=s,
+            destination_vertex=t,
+            origin_forward_resource=resource,
+            destination_backward_resource=resource,
+            cost_function=Base.:+,
+            forward_functions=FF,
+            backward_functions=FB,
+        )
+    else
+        CSPInstance(;
+            graph,
+            origin_vertex=s,
+            destination_vertex=t,
+            origin_forward_resource=resource,
+            cost_function=identity,
+            forward_functions=FF,
+        )
+    end
     return generalized_constrained_shortest_path(instance)
 end