Discrete Optimization
By Prof. Seungchul Lee
http://iai.postech.ac.kr/
Industrial AI Lab at POSTECH
http://iai.postech.ac.kr/
Industrial AI Lab at POSTECH
Table of Contents
1. Optimization¶
An important tool in 1) engineering problem solving and 2) decision science
- peolple optimize
- nature optimizes
3 key components
- objective
- decision variable or unknown
- constraints
Procedures
- The process of identifying objective, variables, and constraints for a given problem is known as "modeling"
- Once the model has been formulated, optimization algorithm can be used to find its solutions.
In mathematical expression
Remarks) equivalent
1.1. Modeling Example¶
Manufacturer produces two products P and Q with machine A and B.
time of producing each unit of P
- in machine A: 50 min
- in machine B: 30 min
time of producing each unit of Q
- in machine A: 24 min
- in machine B: 33 min
working plan for a week $\begin{cases} \text{40 hrs of work on machine A} \\ \text{35 hrs of work on machine B} \end{cases}$
The week starts with $\begin{cases} \text{stock of 30 units of P} \\ \text{stock of 90 units of Q} \\ \text{demand of 75 units of P} \\ \text{demand of 95 units of Q} \\ \end{cases}$
Question: how to plan the production to end the week with the maximum stock?
Define decision variables
- $x$ = unit of P to be produced
- $y$ = unit of Q to be produced
1.2. General or Standard Forms¶
1) Linear programming $$\begin{align} \min_{x} \quad & f^Tx \\ \text{subject to} \quad & Ax \leq b\\ & A_{eq}x \leq b_{eq}\\ & LB \leq x \leq UB \end{align} $$
2) Integer programming $$\begin{align} \min_{x} \quad & f^Tx \\ \text{subject to} \quad & Ax \leq b\\ & A_{eq}x \leq b_{eq}\\ & LB \leq x \leq UB\\ & x \in \mathbb{Z} \text{(integer)} \end{align} $$
3) Mixed integer linear programming $$\begin{align} \min_{x,z} \quad & f^T \begin{bmatrix} x \\ z \end{bmatrix} \\ \text{subject to} \quad & A \begin{bmatrix} x \\ z \end{bmatrix} \leq b \\ & A_{eq} \begin{bmatrix} x \\ z \end{bmatrix} = b_{eq}\\ & LB \leq \begin{bmatrix} x \\ z \end{bmatrix} \leq UB \\ & x \in \mathbb{R} \text{(real)} \\ & z \in \mathbb{Z} \end{align} $$
4) Binary integer programming $$\begin{align} \min_{x} \quad & f^Tx \\ \text{subject to} \quad & Ax \leq b \\ & A_{eq}x = b_{eq} \\ & LB \leq x \leq UB \\ & x \in \{0,1\} \end{align} $$
5) Quadratic programming $$\begin{align} \min_{x} \quad & \frac{1}{2}x^THx + f^Tx \\ \text{subject to} \quad & Ax \leq b \\ & A_{eq}x = b_{eq} \\ & LB \leq x \leq UB \end{align} $$
2. Linear Programming (LP)¶
LP Example 1
In [1]:
f = -[1 1]';
A = [2 1;
1 2];
b = [29,25]';
lb = [2 5];
linprog(f,A,b,[],[],lb,[])
Out[1]:
LP example 2
In [2]:
f = -[1 1 1 1 1]';
Aeq = [1 1 0 0 0;
0 0 1 1 0;
0 1 1 2 5];
beq = [100 70 250]';
lb = [0 0 0 0 0]';
linprog(f,[],[],Aeq,beq,lb,[])
Out[2]:
LP example 3
In [3]:
f = -[1 1]';
A = [2 1];
b = log(0.25) - log(pi) + 2*log(2);
lb = [log(0.06) 0]';
x = linprog(f,A,b,[],[],lb,[])
d = exp(x(1))
n = exp(x(2))
Out[3]:
3. Mixed Integer Linear Programming (MILP)¶
- Matlab command: intlinprog (available from MATLAB R2014b)
In [4]:
help intlinprog
Out[4]:
MILP Example 1
In [5]:
f = [-3 -1]';
intcon = 1:2;
A = [17 11];
b = 83;
LB = [0 0]';
UB = [3 7]';
x = intlinprog(f,intcon,A,b,[],[],LB,UB)
Out[5]:
MIP Example 2
In [6]:
f = [-3; -2; -1];
intcon = 3;
A = [1, 1, 1];
b = 7;
Aeq = [4,2,1];
beq = 12;
lb = zeros(3,1);
ub = [Inf;Inf;1];
x = intlinprog(f,intcon,A,b,Aeq,beq,lb,ub)
Out[6]:
3.1. Knapsack problems revisit¶
- There is a container(knapsack) of capacity C = 20. Further more, there is a set 6 of objects. Each object has a weight and a value.
- Determine the number of each item to include in a collection so that the total size is less than or equal to 20 and the total value is as large as possible.
| items | 1 | 2 | 3 | 4 | 5 | 6 |
|---|---|---|---|---|---|---|
| weights | 10 | 9 | 4 | 2 | 1 | 20 |
| values | 175 | 90 | 20 | 50 | 10 | 200 |
Question: Can we formulate this knapsack problem into a binary integer program?
$\Rightarrow$ decision variables $\quad x_{i} \in \{0,1\} \quad i = 1,\cdots,6$
$\Rightarrow$ intlinprog in Matlab to solve
In [7]:
n = 6; % # of items
weights = [10 9 4 2 1 20];
values = [175 90 20 50 10 200];
maxWeight = 20;
% in a form of MILP
f = -values; % max
A = weights;
b = maxWeight;
intcon = 1:n;
lb = zeros(n,1);
ub = ones(n,1);
x = intlinprog(f,intcon,A,b,[],[],lb,ub)
Out[7]:
3.2. Factory and warehouse¶
- A campany wants to $\begin{cases} \text{1. build new factory in Ulsan, Busan, or both} \\ \text{2. build new warehouse (only one has to be built)} \\ \text{3. warehouse must be built in a city of a new factory} \\ \end{cases} $
- Available capital is 10. We want to maximize total profit.
| Profit | Cost to build | ||
|---|---|---|---|
| 1 | Build a factory in Ulsan | 7 | 5 |
| 2 | Build a factory in Busan | 5 | 3 |
| 3 | Build a warehouse in Ulsan | 6 | 5 |
| 4 | Build a warehouse in Busan | 4 | 2 |
- Form a binary integer programming to find the optimal locations of new factory and warehouse, and solve it using intlinprog in MATLAB.
$\quad\;$ (Hint: define binary decision variables, then use them to define objective function and constraints.)
In [8]:
f = -[7 5 6 4]';
intcon = [1 2 3 4];
A = [5 3 5 2;
-1 -1 0 0;
-1 0 1 0;
0 -1 0 1];
b = [10;
-1;
0;
0];
Aeq = [0 0 1 1];
beq = [1];
lb = [0 0 0 0]';
ub = [1 1 1 1]';
[x fval] = intlinprog(f,intcon,A,b,Aeq,beq,lb,ub)
Out[8]:
3.3. Shortest path¶
- Binary integer programming can be applied to finding the shortest path in the graph.
- Define the binary decision variables and form the binary integer programming. Then you can solve it using intlinprog command in MATLAB
In [9]:
f = [1 4 2 4 8 11 3]';
Aeq = [1 1 0 0 0 0 0
1 0 -1 0 0 -1 0
0 1 1 -1 -1 0 0
0 0 0 1 0 0 -1
0 0 0 0 1 1 1];
beq = [1 0 0 0 1]';
intcon = 1:7;
lb = zeros(7,1);
ub = ones(7,1);
x = intlinprog(f,intcon,[],[],Aeq,beq,lb,ub)
Out[9]:
3.4. Robot Path Design¶
- A circular mobile robot is to be navigated. The robot should take a path from point A to point B, avoiding the polygonal obstacles (hatched regions).
- Design a collision free-path the robot can take with the minimum length.
3.5. Vertex cover problem¶
- Suppose you want to have an ATM in every street (block, district) of a city of Ulsan. However, each additional ATM machine costs more money. Therefore, we want to have the minimum number of ATM machines at the right location while they are accessible from every street.
In [10]:
n = 8;
f = [1 1 1 1 1 1 1 1]';
intcon = 1:n;
A = -[1 0 1 0 0 0 0 0;
1 0 0 1 0 0 0 0;
1 0 0 0 1 0 0 0;
1 0 0 0 0 1 0 0;
1 0 0 0 0 0 1 0;
0 1 0 1 0 0 0 0;
0 1 0 0 1 0 0 0;
0 0 1 0 0 0 0 1;
0 0 0 0 1 0 1 0];
b = -[1 1 1 1 1 1 1 1 1]';
lb = zeros(n,1);
ub = ones(n,1);
x = intlinprog(f,intcon,A,b,[],[],lb,ub)
Out[10]:
3.6. Coloring problems¶
- Find the minimum number of colors that can cover the below Ulsan map without conflicting the adjacent rule.
- Coloring problem: color the vertext of a graph in such a way that adjacent vertices are differently colored, while minimizing the number of colors.
- Define ILP decision variables
In [11]:
%% coloring problem using binary integer programming
% Actually number of colours is at most 5
% Therefore we restrict our number of colours to be 5
% vertex V1:울주군, V2:중구, V3:북구, V4:남구, V5:동구
% edges (or adjacent)
E = [1,2;
1,3;
1,4;
2,3;
2,4;
3,4;
3,5;
4,5];
%% equality constraints
A1 = zeros(5, 30);
for i = 1:5
A1(i,1+(i-1)*5:i*5) = 1;
end
b1 = ones(5,1);
%% inequality
A2 = zeros(length(E)*5, 30);
for k = 1:5
for i = 1:length(E)
A2(i+length(E)*(k-1),5*(E(i,1)-1)+k) = 1;
A2(i+length(E)*(k-1),5*(E(i,2)-1)+k) = 1;
end
end
b2 = ones(length(E)*5, 1);
A3 = zeros(25, 30);
for k = 1:5
for i = 1:5
A3(i+5*(k-1),5*(i-1)+k) = 1; % x_ik
A3(i+5*(k-1),25+k) = -1; % -y_k
end
end
b3 = zeros(25, 1);
%% plug-in to intlinprog standard form
f = [zeros(1,25), ones(1,5)]';
intcon = 1:30;
Aeq = A1;
beq = b1;
A = [A2; A3];
b = [b2; b3];
lb = zeros(30, 1);
ub = ones(30, 1);
solution = intlinprog(f,intcon,A,b,Aeq,beq,lb,ub);
x = solution(1:25);
X = [];
for i = 1:5
X(i,:) = x(1+(i-1)*5:5+(i-1)*5);
end
y = solution(26:30)
Out[11]:
3.7. Flat pattern desing of 3D shapes¶
- A box with the dimensions specified is to be made by bending and welding a flat pattern cut out of sheet metal. Design the at pattern with minimum welding length that can be made to the box.
Constraints
- No disconnected node.
- Only 6-1 edges are required.
In [12]:
n = 12;
f = -[1 1 1 1 2 2 2 2 3 3 3 3]';
intcon = 1:n;
A = -[1 1 0 0 1 0 0 1 0 0 0 0;
0 1 0 1 0 0 0 0 0 1 0 1;
0 0 1 1 0 1 1 0 0 0 0 0;
1 0 1 0 0 0 0 0 1 0 1 0;
0 0 0 0 1 1 0 0 1 1 0 0;
0 0 0 0 0 0 1 1 0 0 1 1];
b = -ones(6,1);
Aeq = ones(1,n);
beq = 5;
lb = zeros(n,1);
ub = ones(n,1);
x = intlinprog(f,intcon,A,b,Aeq,beq,lb,ub)
Out[12]:
3.8. Traveling salesman problem¶
- Starting from his home, a salesman wishes to visit each of (n-1) other cities and return home at minimal cost. He must visit each city exactly once and it cost $c_{ij}$ to travel from city $i$ to city $j$.
- What route should he select?
Starting from node 1, visit 2,3,4 once and return node 1 at minimal cost.
In [13]:
n = 4;
f = [5 10 13 5 8 9 10 8 7 13 9 7]';
intcon = 1:n*(n-1);
lb = zeros(n*(n-1),1);
ub = ones(n*(n-1),1);
Aeq = [1 1 1 0 0 0 0 0 0 0 0 0;
0 0 0 1 1 1 0 0 0 0 0 0;
0 0 0 0 0 0 1 1 1 0 0 0;
0 0 0 0 0 0 0 0 0 1 1 1;
0 0 0 1 0 0 1 0 0 1 0 0;
1 0 0 0 0 0 0 1 0 0 1 0;
0 1 0 0 1 0 0 0 0 0 0 1;
0 0 1 0 0 1 0 0 1 0 0 0];
beq = [1 1 1 1 1 1 1 1]';
x = intlinprog(f,intcon,[],[],Aeq,beq,lb,ub)
Out[13]:
what is wrong?
In [14]:
A = -[0 1 1 0 1 1 0 0 0 0 0 0;
0 0 0 0 0 0 1 1 0 1 1 0];
b = -[1 1]';
x = intlinprog(f,intcon,A,b,Aeq,beq,lb,ub)
Out[14]:
4. Quadratic Programming (QP)¶
4.1. Best location at a concert¶
Question: Find the best location to listen to singer's voice
In [15]:
H = [2 0;
0 2];
f = -[6 6]';
A = [1 1];
b = 5;
lb = [0 0]';
x = quadprog(H,f,A,b,[],[],lb,[])
Out[15]:
4.2. Linear regression¶
Data fitting or approximation
Given $(x_1,y_1), (x_2,y_2),\cdots,(x_m,y_m)$, find $y_i \approx \theta_1 x_i + \theta_0$
such that $ \min_{\theta_1,\theta_2} \; \sum \limits_{i=1}^m(\hat{y}_i-y_i)^{2}$
Optimization problem: find $\theta_1$ and $\theta_0$ such that minimize total sum of $\text{error}^{2} \; ( \text{error}_i = \big \vert \hat{y}_i-y_i \big \vert )$ + no constraints on $\theta_1$ and $\theta_0$
Remark) norm of vector $u = \begin{bmatrix} u_1 \\ u_2 \\ \vdots \\ u_m \end{bmatrix} = \Vert u \Vert = \sqrt{u_1^{2}+u_2^{2}+\cdots+u_m^{2}}$
use quadprog($ 2 X^T X$, $- 2 X^T Y$)
In [16]:
% data
x = [0.1 0.4 0.7 1.2 1.3 1.7 2.2 2.8 3.0 4.0 4.3 4.4 4.9]';
y = [0.5 0.9 1.1 1.5 1.5 2.0 2.2 2.8 2.7 3.0 3.5 3.7 3.9]';
m = length(x);
X = [x ones(m,1)];
H = 2*X'*X;
f = -2*X'*y;
theta = quadprog(H,f);
yhat = X*theta;
figure(1), clf
plot(x,y,'o',x,yhat,'linewidth',2)
Out[16]:
In [17]:
%%javascript
$.getScript('https://kmahelona.github.io/ipython_notebook_goodies/ipython_notebook_toc.js')