# The cannery problem

Original author: Louis Luangkesorn, January 30, 2015.

This tutorial solves the cannery problem from Dantzig, Linear Programming and Extensions, Princeton University Press, Princeton, NJ, 1963. This class of problem is known as a transshipment problem.

The purpose of this tutorial is to demonstrate how to use JSON data in the formulation of a JuMP model.

## Required packages

using JuMP
import HiGHS
import JSON

## Formulation

The cannery problem assumes we are optimizing the shipment of cans from production plants $p \in P$ to markets $m \in M$.

Each production plant $p$ has a capacity, $c_p$, and each market $m$ has a demand $d_m$. The distance from plant to market is $d_{p,m}$.

With a little effort, we can formulate our problem as the following linear program:

\begin{aligned} \min & \sum\limits_{p \in P}\sum\limits_{m \in M} d_{p,m} x_{p,m} \\ \text{s.t.} & \sum\limits_{m \in M} x_{p,m} \le c_p, && \forall p \in P \\ & \sum\limits_{p \in P} x_{p,m} \ge d_m, && \forall m \in M \\ & x_{p,m} \ge 0, && \forall p \in P, m \in M \end{aligned}

## Data

A key feature of the tutorial is to demonstrate how to load data from JSON.

For simplicity, we've hard-coded it below. But if the data was available as a .json file, we could use data = JSON.parsefile(filename) to read in the data.

data = JSON.parse("""
{
"plants": {
"Seattle": {"capacity": 350},
"San-Diego": {"capacity": 600}
},
"markets": {
"New-York": {"demand": 300},
"Chicago": {"demand": 300},
"Topeka": {"demand": 300}
},
"distances": {
"Seattle => New-York": 2.5,
"Seattle => Chicago": 1.7,
"Seattle => Topeka": 1.8,
"San-Diego => New-York": 2.5,
"San-Diego => Chicago": 1.8,
"San-Diego => Topeka": 1.4
}
}
""")
Dict{String, Any} with 3 entries:
"plants"    => Dict{String, Any}("Seattle"=>Dict{String, Any}("capacity"=>350…
"distances" => Dict{String, Any}("San-Diego => New-York"=>2.5, "Seattle => To…
"markets"   => Dict{String, Any}("Chicago"=>Dict{String, Any}("demand"=>300),…

Create the set of plants:

P = keys(data["plants"])
KeySet for a Dict{String, Any} with 2 entries. Keys:
"Seattle"
"San-Diego"

Create the set of markets:

M = keys(data["markets"])
KeySet for a Dict{String, Any} with 3 entries. Keys:
"Chicago"
"Topeka"
"New-York"

We also need a function to compute the distance from plant to market:

distance(p::String, m::String) = data["distances"]["$(p) =>$(m)"]
distance (generic function with 1 method)

## JuMP formulation

Now we're ready to convert our mathematical formulation into a JuMP model.

First, create a new JuMP model. Since we have a linear program, we'll use HiGHS as our optimizer:

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

Our decision variables are indexed over the set of plants and markets:

@variable(model, x[P, M] >= 0)
2-dimensional DenseAxisArray{VariableRef,2,...} with index sets:
Dimension 1, ["Seattle", "San-Diego"]
Dimension 2, ["Chicago", "Topeka", "New-York"]
And data, a 2×3 Matrix{VariableRef}:
x[Seattle,Chicago]    x[Seattle,Topeka]    x[Seattle,New-York]
x[San-Diego,Chicago]  x[San-Diego,Topeka]  x[San-Diego,New-York]

We need a constraint that each plant can ship no more than its capacity:

@constraint(model, [p in P], sum(x[p, :]) <= data["plants"][p]["capacity"])
1-dimensional DenseAxisArray{ConstraintRef{Model, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.LessThan{Float64}}, ScalarShape},1,...} with index sets:
Dimension 1, ["Seattle", "San-Diego"]
And data, a 2-element Vector{ConstraintRef{Model, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.LessThan{Float64}}, ScalarShape}}:
x[Seattle,Chicago] + x[Seattle,Topeka] + x[Seattle,New-York] ≤ 350.0
x[San-Diego,Chicago] + x[San-Diego,Topeka] + x[San-Diego,New-York] ≤ 600.0

and each market must receive at least its demand:

@constraint(model, [m in M], sum(x[:, m]) >= data["markets"][m]["demand"])
1-dimensional DenseAxisArray{ConstraintRef{Model, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.GreaterThan{Float64}}, ScalarShape},1,...} with index sets:
Dimension 1, ["Chicago", "Topeka", "New-York"]
And data, a 3-element Vector{ConstraintRef{Model, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.GreaterThan{Float64}}, ScalarShape}}:
x[Seattle,Chicago] + x[San-Diego,Chicago] ≥ 300.0
x[Seattle,Topeka] + x[San-Diego,Topeka] ≥ 300.0
x[Seattle,New-York] + x[San-Diego,New-York] ≥ 300.0

Finally, our objective is to minimize the transportation distance:

@objective(model, Min, sum(distance(p, m) * x[p, m] for p in P, m in M))
$$$1.7 x_{Seattle,Chicago} + 1.8 x_{Seattle,Topeka} + 2.5 x_{Seattle,New-York} + 1.8 x_{San-Diego,Chicago} + 1.4 x_{San-Diego,Topeka} + 2.5 x_{San-Diego,New-York}$$$

## Solution

optimize!(model)

solution_summary(model)
* Solver : HiGHS

* Status
Termination status : OPTIMAL
Primal status      : FEASIBLE_POINT
Dual status        : FEASIBLE_POINT
Message from the solver:
"kHighsModelStatusOptimal"

* Candidate solution
Objective value      : 1.68000e+03
Objective bound      : 0.00000e+00
Relative gap         : Inf
Dual objective value : 1.68000e+03

* Work counters
Solve time (sec)   : 3.50475e-04
Simplex iterations : 3
Barrier iterations : 0
Node count         : -1


What's the optimal shipment?

println("RESULTS:")
for p in P, m in M
println(p, " => ", m, ": ", value(x[p, m]))
end
RESULTS:
Seattle => Chicago: 300.0
Seattle => Topeka: 0.0
Seattle => New-York: 0.0
San-Diego => Chicago: 0.0
San-Diego => Topeka: 300.0
San-Diego => New-York: 300.0

Tip