using JuMP
using HiGHS
model = Model(HiGHS.Optimizer);

@variable(model, x₁≥0)
@variable(model, x₂≥0)
@objective(model, Min, 50x₁+80x₂)
@constraint(model, c1, 3x₁≥6)
@constraint(model, c2, 2x₁+4x₂≥10)
@constraint(model, c3, 2x₁+5x₂≥8)
model
print(model)
optimize!(model)
@show termination_status(model)
@show objective_value(model)
@show value(x₁), value(x₂)

model = Model(HiGHS.Optimizer);
@variable(model, x[1:2] .≥ 0)
A = [3 0;2 4;2 5]
b = [6,10,8]
c = [50,80]
@objective(model, Min, c'*x)
@constraint(model, (A*x) .≥ b)
print(model)
optimize!(model)
@show termination_status(model)
@show objective_value(model)
@show value.(x)

julia> model = Model(HiGHS.Optimizer)
A JuMP Model
Feasibility problem with:
Variables: 0
Model mode: AUTOMATIC
CachingOptimizer state: EMPTY_OPTIMIZER
Solver name: HiGHS

julia> @variable(model, x[i in 1:3, j in i:3])
JuMP.Containers.SparseAxisArray{VariableRef, 2, Tuple{Int64, Int64}} with 6 entries:
  [1, 1]  =  x[1,1]
  [1, 2]  =  x[1,2]
  [1, 3]  =  x[1,3]
  [2, 2]  =  x[2,2]
  [2, 3]  =  x[2,3]
  [3, 3]  =  x[3,3]

model = Model(HiGHS.Optimizer);
@variable(model, x[1:2] .≥ 0)
@objective(model, Max, x[1])
@constraint(model, x .≤ -10)
optimize!(model)
@show termination_status(model)
@show objective_value(model)
@show value.(x)

model = Model(HiGHS.Optimizer);
@variable(model, x[1:2] .≥ 0)
@objective(model, Max, 0)
@constraint(model, x .≤ 10)
optimize!(model)
@show termination_status(model)
@show objective_value(model)
@show value.(x)

model = Model(HiGHS.Optimizer);
@variable(model, x[1:2] .> 0)

(2023-12-26_julia-optimization) pkg> st JuMP HiGHS
Status `~/Git/julia-training/2023-12-26_julia-optimization/Project.toml`
  [87dc4568] HiGHS v1.7.5
  [4076af6c] JuMP v1.17.0

using JuMP
using HiGHS

model = Model(HiGHS.Optimizer)
@variable(model, x₁ ≥ 0)
@variable(model, x₂ ≥ 0)
@variable(model, x₃ ≥ 0)
@objective(model, Max, 3x₁+5x₂+4x₃)
@constraint(model, c1, 2x₁+x₂+3x₃ ≤ 6)
@constraint(model, c2, 3x₁-2x₂+x₃ ≤ 7)
@constraint(model, c3, x₂+2x₃ ≤ 5)
set_silent(model)
optimize!(model)
@show termination_status(model)
@show objective_value(model)
@show shadow_price(c1)
@show shadow_price(c2)
@show shadow_price(c3)
@show value(x₁), value(x₂), value(x₃)

using Graphs
G = DiGraph(5)
add_edge!(G,1,2)
add_edge!(G,2,3)
add_edge!(G,2,4)
add_edge!(G,3,4)
add_edge!(G,4,5)
nv(G)
ne(G)
outneighbors(G,2)

using GLMakie
using GraphMakie
using Graphs
G = DiGraph(5)
add_edge!(G,1,2)
add_edge!(G,2,3)
add_edge!(G,2,4)
add_edge!(G,3,4)
add_edge!(G,4,5)
nv(G)
ne(G)
outneighbors(G,2)
graphplot(G)

using JuMP
using HiGHS
using LinearAlgebra

n = 5
V = 1:n
s = 1
t = 2
G = [
    0 100  30   0   0
    0   0  20   0   0
    0   0   0  10  60
    0  15   0   0  50
    0   0   0   0   0
]
model = Model(HiGHS.Optimizer)
@variable(model, x[V, V], Bin)
@objective(model, Min, dot(G,x))
@constraint(model, [x[i,j] for i in V, j in V if G[i,j]==0] .== 0)
@constraint(model, [sum(x[i,j] for j in V if j≠i) - sum(x[j,i] for j in V if j≠i) for i in V if i ∉ (s,t)] .== 0)
@constraint(model, sum(x[s,j] for j in V if j≠s) - sum(x[j,s] for j in V if j≠s) .== 1)
@constraint(model, sum(x[t,j] for j in V if j≠t) - sum(x[j,t] for j in V if j≠t) .== -1)
set_silent(model)
optimize!(model)
@show termination_status(model)
@show objective_value(model)
for i in 1:n, j in 1:n
    v = value(x[i,j])
    if v ≠ 0
        println("x[$i,$j] = ", v)
    end
end

using SimpleWeightedGraphs
using JuMP
using HiGHS

n = 5
s = 1
t = 2
V = 1:n
g = SimpleWeightedDiGraph(n)
I = [  1,  1,  2,  3,  3,  4,  4]
J = [  2,  3,  3,  4,  5,  2,  5]
W = [100, 30, 20, 10, 60, 15, 50]
for (i,j,w) in zip(I,J,W)
    add_edge!(g,i,j,w)
end
model = Model(HiGHS.Optimizer)
@variable(model, x[i in 1:n, j in 1:n; has_edge(g, (i,j))], Bin)
@objective(model, Min, sum(v * x[i,j] for (i,j,v) in zip(I,J,W)))
@constraint(model, c_src, sum(x[s,i] for i in outneighbors(g, s)) == 1);
@constraint(model, c_dst, sum(x[i,t] for i in inneighbors(g, t)) == 1);
@constraint(model, c_flow[v in V; v≠s&&v≠t], sum(x[v,i] for i in outneighbors(g,v)) == sum(x[j,v] for j in inneighbors(g,v)))
set_silent(model)
optimize!(model)
@show termination_status(model)
@show objective_value(model)
for x in x
    if value(x) == 1
        println(x)
    end
end

julia> using Graphs, SimpleWeightedGraphs

julia> n = 5
5

julia> s = 1
1

julia> t = 2
2

julia> V = 1:n
1:5

julia> g = SimpleWeightedDiGraph(n)
{5, 0} directed simple Int64 graph with Float64 weights

julia> I = [  1,  1,  2,  3,  3,  4,  4]
7-element Vector{Int64}:
 1
 1
 2
 3
 3
 4
 4

julia> J = [  2,  3,  3,  4,  5,  2,  5]
7-element Vector{Int64}:
 2
 3
 3
 4
 5
 2
 5

julia> W = [100, 30, 20, 10, 60, 15, 50]
7-element Vector{Int64}:
 100
  30
  20
  10
  60
  15
  50

julia> for (i,j,w) in zip(I,J,W)
           add_edge!(g,i,j,w)
       end

julia> r = dijkstra_shortest_paths(g, [1], weights(g))
Graphs.DijkstraState{Float64, Int64}([0, 4, 1, 3, 3], [0.0, 55.0, 30.0, 40.0, 90.0], [Int64[], Int64[], Int64[], Int64[], Int64[]], [1.0, 1.0, 1.0, 1.0, 2.0], Int64[])

julia> r.dists
5-element Vector{Float64}:
  0.0
 55.0
 30.0
 40.0
 90.0

julia> r.parents
5-element Vector{Int64}:
 0
 4
 1
 3
 3

julia> dsts = [2,3,4,5]
4-element Vector{Int64}:
 2
 3
 4
 5

julia> spaths = enumerate_paths(r, dsts)
4-element Vector{Vector{Int64}}:
 [1, 3, 4, 2]
 [1, 3]
 [1, 3, 4]
 [1, 3, 5]

julia> for (i,s) in enumerate(dsts)
           println("1 から $s への最短路: $(spaths[i])")
       end
1 から 2 への最短路: [1, 3, 4, 2]
1 から 3 への最短路: [1, 3]
1 から 4 への最短路: [1, 3, 4]
1 から 5 への最短路: [1, 3, 5]

using Graphs
using SimpleWeightedGraphs
using GraphRecipes
using Plots
using Random

n = 6
g = SimpleWeightedDiGraph(n)
I = [ 1, 1, 2, 3,  3, 4, 4, 5,  5]
J = [ 2, 3, 5, 2,  4, 5, 6, 3,  6]
W = [ 9, 6, 8, 4, 11, 6, 4, 7, 15]
for (i,j,w) in zip(I,J,W)
    add_edge!(g,i,j,w)
end
edgelabel_dict = Dict()
for i in 1:n, j in 1:n
    if has_edge(g,i,j)
        # edgelabel_dict[(i, j)] = g.weights[j,i]
        edgelabel_dict[(i, j)] = get_weight(g,i,j)
    end
end
Random.seed!(42)
GraphRecipes.graphplot(g, names=1:n, edgelabel=edgelabel_dict, curves=false, nodeshape=:rect)

using JuMP
using HiGHS

s = 1
t = 2
V = 1:n
model = Model(HiGHS.Optimizer)
@variable(model, x[i in 1:n, j in 1:n; has_edge(g, (i,j))] .≥ 0)
@objective(model, Max, sum([x[s,i] for i in outneighbors(g,s)]))
@constraint(model, c_flow[v in V; v≠s&&v≠t], sum(x[v,i] for i in outneighbors(g,v)) == sum(x[j,v] for j in inneighbors(g,v)))
@constraint(model, c_max[v1 in V, v2 in V; has_edge(g,v1,v2)], x[v1,v2] ≤ get_weight(g,v1,v2))
set_silent(model)
optimize!(model)
@show termination_status(model)
@show objective_value(model)
for x in x
    println(x, value(x))
end

julia> maximum_flow(g, 1, 2, g.weights')
(13.0, sparse([2, 3, 1, 3, 1, 2, 4, 5, 3, 5, 3, 4, 6, 5], [1, 1, 2, 2, 3, 3, 3, 3, 4, 4, 5, 5, 5, 6], [-9.0, -4.0, 9.0, 4.0, 4.0, -4.0, -2.0, 2.0, 2.0, -2.0, -2.0, 2.0, 0.0, 0.0], 6, 6))

julia> maximum_flow(g, 1, 2)
(2, sparse([2, 3, 1, 3, 1, 2], [1, 1, 2, 2, 3, 3], [-1, -1, 1, 1, 1, -1], 6, 6))