{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"---\n",
"title: Facility Location\n",
"---"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Originally Contributed by**: Mathieu Tanneau and Alexis Montoison"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"using Random\n",
"using LinearAlgebra\n",
"\n",
"using JuMP\n",
"import GLPK\n",
"using Plots"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Uncapacitated facility location"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Problem description\n",
"\n",
"We are given\n",
"* $M=\\{1, \\dots, m\\}$ clients\n",
"* $N=\\{ 1, \\dots, n\\}$ sites where a facility can be built\n",
"\n",
"**Decision variables**\n",
"Decision variables are split into two categories:\n",
"* Binary variable $y_{j}$ indicates whether facility $j$ is built or not\n",
"* Binary variable $x_{i, j}$ indicates whether client $i$ is assigned to facility $j$\n",
"\n",
"**Objective**\n",
"The objective is to minimize the total cost of serving all clients.\n",
"This costs breaks down into two components:\n",
"* Fixed cost of building a facility.\n",
"In this example, this cost is $f_{j} = 1, \\ \\forall j$.\n",
"* Cost of serving clients from the assigned facility.\n",
"In this example, the cost $c_{i, j}$\n",
"of serving client $i$ from facility $j$\n",
"is the Euclidean distance between the two.\n",
"\n",
"**Constraints**\n",
"* Each customer must be served by exactly one facility\n",
"* A facility cannot serve any client unless it is open"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### MILP formulation\n",
"\n",
"The problem can be formulated as the following MILP:\n",
"\n",
"\\begin{align}\n",
"\\min_{x, y} \\ \\ \\ &\n",
"\\sum_{i, j} c_{i, j} x_{i, j} + \n",
"\\sum_{j} f_{j} y_{j}\\\\\n",
"s.t. &\n",
"\\sum_{j} x_{i, j} = 1, && \\forall i \\in M\\\\\n",
"& x_{i, j} \\leq y_{j}, && \\forall i \\in M, j \\in N\\\\\n",
"& x_{i, j}, y_{j} \\in \\{0, 1\\}, && \\forall i \\in M, j \\in N\n",
"\\end{align}\n",
"\n",
"where the first set of constraints ensures\n",
"that each client is served exactly once,\n",
"and the second set of constraints ensures\n",
"that no client is served from an unopened facility."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Problem data"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"Random.seed!(314)\n",
"\n",
"m = 12 # number of clients\n",
"n = 5 # number of facility locations\n",
"\n",
"# Clients' locations\n",
"Xc = rand(m)\n",
"Yc = rand(m)\n",
"\n",
"# Facilities' potential locations\n",
"Xf = rand(n)\n",
"Yf = rand(n)\n",
"\n",
"# Fixed costs\n",
"f = ones(n);\n",
"\n",
"# Distance\n",
"c = zeros(m, n)\n",
"for i in 1:m\n",
" for j in 1:n\n",
" c[i, j] = norm([Xc[i] - Xf[j], Yc[i] - Yf[j]], 2)\n",
" end\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"\n"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Display the data\n",
"scatter(Xc, Yc, label = \"Clients\", markershape=:circle, markercolor=:blue)\n",
"scatter!(Xf, Yf, label=\"Facility\", \n",
" markershape=:square, markercolor=:white, markersize=6,\n",
" markerstrokecolor=:red, markerstrokewidth=2\n",
")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### JuMP implementation"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"# Create a JuMP model\n",
"ufl = Model(GLPK.Optimizer)\n",
"\n",
"# Variables\n",
"@variable(ufl, y[1:n], Bin);\n",
"@variable(ufl, x[1:m, 1:n], Bin);\n",
"\n",
"# Each client is served exactly once\n",
"@constraint(ufl, client_service[i in 1:m],\n",
" sum(x[i, j] for j in 1:n) == 1\n",
");\n",
"\n",
"# A facility must be open to serve a client\n",
"@constraint(ufl, open_facility[i in 1:m, j in 1:n],\n",
" x[i, j] <= y[j]\n",
")\n",
"\n",
"# Objective\n",
"@objective(ufl, Min, f'y + sum(c .* x));"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Optimal value: 5.226417970467934\n"
]
}
],
"source": [
"# Solve the uncapacitated facility location problem with GLPK\n",
"optimize!(ufl)\n",
"println(\"Optimal value: \", objective_value(ufl))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Visualizing the solution"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"# The threshold 1e-5 ensure that edges between clients and facilities are drawn when x[i, j] ≈ 1.\n",
"x_ = value.(x) .> 1 - 1e-5\n",
"y_ = value.(y) .> 1 - 1e-5\n",
"\n",
"# Display clients\n",
"p = scatter(Xc, Yc, markershape=:circle, markercolor=:blue, label=nothing)\n",
"\n",
"# Show open facility\n",
"mc = [(y_[j] ? :red : :white) for j in 1:n]\n",
"scatter!(Xf, Yf, \n",
" markershape=:square, markercolor=mc, markersize=6,\n",
" markerstrokecolor=:red, markerstrokewidth=2,\n",
" label=nothing\n",
")\n",
"\n",
"# Show client-facility assignment\n",
"for i in 1:m\n",
" for j in 1:n\n",
" if x_[i, j] == 1\n",
" plot!([Xc[i], Xf[j]], [Yc[i], Yf[j]], color=:black, label=nothing)\n",
" end\n",
" end\n",
"end\n",
"\n",
"display(p)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Capacitated Facility location"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Problem formulation\n",
"\n",
"The capacitated variant introduces a capacity constraint on each facility, i.e., clients have a certain level of demand to be served, while each facility only has finite capacity which cannot be exceeded.\n",
"\n",
"Specifically, let\n",
"* $a_{i} \\geq 0$ denote the demand of client $i$\n",
"* $q_{j} \\geq 0$ denote the capacity of facility $j$\n",
"\n",
"The capacity constraints then write\n",
"\\begin{align}\n",
"\\sum_{i} a_{i} x_{i, j} &\\leq q_{j} y_{j} && \\forall j \\in N\n",
"\\end{align}\n",
"\n",
"Note that, if $y_{j}$ is set to $0$, the capacity constraint above automatically forces $x_{i, j}$ to $0$."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Thus, the capacitated facility location can be formulated as follows\n",
"\n",
"\\begin{align}\n",
"\\min_{x, y} \\ \\ \\ &\n",
"\\sum_{i, j} c_{i, j} x_{i, j} + \n",
"\\sum_{j} f_{j} y_{j}\\\\\n",
"s.t. &\n",
"\\sum_{j} x_{i, j} = 1, && \\forall i \\in M\\\\\n",
"& \\sum_{i} a_{i} x_{i, j} \\leq q_{j} y_{j}, && \\forall j \\in N\\\\\n",
"& x_{i, j}, y_{j} \\in \\{0, 1\\}, && \\forall i \\in M, j \\in N\n",
"\\end{align}\n",
"\n",
"For simplicity, we will assume that there is enough capacity to serve the demand,\n",
" i.e., there exists at least one feasible solution."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"# Demands\n",
"a = rand(1:3, m);\n",
"\n",
"# Capacities\n",
"q = rand(5:10, n);"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"\n"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Display the data\n",
"scatter(Xc, Yc, label=nothing,\n",
" markershape=:circle, markercolor=:blue, markersize= 2 .*(2 .+ a)\n",
")\n",
"\n",
"scatter!(Xf, Yf, label=nothing, \n",
" markershape=:rect, markercolor=:white, markersize= q,\n",
" markerstrokecolor=:red, markerstrokewidth=2\n",
")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### JuMP implementation"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [],
"source": [
"# Create a JuMP model\n",
"cfl = Model(GLPK.Optimizer)\n",
"\n",
"# Variables\n",
"@variable(cfl, y[1:n], Bin);\n",
"@variable(cfl, x[1:m, 1:n], Bin);\n",
"\n",
"# Each client is served exactly once\n",
"@constraint(cfl, client_service[i in 1:m], sum(x[i, :]) == 1)\n",
"\n",
"# Capacity constraint\n",
"@constraint(cfl, capacity, x'a .<= (q .* y))\n",
"\n",
"# Objective\n",
"@objective(cfl, Min, f'y + sum(c .* x));"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Optimal value: 5.836938427763575\n"
]
}
],
"source": [
"# Solve the problem\n",
"optimize!(cfl)\n",
"println(\"Optimal value: \", objective_value(cfl))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Visualizing the solution"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"\n"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"# The threshold 1e-5 ensure that edges between clients and facilities are drawn when x[i, j] ≈ 1.\n",
"x_ = value.(x) .> 1 - 1e-5\n",
"y_ = value.(y) .> 1 - 1e-5\n",
"\n",
"# Display the solution\n",
"p = scatter(Xc, Yc, label=nothing,\n",
" markershape=:circle, markercolor=:blue, markersize= 2 .*(2 .+ a)\n",
")\n",
"\n",
"mc = [(y_[j] ? :red : :white) for j in 1:n]\n",
"scatter!(Xf, Yf, label=nothing, \n",
" markershape=:rect, markercolor=mc, markersize=q,\n",
" markerstrokecolor=:red, markerstrokewidth=2\n",
")\n",
"\n",
"# Show client-facility assignment\n",
"for i in 1:m\n",
" for j in 1:n\n",
" if x_[i, j] == 1\n",
" plot!([Xc[i], Xf[j]], [Yc[i], Yf[j]], color=:black, label=nothing)\n",
" break\n",
" end\n",
" end\n",
"end\n",
"\n",
"display(p)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Further\n",
"* [Benders decomposition](https://github.com/JuliaOpt/JuMPTutorials.jl/blob/master/script/optimization_concepts/benders_decomposition.jl)\n",
"is a method of choice for solving facility location problems.\n",
"* Benchmark instances can be found\n",
"[here](https://resources.mpi-inf.mpg.de/departments/d1/projects/benchmarks/UflLib/)."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Julia 1.5.1",
"language": "julia",
"name": "julia-1.5"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.5.1"
}
},
"nbformat": 4,
"nbformat_minor": 2
}