implement minimum cost flow

Author: jbytecode

src/OperationsResearchModels.jl

@@ -59,6 +59,7 @@ import .ShortestPath: ShortestPathResult, ShortestPathProblem
 
 import .Network: Connection, nodes
 import .MaximumFlow: MaximumFlowResult, MaximumFlowProblem
+import .MaximumFlow: MinimumCostFlowProblem, MinimumCostFlowResult
 import .Assignment: AssignmentProblem, AssignmentResult
 import .Game: game, GameResult, game_solver
 import .MinimumSpanningTree: hasloop, MstResult, MstProblem
@@ -73,8 +74,8 @@ import .TravelingSalesman: TravelinSalesmenResult, travelingsalesman
 import .Simplex: SimplexProblem, simplexiterations, createsimplexproblem, gaussjordan
 
 export TransportationProblem, TransportationResult, balance, isbalanced, northwestcorner, leastcost
-export Connection, ShortestPathResult, MaximumFlowResult, nodes
-export ShortestPathProblem, MaximumFlowProblem
+export Connection, ShortestPathResult, MaximumFlowResult, MinimumCostFlowResult, nodes
+export ShortestPathProblem, MaximumFlowProblem, MinimumCostFlowProblem
 export AssignmentProblem, AssignmentResult
 export game, GameResult, game_solver
 export hasloop, MstResult, MstProblem

src/maximumflow.jl

@@ -6,8 +6,8 @@ using JuMP, HiGHS
 import ..OperationsResearchModels: solve
 
 
-export MaximumFlowProblem
-export MaximumFlowResult
+export MaximumFlowProblem, MaximumFlowResult
+export MinimumCostFlowProblem, MinimumCostFlowResult
 
 
 struct MaximumFlowResult
@@ -19,6 +19,66 @@ struct MaximumFlowProblem
     connections::Vector{Connection}
 end
 
+struct MinimumCostFlowProblem
+    connections::Vector{Connection}
+    costs::Vector{Connection}
+end
+
+struct MinimumCostFlowResult 
+    path::Vector{Connection}
+    cost::Float64
+end 
+
+
+function leftexpressions(x::Matrix{JuMP.VariableRef}, node::Int64, nodes::Vector{Connection}, model)
+    lst = []
+    for conn in nodes
+        if conn.to == node
+            push!(lst, conn)
+        end
+    end
+    if length(lst) == 0
+        return :f
+    end
+    expr = @expression(model, 0)
+    for i = eachindex(lst)
+        expr += x[lst[i].from, lst[i].to]
+    end
+    return expr
+end
+
+function rightexpressions(x::Matrix{JuMP.VariableRef}, node::Int64, nodes::Vector{Connection}, model)
+    lst = []
+    for conn in nodes
+        if conn.from == node
+            push!(lst, conn)
+        end
+    end
+    if length(lst) == 0
+        return :f
+    end
+    expr = @expression(model, 0)
+    for i = eachindex(lst)
+        expr += x[lst[i].from, lst[i].to]
+    end
+    return expr
+end
+
+
+function hassameorder(a::Vector{Connection}, b::Vector{Connection})::Bool
+
+    if length(a) != length(b)
+        return false
+    end
+
+    for i = 1:length(a)
+        if a[i].from != b[i].from || a[i].to != b[i].to
+            return false
+        end
+    end
+
+    return true
+end
 
 """
 
@@ -65,45 +125,10 @@ julia> result.flow
 7.0
 ```
 """
-function solve(problem::MaximumFlowProblem)
+function solve(problem::MaximumFlowProblem)::MaximumFlowResult
 
     cns = problem.connections
 
-    function leftexpressions(node::Int64, nodes::Vector{Connection}, model)
-        lst = []
-        for conn in nodes
-            if conn.to == node
-                push!(lst, conn)
-            end
-        end
-        if length(lst) == 0
-            return :f
-        end
-        expr = @expression(model, 0)
-        for i = eachindex(lst)
-            expr += x[lst[i].from, lst[i].to]
-        end
-        return expr
-    end
-
-    function rightexpressions(node::Int64, nodes::Vector{Connection}, model)
-        lst = []
-        for conn in nodes
-            if conn.from == node
-                push!(lst, conn)
-            end
-        end
-        if length(lst) == 0
-            return :f
-        end
-        expr = @expression(model, 0)
-        for i = eachindex(lst)
-            expr += x[lst[i].from, lst[i].to]
-        end
-        return expr
-    end
-
-
     model = Model(HiGHS.Optimizer)
     MOI.set(model, MOI.Silent(), true)
 
@@ -122,8 +147,8 @@ function solve(problem::MaximumFlowProblem)
 
     # Constraints 
     for nextnode in mynodes
-        leftexpr = leftexpressions(nextnode, cns, model)
-        rightexpr = rightexpressions(nextnode, cns, model)
+        leftexpr = leftexpressions(x, nextnode, cns, model)
+        rightexpr = rightexpressions(x, nextnode, cns, model)
         if leftexpr == :f
             @constraint(model, rightexpr == f)
         elseif rightexpr == :f
@@ -157,4 +182,101 @@ end
 
 
 
+"""
+    solve(problem, flow)
+
+# Description
+
+This function solves the Minimum Cost Flow problem given a flow value.
+
+# Arguments
+
+- `problem::MinimumCostFlowProblem`: The problem in type of MinimumCostFlowProblem
+- `flow::Float64`: The flow value to be used in the problem.
+
+
+"""
+function solve(problem::MinimumCostFlowProblem, flow::Float64)::MinimumCostFlowResult
+    cns = problem.connections
+    costs = problem.costs
+
+    if !hassameorder(cns, costs)
+        throw(ArgumentError("Connections and costs must have the same order."))
+    end
+
+    model = Model(HiGHS.Optimizer)
+    MOI.set(model, MOI.Silent(), true)
+
+    mynodes = nodes(cns)
+    n = length(mynodes)
+
+    startnode = start(cns)
+    finishnode = finish(cns)
+
+    # Variables 
+    @variable(model, x[1:n, 1:n] .>= 0)
+
+    # Objective Function
+    @objective(model, Max, sum(x[conn.from, conn.to] * conn.value for conn in costs))
+
+    # Constraints 
+    for nextnode in mynodes
+        leftexpr = leftexpressions(x, nextnode, cns, model)
+        rightexpr = rightexpressions(x, nextnode, cns, model)
+        if leftexpr == :f
+            @constraint(model, rightexpr == flow)
+        elseif rightexpr == :f
+            @constraint(model, leftexpr == flow)
+        else
+            @constraint(model, leftexpr - rightexpr == 0)
+        end
+    end
+
+    for nd in cns
+        @constraint(model, x[nd.from, nd.to] <= nd.value)
+    end
+
+    optimize!(model)
+
+    xs = value.(x)
+    cost = JuMP.objective_value(model)
+    solutionnodes = []
+    for i = 1:n
+        for j = 1:n
+            if xs[i, j] > 0
+                push!(solutionnodes, Connection(i, j, xs[i, j]))
+            end
+        end
+    end
+
+    return MinimumCostFlowResult(solutionnodes, cost)
+end 
+
+
+"""
+    solve(problem)
+
+# Description
+
+This function solves the Minimum Cost Flow problem by first solving the Maximum Flow problem and 
+then using the flow value to solve the Minimum Cost Flow problem.
+
+# Arguments
+
+- `problem::MinimumCostFlowProblem`: The problem in type of MinimumCostFlowProblem
+
+"""
+function solve(problem::MinimumCostFlowProblem)::MinimumCostFlowResult
+
+    maximumflowproblem = MaximumFlowProblem(problem.connections)
+
+    maximumflowresult = solve(maximumflowproblem)
+
+    f = maximumflowresult.flow
+
+    result::MinimumCostFlowResult = solve(problem, f)
+    
+    return result
+end 
+
 end # end of module

test/runtests.jl

@@ -12,6 +12,7 @@ include("testjohnsons.jl")
 include("testutility.jl")
 include("testshortestpath.jl")
 include("testmaximumflow.jl")
+include("testminimumcostflow.jl")
 include("testgame.jl")
 include("testmst.jl")
 include("testpmedian.jl")

test/testmaximumflow.jl

@@ -46,3 +46,5 @@
         end
     end
 end
+
+

test/testminimumcostflow.jl

@@ -0,0 +1,50 @@
+@testset "Minimum Cost Flow" verbose = true begin
+    @testset "Basic Minimum Cost Flow example" begin
+        conns = [
+            Connection(1, 2, 20),
+            Connection(1, 3, 30),
+            Connection(2, 4, 10),
+            Connection(3, 4, 5)
+        ]
+
+        costs = [
+            Connection(1, 2, 1),
+            Connection(1, 3, 2),
+            Connection(2, 4, 3),
+            Connection(3, 4, 3)
+        ]
+
+        problem = MinimumCostFlowProblem(conns, costs)
+
+        result = solve(problem)
+
+        @test result isa MinimumCostFlowResult
+        @test result.cost == 65.0
+    end
+
+    @testset "Basic Minimum Cost Flow example -2" begin
+        conns = [
+            Connection(1, 2, 20),
+            Connection(1, 3, 30),
+            Connection(2, 4, 10),
+            Connection(3, 4, 5),
+            Connection(3, 2, 25)
+        ]
+
+        costs = [
+            Connection(1, 2, 1),
+            Connection(1, 3, 2),
+            Connection(2, 4, 3),
+            Connection(3, 4, 3),
+            Connection(3, 2, 2)
+        ]
+
+        problem = MinimumCostFlowProblem(conns, costs)
+
+        result = solve(problem)
+
+        @test result isa MinimumCostFlowResult
+        @test result.cost == 95.0
+
+    end
+end