simplesso rivisto - Università di Bologna

19/05/2004
Programmazione Matematica:
VI – Estensioni dell’algoritmo del
Simplesso
Daniele Vigo
D.E.I.S. – Università di Bologna
[email protected]
rev. 1.0 – Aprile 2004
Algoritmo del Simplesso
• L’algoritmo del Simplesso ad ogni pivot richiede
l’aggiornamento dell’intero tableau Y (m+1)×(n+1)
• Principalmente questo serve per trovare la colonna
su cui fare il pivoting (calcolo costi ridotti)
–c*0
0
c'FT
B–1d
I
B–1F
soluzione (serve)
D. Vigo
costi ridotti (basta
trovarne uno < 0)
di questa parte serve
solo la colonna su cui
si fa il pivot
Pmat.Rev.2
1
19/05/2004
1. Algoritmo del Simplesso Rivisto
• Prevede l’aggiornamento esplicito di una matrice
(m+1)×(m+1), che chiamiamo K
• K(0) matrice iniziale, K(i) matrice all’iterazione i
• I costi ridotti e le colonne aggiornate su cui si fa il
pivot vengono generate solo se e quando
necessario
• tempo per iterazione O(m2) << O(nm)
D. Vigo
Pmat.Rev.3
Informazioni del Tableau
• Il tableau corrente è quello iniziale moltiplicato
per l’inversa della base corrente
0
cBT
cFT
–c*0
0
c'FT
d
B
F
B–1d
I
B–1F
D. Vigo
Pmat.Rev.4
2
19/05/2004
Informazioni del Tableau (2)
• Se il tableau iniziale contiene una matrice identità
con costi ridotti 0, nelle iterazioni successive in
quella parte del tableau c’è B−1
0
0
cFT
–c*0
−cBT B−1
c'FT
d
I
F
B–1d
B−1
B–1F
K(0)
K(i)
• Conoscere K(i) è sufficiente per il pivoting
D. Vigo
Pmat.Rev.5
Pivoting rivisto
• Dato il tableau originale Y = (A, c, d ) ed K(i):
1. Generazione dei costi ridotti (pricing)
•
si calcolano una colonna per volta finchè non se ne
trova uno negativo (s) o si termina senza trovarne
c’j = cj − cB B−1 Aj = cj − cB K(i) Aj
2. Generazione della colonna (column generation)
A’j = B−1Aj = K(i) Aj
3. Determinazione dell’elemento pivot (riga r)
4. Aggiornamento di K(i+1) per la prossima
iterazione
D. Vigo
Pmat.Rev.6
3
19/05/2004
Vantaggi
• Nel pricing ci si può fermare alla prima colonna a
costo negativo (Es. regola di Bland)
• Per la column generation oltre a K(i) basta la
matrice A originale
• se A è sparsa la column generation se ne avvantaggia ed
è molto veloce
• La matrice K(i) può essere memorizzata in forma
prodotto (molto più compatta)
K(i) = Pi … P1 K(0) = Pi … P1
dove ogni Pi è un vettore
D. Vigo
Pmat.Rev.7
2. Problema in forma di massimo
•
•
L’algoritmo del simplesso fa riferimento alla
forma standard (f.ob. di minimo)
Se la funzione obiettivo è di massimo:
1. la si trasforma in forma di minimo
2. si usa la condizione di ottimalità modificata per la
funzione obiettivo di massimo:
“la base corrente è ottima se tutte le variabili non in
base hanno costo ridotto non positivo (≤ 0)”
D. Vigo
Pmat.Rev.8
4
19/05/2004
3. Variabili con upper bound
• E’ molto frequente che le variabili oltre ai lower
bound (xj ≥ 0) abbiano anche upper bound (xj ≤ hj)
• Il problema risultante
min
c Tx
Ax
x
x
=d
≥0
≤h
• Può essere trattato ponendo gli upper bound in forma
standard xj + xja = hj
• si aggiungono fino ad n vincoli e variabili slack!
D. Vigo
Pmat.Rev.9
Variabili con upper bound (2)
• Gli upper bound si possono trattare in modo
implicito (come i lower bound)
• Soluzione Base Ammissibile Estesa (SBAE):
le n–m variabili fuori base possono essere
• al lower bound (xj = 0)
• all’upper bound (xj = hj)
• Ai fini dell’ottimalità i ruoli delle variabili fuori
base sono diversi:
• se al l. b. può solo aumentare e conviene se c’j < 0
• se all’u. b. può solo diminuire e conviene se c’j > 0
D. Vigo
Pmat.Rev.10
5
19/05/2004
Criterio di ottimalità estesa
• Per un problema di minimo una SBAE è ottima se
il costo ridotto di ogni variabile fuori base è:
• c’j ≥ 0 se la variabile è al lower bound (xj = 0)
• c’j ≤ 0 se la variabile è all’upper bound (xj = hj)
• Il metodo del simplesso visto fin qui si estende
facilmente per trattare implicitamente le SBAE
D. Vigo
Pmat.Rev.11
Forma canonica estesa
• Supponiamo di definire due variabili non negative:
xj+ = xj (l.b.) e xj– = hj – xj (u.b.)
• In questo caso il test di ottimalità è sempre c’j ≥ 0
per entrambe
• In una data iterazione il simplesso usa una sola
delle due variabili xj+ = xj e xj– = hj – xj
• Si usa un vettore di flag ej ={+,–} per sapere quale
• Il tableau rappresenta quindi le equazioni:
n
∑y
j =1
D. Vigo
e
ij
x j j = yi 0
i = 1, K, m
Pmat.Rev.12
6
19/05/2004
Pivoting esteso
•
•
Determina una variabile xj fuori base con c’j < 0
determina θjE = min{αj, θj, γj } con:
•
•
•
αj = hj
θj = min { yi0/yij con yij > 0}
γj = min {(yi0 − hβ(i))/yij con yij < 0}
1. Se θjE = αj la variabile xj va al bound opposto
a) sottrai hj volte la colonna j dalla colonna 0
b) moltiplica la colonna j per −1 (cambia anche ej)
c) la base non cambia
D. Vigo
Pmat.Rev.13
Pivoting esteso (2)
2. Se θjE = θj = min { yi0/yij con yij > 0} (su riga k)
a) la variabile xβ(k) attualmente in base esce dalla base al
suo vecchio bound (eβ(k) non cambia)
b) esegui il pivot sull’elemento ykj
3. Se θjE = γj = min {(yi0 − hβ(i))/yij con yij < 0}
(su riga k)
a) la variabile xβ(k) attualmente in base esce dalla base al
bound opposto (eβ(k) cambia)
b) sottrai hβ(k) da yk0 , cambia il segno di ykβ(k) ed eβ(k)
c) esegui il pivot sull’elemento ykj
D. Vigo
Pmat.Rev.14
7
19/05/2004
Esempio
min 2x1
s.t. x1
+x2 +3x3 −2x4 +10x5
+x3 −x4 +2x5
x2 +2x5 +2x4 + x5
0 ≤ x1 ≤ 7, 0 ≤ x2 ≤ 10, 0 ≤ x3
0 ≤ x4 ≤ 5, 0 ≤ x5 ≤ 3
=5
=9
≤ 1,
x1
x2
x3
x4
x5
−z
0
2
1
3
−2
10
x1
5
1
0
1
−1
2
x2
9
0
1
2
2
1
+
+
+
+
+
e
D. Vigo
Pmat.Rev.15
Esempio (2)
x1
x2
x3
x4
x5
−z
−19
0
0
−1
−2
5
x1
5
1
0
1
−1
2
x2
9
0
1
2
2
1
+
+
+
+
+
e
pivot su colonna x4:
αj
θj
γj
=2
D. Vigo
= hj = 5
= min { yi0/yij con yij > 0} = 9/2
= min {(yi0 − hβ(i))/yij con yij < 0}
Í
Pmat.Rev.16
8
19/05/2004
Esempio (3)
Caso 3. Se θjE = min {(yi0 − hβ(i))/yij con yij < 0}
(su riga k =1)
a) la variabile xβ(1)= x1 esce dalla base al bound
opposto (eβ(k) cambia)
b) sottrai hβ(1)= 7 da y10 , cambia il segno di y11 ed e1
c) esegui il pivot sull’elemento y14
x1
x2
x3
x4
x5
−z
−19
0
0
−1
−2
5
x1
−2
−1
0
1
−1
2
x2
9
0
1
2
2
1
−
+
+
+
+
e
D. Vigo
Pmat.Rev.17
Esempio (4)
x1
x2
x3
x4
x5
−z
−15
2
0
−3
0
1
x4
2
1
0
−1
1
−2
x2
5
−2
1
4
0
5
−
+
+
+
+
e
pivot su colonna x3:
αj
θj
γj
=3
D. Vigo
= hj = 1
= min { yi0/yij con yij > 0} = 5/4
= min {(yi0 − hβ(i))/yij con yij < 0}
Í
Pmat.Rev.18
9
19/05/2004
Esempio (5)
Caso 1. Se θjE = αj , x3 va al bound opposto
a) sottrai h3=1 volte la colonna 3 dalla colonna 0
b) moltiplica la colonna 3 per −1 (cambia anche e3)
c) la base non cambia
x1
x2
x3
x4
x5
−z
−12
2
0
3
0
1
x4
3
1
0
1
1
−2
x2
1
−2
1
−4
0
5
−
+
−
+
+
e
Stop !
x1 = 7, x2 = 1, x3 = 1, x4 = 3, x5 = 0
D. Vigo
Pmat.Rev.19
Generazione di colonne
• In molti casi A e c di un problema PL o PLI sono
definiti implicitamente da una regola di
costruzione ed hanno un numero di elementi n >>1
• Il numero di elementi di A e c in base è m << n
• A e c servono solo per calcolare i costi ridotti ed
individuare le variabili da far entrare in base
• Se è possibile valutare i costi ridotti sulla base
della regola implicita non è necessario costruire
esplicitamente tutto il problema ma solo le
colonne che via via entrano in base
D. Vigo
Pmat.Rev.20
10
19/05/2004
Esempio: Bin Packing Problem
• m tipi di oggetti, per ogni tipo (i =1,…,m)
• wi dimensione di un oggetto di tipo i
• bi quantitativo richiesto
• n numero totale di oggetti = Σi=1,…,m bi
• n contenitori (bin), ciascuno di capacità K
• impaccare tutti gli oggetti nel minor numero
possibile di contenitori in modo che la somma
delle dimensioni degli oggetti inseriti in ogni
contenitore non superi la capacità K
• Modelli PLI con numero polinomiale di variabili
D. Vigo
Pmat.Rev.21
Modello di Gilmore–Gomory
• Primo esempio di generazione di colonne (1963)
• Pattern di caricamento: vettore di m elementi che
specifica il numero di oggetti di ciascun tipo
inseriti in un bin (in modo ammissibile)
• Es. m=5, w=(2,3,4,5,7), b=(20,10,30,15,20), K=10
• a1=(5,0,0,0,0), a2=(0,2,1,0,0), …
• aij = oggetti di tipo i nel pattern j
• J = insieme di tutti i pattern ammissibili
• |J| enorme (cresce esponenzialmente con m ed n)
D. Vigo
Pmat.Rev.22
11
19/05/2004
Modello di Gilmore–Gomory (2)
• xj = numero di volte in cui il pattern j è usato
(GG)
min
Σj∈J
Σj∈J
aij xj ≥ bi
xj
(i = 1, …, m)
≥ 0 intero
xj
(j∈J)
• Il vincolo potrebbe essere anche = ma non serve
• Consideriamo il rilassamento continuo (accettabile
se all’ottimo le x>0 sono abbastanza grandi)
D. Vigo
Pmat.Rev.23
Modello di Gilmore–Gomory (3)
• Supponiamo di conoscere J’⊆ J con |J’|<<|J| tale che
le colonne di J’ contengano una soluzione
ammissibile (= base) di (GG)
• Ad es. J’ ricavato da una soluzione euristica
• Risolviamo il problema rispetto a J’:
•
•
•
•
x’ soluzione primale (xj=0 j∈J\J’), u soluzione duale
la soluzione x’ è ottima se i costi ridotti c’j ≥ 0, j∈J
Pricing: trova una colonna da inserire in base (la migliore)
ossia trova un pattern ∉J’ con costo ridotto minimo
(
)
min c' j = min 1 − ua j = 1 − max ua j
j∈J
D. Vigo
j∈J
j∈J
Pmat.Rev.24
12
19/05/2004
Modello di Gilmore–Gomory (4)
• Pricing: trova un pattern ammissibile con il minimo
costo ridotto (se negativo va inserito in J’ altrimenti
stop: base ottima)
• ai = numero di oggetti di tipo i nel pattern
z = max
Σi=1,...,m
uiai
Σi=1,...,m wiai ≤ K
ai ≥ 0 intero (i=1,...,m)
• è un Knapsack intero !
• se z ≤ 0 la soluzione è ottima !
D. Vigo
Pmat.Rev.25
Esempio
•
•
•
•
•
m=5, w=(2,3,4,5,7), b=(20,10,30,15,20), K=10
soluzione/base iniziale: deve usare tutti i tipi
a1=(1,1,1,0,0), a2=(0,0,0,2,0), a3=(0,0,0,0,1)
mancano due colonne: slack di alcuni vincoli
Per soddisfare b3=30 bisogna usare a1 30 volte,
quindi 1 e 2 sono prodotti in eccesso e si possono
usare a1a=(–1,0,0,0,0) ed a2a=(0,–1,0,0,0)
D. Vigo
Pmat.Rev.26
13
19/05/2004
Esempio (2)
• Risoluzione rilassamento iniziale:
• si ricava la soluzione x=(30,15/2,20,0,0), che si arrotonda
in x=(30,8,20,0,0) che richiede 58 bin, la soluzione duale è
u=(0,0,1,1/2,1)
• Il pricing genera la colonna a4=(0,0,2,0,0) e z = 2
• Pivot sulla nuova colonna: base=(a1,a2,a3,a4,a2a)
• si ottiene x=(20,15/2,20,5,10), arrotondata in
x=(20,8,20,5,10) che richiede 53 bin ed u=(1/2,0,1/2,1/2,1)
D. Vigo
Pmat.Rev.27
Esempio (3)
• Il pricing genera la colonna a5=(5,0,0,0,0) e z = 5/2
• Pivot sulla nuova colonna: base=(a1,a2,a3,a4,a5)
• (nessuno slack sui vincoli)
• si ottiene x=(10,15/2,20,10,2), arrotondata in
x=(10,8,20,10,2) che richiede 50 bin ed
u=(1/5,3/10,1/2,1/2,1)
• Il pricing genera la colonna a6=(0,1,0,0,1) e z =
13/10
• Pivot sulla nuova colonna: base=(a6,a2,a3,a4,a5)
• si ottiene x=(10,15/2,10,15,4), arrotondata in
x=(10,8,10,15,4) che richiede 47 bin ed u=(1/5,0,1/2,1/2,1)
D. Vigo
Pmat.Rev.28
14
19/05/2004
Esempio (4)
• Il pricing genera la colonna a7=(1,0,0,0,1) e z = 6/5
• Pivot sulla nuova colonna: base=(a6,a2,a7,a4,a5)
• si ottiene x=(10,15/2,10,15,2), arrotondata in
x=(10,8,10,15,2) che richiede 45 bin ed
u=(1/5,1/5,1/2,1/2,4/5)
• Il pricing genera la colonna a8=(1,0,2,0,0) e z = 6/5
• Pivot sulla nuova colonna: base=(a6,a2,a3,a4,a8)
• si ottiene x=(10,15/2,10,5,10), arrotondata in
x=(10,8,10,5,10) che richiede 43 bin ed u=(0,0,1/2,1/2,1)
D. Vigo
Pmat.Rev.29
Esempio (5)
• Il pricing genera nuovamente la colonna
a3=(0,0,0,0,1) e z = 1
• La colonna fuori base a costo ridotto minimo ha costo
ridotto nullo: la base corrente è ottima !
• La soluzione ottima frazionaria è:
• base=(a6,a2,a3,a4,a8), x=(10,15/2,10,5,10), e usa 42,5 bin
• La soluzione intera arrotondata è:
• base=(a6,a2,a3,a4,a8), x=(10,8,10,5,10), e 43= ⎡42,5⎤ bin
Î soluzione ottima per il problema intero !
D. Vigo
Pmat.Rev.30
15