# MS-E2121 Exercise session 9
### Problem 9.1: Uncapacitated Facility Location (UFL)

**a)**  
Let $N = \{1,\dots,n\}$ be a set of potential facilities and $M = \{1,\dots,m\}$ a set of clients. Let $y_j = 1$ if facility $j$ is opened, and $y_j = 0$ otherwise. Moreover, let $x_{ij}$ be the fraction of client $i$'s demand satisfied from facility $j$. The UFL can be formulated as the mixed-integer problem (MIP): 

$$\begin{align}
 \text{(UFL-W)} : \quad &\min_{x,y} \sum_{j\in N} f_jy_j + \sum_{i\in M}\sum_{j\in N} c_{ij}x_{ij} \\
       &\text{s.t.} \\
       &\quad \sum_{j\in N}x_{ij} = 1, &\forall i \in M,\\
       &\quad \sum_{i\in M}x_{ij} \leq my_j, &\forall  j \in N,\\
       &\quad x_{ij} \geq 0, &\forall i \in M, \forall j \in N,\\
       &\quad y_j \in \{0,1\}, &\forall j\in N,
\end{align}$$

where $f_j$ is the cost of opening facility $j$, and $c_{ij}$ is the cost of satisfying client $i$'s demand from facility $j$. Consider an instance of the UFL with opening costs $f=(4,3,4,4,7)$ and client costs

$$\begin{align*}
 (c_{ij}) = \left(
	\begin{array}{ccccc}
		12 & 13 & 6 & 0  & 1 \\
		8  & 4  & 9 & 1  & 2 \\
		2  & 6  & 6 & 0  & 1 \\
		3  & 5  & 2 & 1  & 8 \\
		8  & 0  & 5 & 10 & 8 \\
		2  & 0  & 3 & 4  & 1
	\end{array}
 \right)
\end{align*}$$

Implement (the model) and solve the problem with Julia using JuMP.

**b)**  
An alternative formulation of the UFL is of the form

$$\begin{align}
 \text{(UFL-S)} : \quad &\min_{x,y} \sum_{j\in N}f_jy_j + \sum_{i\in M}\sum_{j\in N}c_{ij}x_{ij}\\
       &\text{s.t.}  \\
       &\quad \sum_{j\in N}x_{ij} = 1, &\forall i \in M,\\
       &\quad x_{ij} \leq y_j, &\forall  i\in M, \forall j \in N,\\
       &\quad x_{ij} \geq 0, &\forall i \in M, \forall j \in N,\\
       &\quad y_j \in \{0,1\}, &\forall j\in N.
\end{align}$$


Linear programming (LP) relaxations of these problems can be obtained by relaxing the binary constraints $y_j\in \{0,1\}$ to $0 \leq y_j \leq 1$ for all $j \in N$. For the same instance as in part (a), solve the LP relaxations of UFL-W and UFL-S with Julia using JuMP, and compare the optimal costs of the LP relaxations against the optimal integer cost obtained in part (a).


```julia
using JuMP, Cbc
```

Write down the problem data


```julia
f = [4 3 4 4 7] # Facility opening costs
c = [12 13 6 0 1; 8 4 9 1 2; 2 6 6 0 1; 3 5 2 1 8; 8 0 5 10 8; 2 0 3 4 1] # Cost of satisfying demand
(m, n) = size(c)
M = 1:m # Set of facilities
N = 1:n;# Set of clients
```

Implement the problem in JuMP


```julia
ufl_w = Model(Cbc.Optimizer)

@variable(ufl_w, x[M,N] >= 0) # Fraction of demand (client i) satisfied by facility j
@variable(ufl_w, y[N], Bin)   # Facility location

# Minimize total cost
@objective(ufl_w, Min, sum(f[j]*y[j] for j in N) + sum(c[i,j]*x[i,j] for i in M, j in N)) 

# For each client, the demand must be fulfilled
@constraint(ufl_w, demand[i in M], sum(x[i,j] for j in N) == 1)
# A big-M style constraint stating that facility j can't send out anything if y[j]==0
@constraint(ufl_w, supply[j in N], sum(x[i,j] for i in M) <= m*y[j])

optimize!(ufl_w)
```


```julia
println("UFL-W MILP:")
println("Optimal value $(objective_value(ufl_w))")
println("with y = $(value.(y).data)")
```


```julia
ufl_w_rel = Model(Cbc.Optimizer)

@variable(ufl_w_rel, x[M,N] >= 0) # Fraction of demand (client i) satisfied by facility j
@variable(ufl_w_rel, 0<=y[N]<=1)  # Facility location

# Minimize total cost
@objective(ufl_w_rel, Min, sum(f[j]*y[j] for j in N) + sum(c[i,j]*x[i,j] for i in M, j in N)) 

# For each client, the demand must be fulfilled
@constraint(ufl_w_rel, demand[i in M], sum(x[i,j] for j in N) == 1)
# A big-M style constraint stating that facility j can't send out anything if y[j]==0
@constraint(ufl_w_rel, supply[j in N], sum(x[i,j] for i in M) <= m*y[j])

optimize!(ufl_w_rel)
```


```julia
println("UFL-W LP:")
println("Optimal value $(objective_value(ufl_w_rel))")
println("with y = $(value.(y).data)")
```


```julia
ufl_s_rel = Model(Cbc.Optimizer)

@variable(ufl_s_rel, x[M,N] >= 0)
@variable(ufl_s_rel, 0<=y[N]<=1)

@objective(ufl_s_rel, Min, sum(f[j]*y[j] for j in N) + sum(c[i,j]*x[i,j] for i in M, j in N))

@constraint(ufl_s_rel, demand[i in M], sum(x[i,j] for j in N) == 1)
# The difference between the models is that UFL-S has m constraints telling that nothing can be sent to client i from facility j if y[j]==0
# In UFL-W, there is a single constraint telling that nothing can be sent from facility j if y[j]==0
@constraint(ufl_s_rel, supply[i in M, j in N], x[i,j] <= y[j])

optimize!(ufl_s_rel)
```


```julia
println("UFL-S LP:")
println("Optimal value $(objective_value(ufl_s_rel))")
println("with y = $(value.(y).data)")
```

#### Branching
We see that the UFL-S relaxation produces an integer solution, meaning that we have an integer optimal solution and no branching needs to be done. However, if we used UFL-W instead, we would need to do B&B or something else to obtain the integer optimum. In the UFL-W LP relaxation solution (0, 1/3, 0, 2/3, 0), we have two fractional variables $y_2$ and $y_4$, and we can branch on one of them. Let's choose $y_2$ and see what happens if we set it to 0 or 1.


```julia
ufl_w_rel_y2_0 = Model(Cbc.Optimizer)

@variable(ufl_w_rel_y2_0, x[M,N] >= 0)
@variable(ufl_w_rel_y2_0, 0<=y[N]<=1)

@objective(ufl_w_rel_y2_0, Min, sum(f[j]*y[j] for j in N) + sum(c[i,j]*x[i,j] for i in M, j in N))

@constraint(ufl_w_rel_y2_0, demand[i in M], sum(x[i,j] for j in N) == 1)
@constraint(ufl_w_rel_y2_0, supply[j in N], sum(x[i,j] for i in M) <= m*y[j])
@constraint(ufl_w_rel_y2_0, y[2] == 0)

optimize!(ufl_w_rel_y2_0)
```


```julia
ufl_w_rel_y2_1 = Model(Cbc.Optimizer)

@variable(ufl_w_rel_y2_1, x[M,N] >= 0)
@variable(ufl_w_rel_y2_1, 0<=y[N]<=1)

@objective(ufl_w_rel_y2_1, Min, sum(f[j]*y[j] for j in N) + sum(c[i,j]*x[i,j] for i in M, j in N))

@constraint(ufl_w_rel_y2_1, demand[i in M], sum(x[i,j] for j in N) == 1)
@constraint(ufl_w_rel_y2_1, supply[j in N], sum(x[i,j] for i in M) <= m*y[j])
@constraint(ufl_w_rel_y2_1, y[2] == 1)

optimize!(ufl_w_rel_y2_1)
```


```julia
println("UFL-W LP with y2=0:")
println("Optimal value $(objective_value(ufl_w_rel_y2_0))")
println("with y = $(value.(ufl_w_rel_y2_0[:y]).data)")

println()

println("UFL-W LP with y2=1:")
println("Optimal value $(objective_value(ufl_w_rel_y2_1))")
println("with y = $(value.(ufl_w_rel_y2_1[:y]).data)")
```

Both branches have fractional solutions, and more branching is thus needed. You can practice that in the next exercise.
#' 
### Problem 9.3: Solving Branch & Bound (B&B) graphically

*You can do this with pen and paper if you want to do it graphically, or solve the problems using JuMP instead if you don't feel like drawing.*

Consider the following integer programming problem $IP$:

$$\begin{matrix}
\text{max} &x_{1} &+&2x_{2} & \\
\text{s.t.}&-3x_{1} &+&4x_{2} &\le 4 \\
&3x_{1} &+&2x_{2} &\le 11 \\
&2x_{1} &-&x_{2} &\le 5 \\
&x_{1}, &x_{2} & \text{integer} &\\
\end{matrix}$$

Plot (or draw) the feasible region of the linear programming (LP) relaxation of the problem $IP$, then solve the problems using the figure. Recall that the LP relaxation of $IP$ is obtained by replacing the integrality constraints $x_1,x_2\in \mathbb{Z}_+$ by linear nonnegativity $x_1,x_2\geq 0$ and including (possible) upper bounds corresponding to the upper bounds of the integer variables ($x_1,x_2\leq 1$ for binary variables). 

(a) What is the optimal cost $z_{LP}$ of the LP relaxation of the problem $IP$? What is the optimal cost $z$ of the problem $IP$?

(b) Draw the border of the convex hull of the feasible solutions of the problem $IP$. Recall that the convex hull represents the *ideal* formulation for the problem $IP$.

(c) Solve the problem $IP$ by LP-relaxation based Branch \& Bound (B\&B). You can solve the LP relaxations at each node of the B\&B tree graphically. Start the B\&B procedure without any primal bound.

Check your solutions using JuMP. Make sure to point out the optimal solutions in the figure, as well as giving their numerical values.


```julia
# TODO: add your code here
```


```julia
# TODO: add your code here
```


```julia
# TODO: add your code here
```


```julia
# TODO: add your code here
```


```julia
# TODO: add your code here
```


```julia
# TODO: add your code here
```


```julia
# TODO: add your code here
```


```julia
# TODO: add your code here
```


```julia
# TODO: add your code here
```


```julia

```

96c41b93117a5f4ccb0c73fbb68ab3c96f547684
ipynb
Jupyter Notebook
exercises/ex-09/exercise9_skeleton.ipynb
nnhjy/lp-julia
06c63cf06b6911bc4ab99a8c18b33bca2b61083c
MIT
exercises/ex-09/exercise9_skeleton.ipynb
nnhjy/lp-julia
06c63cf06b6911bc4ab99a8c18b33bca2b61083c
MIT
exercises/ex-09/exercise9_skeleton.ipynb
nnhjy/lp-julia
06c63cf06b6911bc4ab99a8c18b33bca2b61083c
MIT
Qwen/Qwen-72B
1. YES
2. YES