1. No category

Universidad Politécnica de Madrid–Escuela Técnica Superior de Ingenieros Industriales
Grado en Ingeniería en Tecnologías Industriales. Curso 2015-2016-3º
Matemáticas de Especialidad–Ingeniería Eléctrica
Programación Entera
José Luis de la Fuente O’Connor
[email protected]
[email protected]
Clase_Entera_2017.pdf
1/65
2/65
Índice
Introducción, formulación y ejemplos
Resolución gráfica del problema
Relajaciones en la formulación
El algoritmo de los cortes de Gomory
Algoritmos Branch and Bound
Programación no lineal en variables enteras
Introducción
3/65
La programación entera se ocupa de los problemas de optimizar una función de
diversas variables sujeta a condiciones de igualdad y/o desigualdad,
restringiéndose todas o alguna de esas variables a tomar valores enteros.
Sus orígenes se remontan a los años 50 del siglo XX. El pionero fue Ralph
Gomory, EE.UU. 1929.
Las áreas de aplicación práctica son muchas donde hay que asignar recursos sólo
disponibles en cantidades discretas: distribución de mercancías, programación de
la producción en factorías, secuenciación de maquinaria en procesos industriales
y productivos, asignación de grupos generadores de energía eléctrica, cadenas de
suministro, logística, programación de rutas de vuelos, etc.
4/65
La formulación general del problema lineal de Programación Entera es
T
T
maximizar
c
x
C
h
y
x 2Zn y 2Rp
sujeta a
Ax C G y b
x; y 0:
Este problema es un programa entero mixto.
Al conjunto S D fx 2 Zn; y 2 Rp W Ax C G y b; x 0; y 0g se le
denomina región factible.
Un punto Œx T ; y T T 2 S se denomina factible. Al punto Œx T ; y T T 2 S, tal
que
c T x C bT y c T x C hT y; para todo Œx T ; y T T 2 S ;
se le denomina solución óptima del problema.
Al problema en el que todas las variable son enteras,
5/65
max.n c T x
x 2Z
s. a Ax b
x 0;
se le conoce como programa entero puro, programa lineal en variables enteras o,
simplemente, programa entero. Un caso particular de programa entero es aquel
en el que las variables sólo pueden tomar valores 0 ó 1.
En programación entera no existe un método universal, como el simplex de
programación lineal, que obtiene la solución sirviéndose de las propiedades de
convexidad del problema a resolver.
Esto es así porque en programación entera la convexidad desaparece, no
pudiéndose utilizar por tanto la noción de gradiente para caracterizar y buscar el
óptimo de un problema, haciéndose necesario emplear métodos de resolución
específicos del aspecto combinatorio de las variables enteras.
In all these cases, the decision variable must be an integer, since we cannot stop 1.3 times or
build 2.7 distillation columns, or in the latter existence/non-existence example, the variable
will be a binary decision variable. Constrained optimisation problems where the decision
variables are integers are known as integer programming, or a binary linear program if the
variables are constrained to be either zero or one, or finally if only some are required to be
integers, then it is known as a mixed-integer linear programming or MILP problem.
6/65
Es más, los programas enteros,
como muestran los dos primeros ejemplos de la
Surprisingly binary or integer programming problems are much harder to solve than optimisation problems with real number solutions. The reason for this is that being discrete, we
figura, si se quiere minimizar
y, no tienen por qué tener el CHAPTER
óptimo
cerca
del
of possible 8.
solutions
406cannot use gradients or search directions, and therefore the number
OPTIMISATION
combinatorially. Explicit enumeration therefore is quickly out of the question. One
óptimo del programa congrows
variables
ni first
siquiera
óptimo
es un
may
be tempted incontinuas,
this instance to simply
solve the LPel
problem
in the normal
man- entero
ner, andvariables,
then subsequently
this optimum
nearest
integers.
decision
(circlesround
in Fig.
8.57), weresult
see tointhe
fact
there
is no However,
feasible as
region at all,
8.10. próximo
OPTIMISATION
405
a WITH
ese INTEGER
óptimoCONSTRAINTS
continuo
segunda
figura).
illustrated
in(la
Fig.
8.56
thatrelaxed
solutionlinear
may not
be reasonable,
or even
feasible. feasible region.
despite
the fact
that
the
problem
did have
a non-zero
Figure 8.55: A linear programming
problem where we restrict our attenFigure 8.57: Depending on the spaction to integer values of the decision
ing,
a perfectly well-formed linear provector x meaning that the allowable
gram
may
notfallhave
feasiblegrid.
integer
values
of x
on aa regular
solution
at
all.
Note
in
this
example,
Here the optimum of the linear prothere
the feagramare
andno
theintegers,
optimum◦,ofinside
the integer
sible
region.
linear
program are not the same.
Figure 8.56: Note that integer optimum is not simply the truncated or
rounded version of the real-valued
optimum, but could in fact be ‘arbitrarily’ far away.
be other
complications
as well. problems
Consider the
Fig. 8.57 which
ToThere
solvecan
integer
linear
programming
weexample
could shown
do an in
exhaustive
search over all
a perfectly
linear
optimisation
problem
a convex feasiblelarge,
region.soClearly
and worth to place in our small suitcase to maximise
value,
(the valid
’knapsack’
problem),
and
theshows
feasible
integers.
For
large
problems,
this iswith
combinatorially
a better solution
somewhere
in two
that storage
will
be an
optimum, irrespective of what our particular
whether we should or should not install a newstrategy,
pump
between
tanks,
(you
can’t
in fact
oneregion
of thethere
only
known
practical
solution strategies, is to use what is known
objective function happens to be. However when we restrict our attention to just integer
install half a pump!).
as a branch and bound technique, [45, p411–415]. Here we start by solving a relaxed ver-
También puede que la región
factible no registre ni un sólo entero factible,
sion of the original problem, namely we simply ignore the integer constraints, and solve
existiendo dicho espacio en el sentido tradicional.
In all these cases, the decision variable must be an integer, since we cannot stop 1.3 times or
the linear program asexample,
usual. The
term ‘relaxed’ means that we have deliberately relaxed
build 2.7 distillation columns, or in the latter existence/non-existence
the variable
some
of
the
constraints
in
order
to
achieve a solution, albeit perhaps an infeasible one. If
will be a binary decision variable. Constrained optimisation problems where the decision
we happen or
fortuitously
to stumble
variables are integers are known as integer programming,
a binary linear
program across
if the an integer solution, we are done. However in
mostorcases,
will some
not get
integertosolution,
but we do have an under-estimator of the
variables are constrained to be either zero or one,
finallywe
if only
arean
required
be
minimum
cost function.
In the case of Fig. 8.56, we can see that the real-valued optimum
integers, then it is known as a mixed-integer linear
programming
or MILP problem.
7/65
Ejemplos
Gestión de un servicio hospitalario
En un determinado servicio de asistencia hospitalaria hay i enfermos a la espera
de una operación. Cada enfermo necesita una operación de duración Di .
Dada una disponibilidad de cirujanos, la suma de las duraciones de las
operaciones que pueden hacerse cada día j del período de estudio es igual a Tj .
Se trata de minimizar una función económica que refleja la suma de las
penalizaciones por esperas de los diferentes enfermos (esta penalización es una
función lineal creciente de la espera).
8/65
El problema se formula de la siguiente manera:
XX
max.
cij xij
i
s. a
j
X
i
Di xij Tj
X
xij D 1
para todo j
para todo i
j
xij D 0 ó 1
para todo i y j .
La variable de decisión xij es 1 si el enfermo i se le opera el día j y 0 en el
caso contrario.
9/65
El problema de la mochila
Este problema es el más clásico de la programación entera.
Trata de un transportista que dispone de un vehículo con capacidad volumétrica
de transporte b, en el que se pueden alojar diferentes objetos j de volúmenes
aj y valores económicos cj , y que quiere maximizar el valor de lo transportado
en cada viaje.
Su formulación es:
max.
X
s. a
X
cj xj
j
aj xj b
j
xj D 0 ó 1:
Resolución gráfica del problema
10/65
Recordemos de programación lineal que las condiciones de
T
T
maximizar
c
x
C
h
y
x 2Zn y 2Rp
sujeta a Ax C G y b
x; y 0
en el caso en que x 2 Rn y y 2 Rp , definen un politopo (poliedro si es
acotado) en uno de cuyos puntos extremos la función objetivo c T x C hT y
alcanza el máximo.
En general, desafortunadamente, las variables xj es ese punto extremo no
tienen por qué ser enteras.
El politopo contiene todas las soluciones en las que x es entero. Gráficamente
se puede determinar la solución mediante el punto extremo óptimo considerando
las variables como continuas y luego, desplazando c T x C hT y hacia dentro del
politopo, encontrando el primer punto Œx T ; y T T en el que x sea entero.
11/65
Ejemplo Consideremos el problema de programación entera
x 1 + x 2 ≤ 10
x2
−x 1 + 4x 2 ≤ 16
B
max. 23 x1 C x2
s. a
x1 C 4x2 16
x1 C x2 10
5x1 C x2 36
x1 ; x2 0
x1 y x2 enteras.
A
x∗
C
5x 1 + x 2 ≤ 36
f. O.
D
0
x1
Sin tener en cuenta la condición de que las variables han de ser enteras, la
región factible es el poliedro convexo que generan los vértices 0; A; B; C y D
de la figura.
El valor máximo de la función objetivo en este poliedro se alcanza en C .
Como la auténtica región factible la constituyen los puntos interiores enteros de
u
ese poliedro, el máximo del problema se alcanza en el punto x D Œ6; 4T .
12/65
Índice
Introducción, formulación y ejemplos
Resolución gráfica del problema
Relajaciones en la formulación
El algoritmo de los cortes de Gomory
Algoritmos Branch and Bound
Programación no lineal en variables enteras
13/65
Relajaciones en la formulación
La idea es tratar de sustituir el programa entero original por otro más fácil de
resolver en el que se relajan algunas de las condiciones.
La solución de éste proporciona una cota de la función objetivo del problema
original.
Mediante un proceso iterativo que integre estas relajaciones y la resolución de
los subproblemas se puede acotar cada vez más la función objetivo hasta
obtener la solución final.
14/65
Una relajación del programa entero
˚
max. c T x W x 2 S
.PE/
la constituye cualquier subproblema de la forma
max. fzPR .x/ W x 2 SPR g
.PR/
con las propiedades siguientes:
SPR S
y c T x zPR .x/ para x 2 S:
Proposición Si el problema (PR) no es factible, tampoco lo es (PE). Si (PR) es
factible, los valores de las funciones objetivo cumplen que zPE zPR .
15/65
Relajación lineal
La relajación más inmediata resulta de omitir la condición de que las variables
han de ser enteras.
Si el conjunto o región factible del programa entero es S D P \ Zn, la
relajación lineal es
.PL/
max. c x W x 2 P ; P D fx 2 Rn W Ax b; x 0g
˚
T
Proposición Si (PL) es no acotado, el programa entero (PE) es no factible o no
acotado.
16/65
Algoritmo tipo
.1/
Paso 0 – Inicialización Hacer i D 1, w D 1 y z D 1. Escoger un SR S tal que
.1/
zR c T x para x 2 S .
Paso 1 – Relajación Resolver el problema PE relajado:
n
o
.i /
.i /
.i /
R
max. zR .x/ W x 2 SR :
Paso 2 – Comprobación de óptimo Sea x .i / la solución del problema anterior. Si x .i / 2 S:
PARAR ; ésta es la solución óptima de .PE/ con función objetivo w D c T x .i / D
z ; si no, seguir.
.i /
Paso 3 – Mejora de la solución Hacer w D zR , z D c T x .i / si x .i / 2 S y i
i C 1.
.i C1/
.i C1/
.i /
.i C1/
.iC1/
Escoger SR
tal que S SR
SR y zR .x/ tal que c T x zR .x/ .i /
.i C1/
.i /
.i C1/
.i /
zR .x/ para x 2 S con SR
¤ SR o zR .x/ ¤ zR .x/.
Ir al paso 1
17/65
Índice
Introducción, formulación y ejemplos
Resolución gráfica del problema
Relajaciones en la formulación
El algoritmo de los cortes de Gomory
Algoritmos Branch and Bound
Programación no lineal en variables enteras
Algoritmo de planos o cortes de Gomory
Este algoritmo, formulado por Gomory en 1960, es uno de los primeros que se
emplearon para resolver el problema general de programación entera.
Consideraremos que el problema que resuelve –notación original– es
˚
max. c T x W x 2 S .e/ ; donde S .e/ D fx 2 Zn W Ax D b; x 0g;
y está escrito de la siguiente forma:
n
o
T T
.0/
max. x0 W x0; x
2S
;
.PE/
donde
S
.0/
˚
D x0 2 Z; x 2 Zn W x0
T
c x D 0; Ax D b; x 0 :
Supondremos que se dispone de una solución y base óptimas de la relajación
lineal del programa entero.
18/65
19/65
El problema se puede escribir así:
maximizar x0
xBi C
s. a
X
aN ij xj D aN i 0
para i D 0; 1; : : : ; m
(1)
j 2H
xB0 2 Z; xBi 2 ZC ;
xj 2 ZC para j 2 H;
x0 D xB0 , xBi ; i D 1; : : : ; m,
no básicas (xB D B
1b
B
para i D 1; : : : ; m
son las variables básicas y xj ,
1N x
N
y z D cB B
1b
C .cN
cB B
j 2 H N D f1; : : : ; ng,
1 N /x
N.
las
).
Como esta base es factible en el programa primal y dual se cumple que aN i0 0
para i D 1; : : : ; m, y que aN 0j 0, para j 2 H .
Suponiendo existe un i tal que aN i0 … Z, se tiene el siguiente resultado.
Proposición
P
Corte Fraccionario de Gomory Si aN i 0 … Z, j 2H fij xj D fi 0 C xnC1 ; xnC1 2 ZC ,
es una desigualdad válida en S .0/ , donde fij D aN ij baN ij c para j 2 H y fi 0 D aN i 0 baN i 0 c.
20/65
Ejemplo Considérese el programa entero
maximizar 7x1 C 2x2
.PE/
s. a
x1 C 2x2 4
5x1 C x2 20
2x1 2x2 7
x 2 Z2C:
Introduzcamos las variables de holgura x3, x4 y x5. El programa1 queda:
maximizar 7x1 C 2x2
s. a
1
x1 C 2x2 C x3
D
5x1 C x2
C x4
D
2x1 2x2
C x5 D
x1; x2 2 ZCI x3; x4; x5 4
20
7
0:
Aunque las variables de holgura serán enteras, no hay ninguna necesidad de restringirlas a tomar valores
enteros.
Restringiendo las variables x3, x4 y x5 a tomar valores enteros el problema es:
maximizar x0
7x1 2x2
D 0
x1 C 2x2 C x3
D 4
5x1 C x2
C x4
D 20
2x1 2x2
C x5 D 7
x0 2 Z; xj 2 ZC para j D 1; : : : ; 5:
x0
s. a
La solución de la relajación lineal de este programa entero es
C
x0
x1
x2 C
3
x
11 3
1
x
11 3
5
x
11 3
8
x
11 3
C
C
C
C
16
x
11 4
2
x
11 4
1
x
11 4
6
x
11 4
D
D
D
C x5 D
332
11
36
11
40
11
75
:
11
El corte de Gomory que se puede obtener de la primera ecuación
3
5
2
x3 C x4 D
C x6 ;
11
11
11
x6 2 ZC:
21/65
22/65
Paso 0 – Inicialización Hacer
˚
SR.1/ D x0 2 R; x 2 Rn W x0
c T x D 0; Ax D b; x 0 ;
.1/
i D 1 y zR
D x0 .
Paso 1 – Resolución de la relajación lineal Resolver el programa entero relajado
n
o
max. x0 W Œx0 ; x T T 2 SR.i / :
.PR.i/ /
T
Si PR.i / es factible y tiene solución óptima, Œx0.i / ; x .i / T , continuar.
Paso 2 – Comprobación de óptimo Si x .i / 2 ZnC , ésta es la solución óptima del programa entero.
Paso 3 – Comprobación de no factibilidad Si PR.i / no es factible, el programa entero original no es factible.
P
.i /
Paso 4 – Adición de un corte de Gomory Escoger una fila xBi C j 2H .i / aN ij
xj D aN i.i0/ con aN i.i0/ … Z. Hacer
X
fkj xj
xnCi D fk0 ;
k D 1; : : : ; m; xnCi 2 ZC ;
j 2H .i/
el corte de Gomory de esa fila. Hacer
8
<
SR.i C1/ D SR.i / \ x0 2 R; x 2 RnCi W
:
Hacer i
X
j 2H .i /
fkj xj
xnCi
9
=
D fk0 ; x 0 :
;
i C 1 e ir al paso 1.
Algoritmo de los planos o cortes de Gomory
23/65
Cuando se añade un corte, la nueva base, que incluye xnCi como variable
básica, es factible del dual.
La factibilidad del primal se incumple sólo en que xnCi es negativa.
Para reoptimizar el problema lo más adecuado es usar el método dual del
simplex.
Si PR.i / es no acotado, el programa entero, PE, es no acotado o no factible.
Ejemplo, continuación. La última solución de PR.1/ a la que hacíamos
referencia era
16
332
3
C
x0
x1
x2 C
x
11 3
C
x
11 4
D
11
1
x
11 3
C
2
x
11 4
D
36
11
5
x
11 3
C
1
x
11 4
D
40
11
8
x
11 3
C
6
x
11 4
C x5 D
75
11
;
donde x3hD x4 D 0. i
T
.1/
.1/
.1/
.1/
D Œx0 ; x1 ; : : : ; x5  D 332
Es decir x0 ; x .1/
;
11
36 40
; ;
11 11
0; 0;
El corte de Gomory respecto de la tercera ecuación es 115 x3 C 111 x4 D
x6 2 ZC. Añadiéndolo a PR.1/ se obtiene PR.2/ cuya solución es
C 75 x4
x0
C
C 15 x4
x1
C
x2
x3 C
2
x
5 4
1
x
5 4
C x5 C
3
x
5 6
1
x
5 6
D
D
149
5
17
5
x6 D 3
8
x
5 6
11
x
5 6
D
D
29
5
7
5
:
75
11
7
11
.
C x6 ,
24/65
25/65
El corte de Gomory correspondiente a la primera de estas ecuaciones es
2
x C 35 x6 D 45 C x7 , x7 2 ZC . Añadiéndolo a PR.2/ se obtiene PR.3/ , de solución:
5 4
C x4
x0
C 13 x4
x1
x2
x3
C
2
x
3 4
5
x
3 4
2
x
3 4
2
x
3 4
C
C x5
El siguiente corte se elige de la fila x2
C
C x6
2
x
3 4
2
2
1
x4 C x7 D C x8;
3
3
3
x7 D 29
1
x
3 7
5
x
3 7
11
x
3 7
8
x
3 7
5
x
3 7
D
D
D
D
D
11
3
5
3
13
3
11
3
4
:
3
C 53 x7 D 53 . Es2
x8 2 ZC:
La solución óptima de PR.4/ –también del problema original–, es
h
i
.4/
.4/
.4/
.4/
.4/
x0 ; x1 ; : : : ; x6 ; x7 ; x8 D Œ28; 4; 0; 8; 0; 1; 3; 1; 0:
2
Recordemos
2
3
˘
D
2
3
.
3
/
3
D 13 .
26/65
En términos de las variables originales, los tres cortes añadidos al problema son
x2 3, 2x1 C x2 9 y 3x1 C x2 12.
x2
x (1)
x (2)
3
Corte 1
2
x (3)
Corte 2
Corte 3
1
x (4)
1
2
3
4
x1
27/65
Índice
Introducción, formulación y ejemplos
Resolución gráfica del problema
Relajaciones en la formulación
El algoritmo de los cortes de Gomory
Algoritmos Branch and Bound
Programación no lineal en variables enteras
28/65
Algoritmos de ramificación y acotamiento
o Branch and Bound
La idea que los anima es dividir la región factible, S, como si las variables que la
definen fuesen las ramas de un árbol e ir acotando la solución del estudio de
esas ramas.
El conjunto fS .i/ W i D 1; : : : ; kg se denomina división de la región factible, S,
de un programa entero, si [kiD1S .i / D S.
Una división se denomina partición si S .i/ \ S .j / D ; para i; j D 1; : : : ; k,
i ¤ j.
29/65
La división se suele hacer de forma recursiva así:
S = S (0)
S (1)
S (11)
S (12)
S (121)
S (2)
S (13)
S (21)
S (122)
S (22)
30/65
Cuando S Zn y las variables han de ser 0 ó 1, la forma más simple de hacer
la división recursiva apuntada es la indicada en la figura.
S
x1 = 0
x1 = 1
S (0)
x2 = 0
S (1)
x2 = 1
S (00)
x3 = 0
S (000)
x2 = 0
S (01)
x3 = 1
S (001)
x3 = 0
S (010)
x3 = 1
S (011)
x2 = 1
S (11)
S (10)
x3 = 0
S (100)
x3 = 1
S (101)
x3 = 0
S (110)
x3 = 1
S (111)
En la práctica estos métodos necesitan algún procedimiento para eliminar ramas
del árbol de búsqueda que se forma –podarlo–.
31/65
Proposición El árbol de búsqueda se puede podar a partir del nudo correspondiente
a S .i/ si se cumplen cualquiera de las siguientes condiciones:
1. No factibilidad: S .i / D ;.
2. Optimalidad: se ha llegado al óptimo de PE .i / .
.i /
3. Existencia de un valor dominante zPE zPE .
Para no tener que resolver PE .i/ se resuelve su relajación lineal; esto es PL.i/
.i /
.i/
con S .i / SPL y zPL.x/ c T x para x 2 S .i/.
Proposición El árbol de búsqueda puede podarse a partir del nudo correspondiente a
S .i/ si se cumple cualquiera de las condiciones siguientes:
1. El programa PL.i / no es factible.
.i /
.i /
.i /
.i /
2. El óptimo xPL de PL.i / cumple que xPL 2 S .i / y que zPL D c T xPL .
.i /
3. zPL zPE , donde zPE es alguna solución factible de PE.
Ejemplo Considérese el programa entero
32/65
S
x1 = 0
.ZPE /
maximizar
100x1 C 72x2 C 36x3
s. a
2x1 C x2
4x1
C
x1 C x2 C
x 2 Z3 .0; 1/:
0
x3 0
x3 0
x1 = 1
S (0)
S (1)
x2 = 0
x2 = 1
S (10)
S (11)
x3 = 0
x3 = 1
S (110)
La división de este ejemplo y su árbol de búsqueda son los de la figura.
La solución de la relajación lineal es Œ1=2; 1; 1T ; la f.o. 58.
Usando la última proposición para podar el árbol de búsqueda,
S .0/ D fx 0 W x1 D 0; x2 0; x3 0; x2 C x3 1g D ;;
por lo que el nudo S .0/ se desecha.
S (111)
33/65
La condición de óptimo de los nudos S .110/ y S .111/ se cumple pues estos
conjuntos contienen las soluciones Œ1; 1; 0T y Œ1; 1; 1T , respectivamente.
Como, zPE < zPE D 8, se tendrá que zPE D zPE D 8.
Aplicando la tercera condición de la mencionada proposición al nudo S .10/ se ve
.10/
que zPL D 64 < zPE , por lo que se desecha.
Como resultado final se obtiene que x D Œ1; 1; 1T es la solución óptima del
programa entero.
El valor de su función objetivo es ZPE D 8.
.110/
.111/
.111/
u
34/65
Paso 0 – Inicialización Hacer L D fPEg, S .0/ D S , z .0/ D 1 y zPE D 1.
Paso 1 – Comprobación de final Si L D ;, la solución x que daba la función objetivo zPE D c T x es la
óptima.
Paso 2 – Selección de problema, relajación lineal y resolución Seleccionar y borrar de L un problema PE .i / .
.i /
Resolver su relajación lineal PL.i/ . Sea zPL el valor óptimo de la función objetivo de este problema
.i/
y xPL su solución óptima (si existe).
.i/
Paso 3 – Poda Si zPL zPE ir al paso 1 (nótese que si la relajación se resuelve por un algoritmo dual, este
paso es aplicable tan pronto como el valor de la función objetivo del programa dual se hace
inferior a zPE ).
.i/
Si xPL … S .i/ ir al paso 4.
.i/
.i/
.i /
Si xPL 2 S .i/ y c T xPL > zPE , hacer zPE D c T xPL . Borrar de L todos los problemas tales
.i/
.i /
que z .i/ zPE . Si c T xPL D zPL ir al paso 1; si no, al paso 4.
Paso 4 – División Sea fS .ij / gjkD1 una división de S . Añadir los problemas fPE .ij / gjkD1 a L, donde z .ij / D
.i/
zPL para j D 1; : : : ; k .
Ir al paso 1.
Algoritmo Branch and Bound
L designa una colección de programas enteros, PE .i/, cada uno tiene la forma
.i /
zPE D max. fc T x W x 2 S .i/g, donde S .i/ S .
35/65
Branch-and-bound con relajación lineal
Consideraremos los programas enteros generales
˚
.PE/
max. c T x W x 2 S ; con S D fx 2 Zn W Ax b; x 0g ;
y la forma de resolverlos
Mediante esta relajación, el conjunto factible, S, se reemplaza por
.0/
SPL D fx 2 Rn W x 0; Ax bg,
haciéndose en cada relajación zPL.x/
D c T x.
a i = 1, 2, zP L = max. {cT x : x ∈ SP L } < zP L .
EnDivisión
la práctica sólo
usan vectores [dT , d0 ]T muy concretos:
de seramas
36/65
• Dicotomı́as de variables. En este caso d = ej para algún j ∈ N . El punto x(0) será n
(0)
factible
en
las
relajaciones
resultantes
si
x
∈
/ Z y d0 = xj (ver figura 13.5). Nótes
j
Como se usa relajación lineal, la división se debe hacer
con condiciones lineales.
que si xj ∈ B la rama izquierda hace xj = 0 y la derecha xj = 1.
La forma más conveniente de hacerlo es como indica la figura.
xj ≤ xj xj ≥ xj + 1
Es decir, desde el nudo que define
la solución
Figura
13.5 correspondiente, si esta no tiene
algunoDivisión,
de sus elementos,
j , que
havariable
de ser entero,
entero,
hace partir dos
por dicotomı́a
de la
xj , en un
árbol se
enumerativo
ramas: una en la que ese xj se acote superiormente por su entero inmediato
Una
importante
ventaja
práctica
de esta división
es que superior.
las condiciones que se añade
inferior;
otra que
se acote
inferiormente
por su entero
a la relajación lineal correspondiente son simples cotas inferiores o superiores. Sólo ser
Usaremos indistintamente estos términos.
Los dos nuevos programas se pueden resolver mediante el método dual del
simplex; el tamaño de la base seguirá igual.
37/65
Selección del nudo a estudiar
Proposición Si P D fx 2 Rn W Ax b; x 0g es acotado, el árbol de búsqueda
desarrollado a partir de dicotomías de las variables es finito. En particular, si
P!j D dmKaxfxj W
x 2 P ge, ningún camino del árbol de búsqueda puede contener más de j 2N !j ramas.
Proposición Si el nudo t del árbol de búsqueda con conjunto de condiciones S .t / es
.t /
tal que max. fc T x W x 2 SPL g > zPE , ese nudo no se puede podar o rechazar.
Dada una lista L de subárboles del árbol de búsqueda con nudos no rechazados,
qué nudo elegir como el próximo a examinar.
Dos opciones:
8
ˆ
< Analizar uno siguiendo unas reglas preestablecidas; o
ˆ
: Escogerlo según la información disponible en ese punto sobre
cotas de las variables, número de nudos activos, etc.
Existen varias reglas preestablecidas que se pueden utilizar. La más extendida es
la denominada LIFO (del inglés Last In First Out).
LIFO: Primero analizar todos los nudos “hijos"de uno determinado y luego
volver hacia atrás y escoger otro del mismo nivel del nudo “padre” que no haya
sido analizado.
38/65
1
2
3
5
4
6
7
9
8
Los nudos subrayados son los podados o rechazados.
Ventajas:
8 La relajación lineal de un nudo hijo se obtiene de la lineal del nudo padre sin más que
ˆ
ˆ
ˆ
< añadir una cota en una variable. El método dual del simplex, sin reinvertir la base ni
modificar las estructuras de datos sustancialmente, obtiene muy rápidamente la solución
que define el nudo hijo.
ˆ
ˆ
ˆ
: La experiencia dice que las soluciones factibles se encuentran más bien en puntos profundos del árbol de búsqueda.
Otra regla: Búsqueda en amplitud o anchura. Analizar todos los nudos de un
determinado nivel antes de pasar a otro más profundo. En códigos
especializados en variables binarias.
39/65
Un código “sencillo” recursivo de Matlab para resolver un Programa Entero:
function [xstar,Jmin] = ilp_bb_1(c,A,b,xlb,xub,xstar,Jmin)
% Integer linear programming by (recursive) branch & bound
% This algorithm uses linprog.m from the optimisation toolbox
tol = 1e-5;
% anything less is considered zero
if ~exist(’xlb’,’var’), xlb = 0*c; end % canonical form (for the moment)
if ~exist(’xub’,’var’), xub = Inf*ones(size(c)); end
% Upper bound defaults
if ~exist(’xstar’,’var’), xstar = []; end
if ~exist(’Jmin’,’var’), Jmin = Inf; end
% Try to minimise this
% Use the following two line if you prefer OPTI Optimization toolbox
% LPrelax = opti(’f’,c,’ineq’,A,b, ’bounds’, xlb, xub);
% [x,j,exitflag] = solve(LPrelax);
optns = optimset(’display’,’off’);
% turn off diagnostics in Linprog
[x,j,exitflag] = linprog(c,A,b,[],[],xlb,xub,[],optns);
if exitflag~=1, return, end
% if infeasible, branch ended
if j>Jmin, return, end
% if current cost J=cT x is worse, drop
idx = find(abs(x-round(x)) > tol);
% which are non-integers ?
if isempty(idx)
% All integer solutions
if j<Jmin, xstar=round(x); Jmin=j; end % New integer optimum
return
% problem solved
end
% Depth first recursion
xnlb = xlb; xnlb(idx(1)) = ceil(x(idx(1)));
% Take x<--int. greater than x
[xstar,Jmin] = ilp_bb_1(c,A,b,xnlb,xub,xstar,Jmin);
xnub = xub; xnub(idx(1)) = floor(x(idx(1)));
% Take x<--int. less than x
[xstar,Jmin] = ilp_bb_1(c,A,b,xlb,xnub,xstar,Jmin);
end
40/65
Ejemplo Resolvamos este problema de programación entera
minimizar
s. a
x1 C 10x2
66x1 C 14x2 1430
82x1 C 28x2 1306
x1 ; x2 enteras; 0:
Con el programa presentado
>> c3=[1;10]
>> A3=[-66 -14;82 -28]
>> b3=[-1430;-1306]
>> [xreal,jreal]=linprog(c3,A3,b3)
Optimization terminated.
xreal =
7.261682243025319
67.909212283785493
jreal =
6.863538050808803e+002
>> [xint,j]=ilp_bb_1(c3,A3,b3)
xint =
7
70
j =
707
41/65
Con la herramienta OPTI Toolbox
>> f1=[1 10]’;
>> A1=[-66 -14; 82 -28];
>> b1=[-1430 -1306]’;
>> lb = [0;0]; ub = [];
>> xtype = ’II’;
>> Opt=opti(’f’,f1,’ineq’,A1,b1,’lb’,lb,’xtype’,xtype)
-----------------------------------------------------Mixed Integer Linear Program (MILP) Optimization
min f’x
s.t. rl <= Ax <= ru
lb <= x <= ub
xi = Integer / Binary
-----------------------------------------------------Problem Properties:
# Decision Variables:
2
# Constraints:
6
# Linear Inequality:
2
# Bounds:
2
# Integer Variables:
2
-----------------------------------------------------Solver Parameters:
Solver:
CBC
----------------------------------------------------->> [xi,fval]=solve(Opt)
xi =
7
70
fval =
707
>> plot(Opt)
75
MILP Optimization Contours + Constraints - Min: [7; 70] Fval: 707
760
750
750
74
740
740
73
730
730
72
720
720
71
x2
710
710
70
700
700
69
690
690
68
680
680
67
670
670
66
660
65
2
3
4
5
6
7
x1
8
660
9
10
11
12
42/65
Resolvámoslo ahora con el software BBMI de programación entera y lineal. El
fichero de datos que requeriría este problema es el que sigue
PROCESO ITERATIVO
NAME
Ej12.2.3
ROWS
N OBJ
G FILA1
G FILA2
COLUMNS
INT
X1
OBJ
X1
FILA2
X2
OBJ
X2
FILA2
RHS
RHS
FILA1
RHS
FILA2
ENDATA
1
-82
10
28
1430
1306
FILA1
66
FILA1
14
665
Ejercicios
La gráfica de cómo opera BBMI es esta:
43/65
66x1 + 14x2 ≥ 1430
x2 = 74
−82x1 + 28x2 ≥ 1306
x2 = 73
c
x2 = 72
4
x2 = 71
x2 ≥ 68
∗
8
x2 = 70
1
x2 ≤ 67
3
2
11
x2 = 69
x2 = 68
6=7
2
1
x2 ≥ 71
x1 ≥ 8
x1 ≤ 7
3
6
x2 ≤ 70
x1 ≤ 6
4
x2 = 67
x1 = 7
x1 = 8
x1 = 9
x2 ≥ 70
7
x2 ≤ 69
8
5
x1 = 6
x1 = 7
10
9
El resultado impreso de BBMI se lista a continuación.
Problema Ej12.2.3
*** Estadísticas del Problema
3 Fila(s)
2 Restricción(es)
2 Variable(s) de decisión
2 Variable(s) de holgura/artificiales
6 Elementos no cero
2 Variable(s) entera(s)
Densidad de la matriz de coeficientes A: 100.000%
*** Estadísticas de INVERT
3 elementos no cero en la base
0 columnas estructurales en la base
0 vectores columna antes del "bump"
3 vectores columna después del "bump"
L:
0 elementos no cero;
0 vectores ETA
U:
2 elementos no cero;
2 vectores ETA
Total:
0 elementos no en la diagonal;
2 vectores ETA
Máximo de transformaciones ETA:
2; máximo número de elementos ETA:
Error relativo en x:
.000000D+00
Ite. Inf. Nopt. Sum.Inf/F.Obj Entra de-a Cost.Red Sale de-a
Paso
Pivote El.Eta
1
2
1
.2736000D+04
2 L->B -.42D+02
3 B->L
2
0
1
.6863538D+03
1 L->B -.46D+02
4 B->L
9
.10D+03 -.14D+02
.73D+01 .21D+03
2
5
--- SOLUCIÓN INICIAL PROGRAMA LINEAL ---------------------------------Nombre del problema: Ej12.2.3
No. de iteraciones:
3
Valor de la función objetivo:
686.353805073431
*** FILAS
No. ..Fila.. en ....Valor.... ...Holgura... .Lí.Inferior. .Lí.Superior. Val.Dual.
1 OBJ
2 FILA1
3 FILA2
BS 686.35381
LI 1430.0000
LI 1306.0000
-686.35381
.00000000
.00000000
Ninguno
1430.0000
1306.0000
Ninguno
Ninguno
Ninguno
1.000
.2830
.2156
*** COLUMNAS
No. .Columna en ....Valor.... Coste en F.O. .Lí.Inferior. .Lí.Superior. Cos.Red.
1 X1
2 X2
BS 7.2616822
BS 67.909212
1.0000000
10.000000
.00000000
.00000000
Ninguno
Ninguno
.000
.000
44/65
45/65
Variable
Separación
Iteración
2
1
2
Nivel
Dirección
Func. Obj.
1
X> 68
2
X>
8
3
X> 71
Nudos
en Lista
1
2
3
* Nueva solución entera; z(PE)=
Variables
Ent. No Ent.
2
1
1
Valor
4 D
6 D
7 D
718.00000; Tiempo desde última:
6.8724242D+02
7.0871429D+02
7.1800000D+02
.0001 seg.
--- SOLUCIÓN ENTERA ----------------Tiempo desde última:
.0001 seg.
Nombre del problema: Ej12.2.3
No. de iteraciones:
7
Valor de la función objetivo:
718.000000000000
*** FILAS
No. ..Fila.. en ....Valor.... ...Holgura... .Lí.Inferior. .Lí.Superior. Val.Dual.
1 OBJ
2 FILA1
3 FILA2
BS 718.00000
BS 1522.0000
BS 1332.0000
-718.00000
92.000000
26.000000
Ninguno
1430.0000
1306.0000
Ninguno
Ninguno
Ninguno
1.000
.1269E-15
.1349E-16
*** COLUMNAS
No. .Columna en ....Valor.... Coste en F.O. .Lí.Inferior. .Lí.Superior. Cos.Red.
1 X1
2 X2
Variable
Separación
Iteración
2
1
2FIJ
LIB 8.0000000
LIB 71.000000
1.0000000
10.000000
Nivel
Dirección
Func. Obj.
3
X< 70
2
X<
7
1
3
X> 70
Nudos
en Lista
2
1
3
* Nueva solución entera; z(PE)=
8.0000000
71.000000
Variables
Ent. No Ent.
707.00000; Tiempo desde última:
Tiempo desde última:
.0001 seg.
Nombre del problema: Ej12.2.3
No. de iteraciones:
10
Valor de la función objetivo:
1.00
10.0
Valor
-Nudo desechado en BKTRAK- 1.0000000D+20
0
9 P
6.9842857D+02
1
10 D
7.0700000D+02
--- SOLUCIÓN ENTERA -----------------
*** FILAS
Ninguno
Ninguno
707.000000000000
.0001 seg.
No. ..Fila.. en ....Valor.... ...Holgura... .Lí.Inferior. .Lí.Superior. Val.Dual.
1 OBJ
2 FILA1
3 FILA2
BS 707.00000
BS 1442.0000
BS 1386.0000
-707.00000
12.000000
80.000000
Ninguno
1430.0000
1306.0000
Ninguno
Ninguno
Ninguno
1.000
.1269E-15
.1349E-16
*** COLUMNAS
No. .Columna en ....Valor.... Coste en F.O. .Lí.Inferior. .Lí.Superior. Cos.Red.
1 X1
2 X2
Variable
Separación
Iteración
2
1
2
EQB 7.0000000
LIB 70.000000
1.0000000
10.000000
Nivel
Dirección
Func. Obj.
3
X< 69
3
X<
6
1
X< 67
7.0000000
70.000000
Nudos
en Lista
2
1
0
Ninguno
Ninguno
Variables
Ent. No Ent.
1.00
10.0
Valor
-Nudo desechado en BKTRAK- 1.0000000D+20
-Nudo desechado en BKTRAK- 7.4457143D+02
-Nudo desechado en BKTRAK- 1.0000000D+20
--- SOLUCIÓN ENTERA OPTIMA -----------------------Nombre del problema: Ej12.2.3
No. de iteraciones:
10
Valor de la función objetivo:
707.000000000000
*** FILAS
No. ..Fila.. en ....Valor.... ...Holgura... .Lí.Inferior. .Lí.Superior. Val.Dual.
1 OBJ
2 FILA1
3 FILA2
BS 707.00000
BS 1442.0000
BS 1386.0000
-707.00000
12.000000
80.000000
Ninguno
1430.0000
1306.0000
Ninguno
Ninguno
Ninguno
1.000
.1269E-15
.1349E-16
*** COLUMNAS
No. .Columna en ....Valor.... Coste en F.O. .Lí.Inferior. .Lí.Superior. Cos.Red.
1 X1
2 X2
EQB 7.0000000
LIB 70.000000
1.0000000
10.000000
Tiempo total de CPU en cálculos:
7.0000000
70.000000
.0603 segundos
Ninguno
Ninguno
.000
10.0
46/65
47/65
Para darnos una idea de las prestaciones de estos algoritmos, resolvamos un
problema relativamente difícil de Programación Entera con 50 variables, todas
enteras, y cuya matriz de condiciones tiene este esquema:
0
5
10
15
20
25
30
35
40
0
10
20
30
40
50
nz = 221
Lo vamos a resolver con el solver de PE de OPTI Toolbox y con ilp_bb_1
Con OPTI Toolbox:
>> prob_1 = coinRead(’DEMOS5_noint.mps’);
>> spy(prob_1.A);
>> Opt_1 = opti(’f’,-prob_1.f,’ineq’,prob_1.A,prob_1.ru,’bounds’,prob_1.lb,prob_1.ub,’xtype’,xtype)
-----------------------------------------------------Mixed Integer Linear Program (MILP) Optimization
min f’x
s.t. rl <= Ax <= ru
lb <= x <= ub
xi = Integer / Binary
-----------------------------------------------------Problem Properties:
# Decision Variables:
50
# Constraints:
191
# Linear Inequality:
41
# Bounds:
100
# Integer Variables:
50
-----------------------------------------------------Solver Parameters:
Solver:
CBC
----------------------------------------------------->> tic,[x,fval,exitflag,info] = solve(Opt_1),toc
x =
1.0000
1.0000
1.0000
1.0000
0
1.0000
1.0000
0
1.0000
0
0
1.0000
1.0000
1.0000
1.0000
0
1.0000
1.0000
1.0000
0
0
1.0000
1.0000
0
1.0000
1.0000
0
48/65
1.0000
1.0000
1.0000
1.0000
1.0000
0
1.0000
0
0
0
0
0
1.0000
0
0
1.0000
0
0
0
1.0000
0
1.0000
1.0000
fval =
-1046
exitflag =
1
info =
Nodes:
AbsGap:
RelGap:
Time:
Algorithm:
Status:
2
24.9233
0.0238
0.0581
’CBC: Branch and Cut using CLP’
’Integer Optimal’
Elapsed time is 0.060537 seconds.
>>
Con ilp_bb_1:
>> tic,[xint,j]=ilp_bb_1(-prob.f,A,prob_1.ru,prob_1.lb,prob_1.ub),toc
xint =
1
1
1
1
0
1
1
0
1
0
0
1
1
1
1
0
1
1
1
0
0
1
1
0
1
1
0
1
1
1
1
1
0
1
0
0
0
0
0
1
0
0
1
0
0
0
1
0
1
1
j =
-1.0460e+03
Elapsed time is 7295.551641 seconds.
>>
49/65
El tiempo empleado es inaceptable.
La estrategia es mala al no considerar ninguna forma de ahorrarse visitar todos los nudos del árbol de
búsqueda.
Selección de la variable de ramificación
50/65
Es evidente que para mejorar las prestaciones de un algoritmo Branch and
Bound es necesario hacer algo que evite visitar todos los nudos del árbol de
búsqueda de soluciones enteras. Veamos lo que hace BBMI.
Supongamos que se ha elegido i como el nudo que se va a analizar. Asociado a
ese nudo habrá un programa lineal de solución x .i/,
n
X
xBi D aN i0 C
aN ij . xj /
para i D 1; : : : ; mI
j DmC1
P
la función objetivo es x0 D aN 00 C jnDmC1 aN 0j . xj /: Los aN 0j ,
j D m C 1; : : : ; n, son los costes reducidos de las variables no básicas.
El paso siguiente de escoger la variable que va a definir la ramificación o división
que parta de i es crítico para las prestaciones del procedimiento. Apuntemos
alguna estrategia de mejora.
51/65
Penalizaciones de variables básicas
Supongamos que en el nudo i alguna variable básica xp no es entera debiendo
serlo, es decir
n
X
xp D aNp0 C
aNpj . xj / D np0 C fp0;
j DmC1
donde np0 designa la parte entera de xp , bxp c, y fp0 la parte fraccionaria.3
La división que definirá la variable xp es
xp np0 C 1 y xp np0:
3
Recordemos b3; 22c D 3; b 4; 12c D
5.
52/65
Del algoritmo dual del simplex podemos deducir que la imposición de la nueva
cota xp np0 C 1 a la variable xp implicará un empeoramiento (reducción) del
valor de la función objetivo.
Como la variable básica xp , realizando una iteración del método dual del
simplex, pasaría como mínimo a ser no básica, otra variable no básica, xj ,
que estuviese en uno de sus límites, se incrementaría o decrementaría según
estuviese en el inferior o superior.
Es decir, se reduciría el valor de la función objetivo en una cantidad o
penalización que vendría dada por la expresión
PU D
8
ˆ
< mKınj;aNpj <0 .1
fp0/
aN 0j
;
aNpj
si xj D lj ; o por
ˆ
: mKınj;aNpj >0 .1
fp0/
aN 0j
;
aNpj
si xj D uj :
Razonando de manera similar, la imposición de la nueva cota xp np0 a la
variable xp implicará un empeoramiento (reducción) del valor de la función
objetivo que vendrá dado por la expresión
PD D
8
< mKınj;aNpj >0 fp0 aN 0j ;
aNpj
aN 0j
: mKın
j;aNpj <0 fp0 aNpj ;
si xj D uj :
Cualquier solución entera que se pudiese obtener partiendo del nudo i estaría
por consiguiente acotada superiormente por
mKaxfx0
si xj D lj ;
PU ; x0
PD g:
Se podrá conseguir una mejor solución desde el nudo i si y sólo si
mKınfPU ; PD g < x0
zPE
para cada variable básica que no sea entera debiendo serlo.
Si esta condición no se cumple se abandona el nudo i .
53/65
54/65
Si todas las variables xp cumplen esa condición, la cuestión siguiente a solventar
es qué variable de ramificación se elige.
Se puede elegir aquella con la penalización asociada más pequeña.
O escoger aquella variable con una penalización asociada más grande y
comenzar imponiendo el límite opuesto al que determina esa penalización.
Si, por ejemplo, la penalización más grande asociada a una variable xp es PU ,
comenzar analizando el problema que resulte de imponer a esa variable la nueva
cota np0 e incluir en la lista de nudos a analizar más adelante el que resulte de
imponer la cota np0 C 1.
Si se almacenan las peores soluciones, una vez que se encuentre una buena, se
rechazarán rápidamente buena parte de los nudos que queden en la lista.
Penalizaciones de variables no básicas
55/65
Posibles cambios en variables no básicas que impondrían nuevas cotas en el
nudo i . Cualquier solución entera estaría acotada por
mKaxfx0
PU ; x0
PD g;
donde ahora
PU
D
(
8
aN 0j .1 fp0 /=. aNpj /;
ˆ
ˆ
ˆ
mK
ı
n
ˆ
j;
a
N
<0
pj
<
mKaxfaN 0j ; aN 0j .1 fp0 /=. aNpj /g;
(
ˆ
aN 0j .1 fp0 /=. aNpj /;
ˆ
ˆ
ˆ
: mKınj; aNpj >0
mKaxfaN 0j ; aN 0j .1 fp0 /=. aNpj /g;
j 2M
)
j 2J
)
j 2M
j 2J
si xj D lj
si xj D uj
y
PD
D
(
8
aN 0j fp0 =aNpj ;
ˆ
ˆ
ˆ
ˆ
< mKınj; aNpj >0 mKaxfaN ; aN f =aN g;
0j
0j p0
pj
(
ˆ
aN 0j fp0 =aNpj ;
ˆ
ˆ
ˆ
: mKınj; aNpj <0
mKaxfaN 0j ; aN 0j fp0 =aNpj g;
j 2M
)
j 2J
)
j 2M
j 2J
si xj D lj
si xj D uj :
J es el conjunto de índices de las variables no básicas que han de ser enteras y
M el de las no básicas que no han de ser enteras.
56/65
El nudo i se podrá desechar si, para cualquier xp a la que se introducen nuevas
cotas,
mKınfPU ; PD g x0 zPE :
Este criterio no es el único; recordemos el corte de Gomory asociado a xp ,
fp0 C
n
X
fpj xj 0;
j DmC1
que se habría de satisfacer en cualquier solución que se obtuviese partiendo del
problema actual.
Asociado a ese corte de Gomory se puede derivar la penalización
8
ˆ
ˆ
ˆ
ˆ
ˆ
ˆ
ˆ
ˆ
ˆ
mKınj
ˆ
ˆ
ˆ
ˆ
ˆ
ˆ
ˆ
<
PG D ˆ
ˆ
ˆ
ˆ
ˆ
ˆ
ˆ
ˆ
ˆ
mKınj
ˆ
ˆ
ˆ
ˆ
ˆ
ˆ
ˆ
:
8
ˆ
fp0 aN 0j =aNpj
ˆ
ˆ
ˆ
ˆ
< .1 fp0 /aN 0j =. aNpj /
ˆ
fp0 aN 0j =fpj
ˆ
ˆ
ˆ
ˆ
: .1 f /aN =.1 f /
p0 0j
pj
8
ˆ
fp0 aN 0j =aNpj
ˆ
ˆ
ˆ
ˆ
< .1 fp0 /aN 0j =. aNpj /
ˆ
fp0 . aN 0j /=fpj
ˆ
ˆ
ˆ
ˆ
: .1 f /. aN /=.1
p0
0j
9
j 2M>
>
>
>
>
j 2M=
si aNpj 0
si aNpj < 0
si fpj fp0 j 2 J >
>
>
>
>
si f > f j 2 J ;
pj
si xj D lj I
p0
si aNpj 0
j
si aNpj > 0
j
si fpj fp0 j
fpj / si fpj > fp0 j
9
2M>
>
>
>
>
2M=
2J >
>
>
>
>
2J ;
si xj D uj :
Se puede comprobar que
PG mKın.PU ; PD /;
lo que trae como consecuencia que el nudo i se podrá desechar si, para
cualquier xp a la que se introducen nuevas cotas,
PG x0
zPE :
57/65
Ejemplo Volvamos al ejemplo: maximizar 7x1 C 2x2
x1 C 2x2 4
5x1 C x2 20
2x1 2x2 7
x 2 Z2C :
s. a
58/65
Con las variables de holgura x3, x4 y x5, la relajación lineal da la solución
C
x0
x1
x2 C
3
x
11 3
1
x
11 3
5
x
11 3
8
x
11 3
C
C
C
C
16
x
11 4
2
x
11 4
1
x
11 4
6
x
11 4
D
D
D
C x5 D
332
11
36
11
40
11
75
;
11
donde x3 D x4 D 0.
.0/
De aquí que zPL D x0 D 332=11 y x .0/ D Œ36=11; 40=11; 0; 0; 75=11T .
Trabajemos con las penalizaciones que incluye BBMI. Las variables básicas son
x1, x2 y x5. Esta última no se requiere que sea entera. Las variables no básicas
de la solución obtenida están todas en su límite inferior.
59/65
Para empezar, como todavía no tenemos ninguna solución factible de zPE ,
determinemos las penalizaciones por bifurcar a:
x1 3 W
0
P1D
x1 4 W
0C
P1U
x2 3 W
0
P2D
D mKınj; aNpj >0 faN 0j fp0 =aNpj ; j 2 M g D
mKın
D mKınj; aNpj <0 faN 0j .1
fp0 /=. aNpj /; j 2 M g D mKın
D mKınj; aNpj >0 faN 0j fp0 =aNpj ; j 2 M g D
mKın
16
11
3
11
.1
16
11
7
11
3
11
2
11
3
/
11
D 24
I
11
1
; 3
11 11
1
11
7
11
24
D 11
I
5
11
D
21
:
55
0C
La penalización P2D
la podemos considerar 1 pues de la solución del programa
lineal anterior se observa que al incrementar la variable x2, x3 y x4 decrecen.
De las penalizaciones calculadas, la mejor (menos mala) se obtiene al hacer
5
1
40
x2 3. Añadamos esa condición, es decir, x2 D 11
x
x 3.
3
11
11 4
60/65
La relajación lineal resultante anterior más esta condición queda:
C
x0
x1
x2 C
3
x
11 3
1
x
11 3
5
x
11 3
8
x
11 3
5
x
11 3
C
C
C
C
16
x
11 4
2
x
11 4
1
x
11 4
6
x
11 4
1
x
11 4
D
D
D
C x5
D
C x6 D
332
11
36
11
40
11
75
11
7
:
11
Después de una iteración del método dual del simplex se llega a la solución
óptima siguiente:
x0
.1/
C 75 x4
C 15 x4
C 35 x6
1
x
x1
5 6
x2
C x6
2
x C x5 C 85 x6
5 4
11
x3 C 15 x4
x
5 6
D 149
5
17
D 5
D 3
D 29
5
D 75 :
De aquí que zPL D x0 D 149=5 y x .1/ D Œ17=5; 3; 7=5; 0; 29=5T , con x6 D 0.
La única variable que debe ser entera en esta solución y no lo es x1.
61/65
Calculemos ahora las penalizaciones por bifurcar a:
x1 3 W
1
P1D
x1 4 W
1C
P1U
D mKınj; aNpj >0 faN 0j fp0 =aNpj ; j 2 M g D
mKın
D mKınj; aNpj <0 faN 0j .1
fp0 /=. aNpj /; j 2 M g D mKın
2
5
7
5
3
5
.1
1
5
D 14
5
2
/ 15 D 95 :
5
De estas penalizaciones la mejor se obtiene al hacer x1 4. Añadamos la
1
condición x1 D 17
x C 15 x6 4.
5
5 4
La relajación lineal resultante anterior más esta condición queda:
C 75 x4
C 15 x4
x0
x1
x2
C
3
x
5 6
1
x
5 6
C
2
x C x5 C
5 4
x3 C 15 x4
1
x
C
5 4
x6
8
x
5 6
11
x
5 6
1
x
5 6
D 149
5
17
D 5
D 3
D 29
5
D 57
x7 D 35 :
Otra iteración del dual del simplex obtiene la primera solución factible entera
x .2/ D Œ4; 0; 8; 0; 1T , con zPE D 28.
62/65
Después se analizaría el nudo 3, rechazándose inmediatamente pues del cálculo
anterior de penalizaciones se vería que la función objetivo obtenible sería
1
149=5 P1D
D 135=5 D 27 < 28. Posteriormente también se rechazaría el
nudo 4.
0
x2 ≤ 3
x2 ≥ 4
1
4
x1 ≤ 3
x1 ≥ 4
3
2
Solución entera.
x(2) = [4, 0, 8, 0, 1]T
(2)
zP L = z P E = zP E = 28
Programación No lineal en variables Enteras
63/65
La dificultad de estos problemas aumenta considerablemente, pero nada impide
utilizar la estrategia Branch and Bound utilizando un solver no lineal –fmincon
por ejemplo– en vez del lineal que hemos presentado. El código de antes
adaptado para resolver un Programa Entero no lineal es éste:
function [xstar,Jmin] = inlp_bb(fun,x0,A,b,xlb,xub,NLCon,xstar,Jmin)
% [xstar,Jmin] = ilp_bb(fun,x0,A,b,xlb,xub,NLCon,xstar,Jmin)
% Integer Non-linear programming by (recursive) branch & bound
% Uses fmincon.m from the optimisation toolbox to solve the relaxed problem
optns = optimset(’display’,’off’); tol = 1e-3; % anything less is zero
if ~exist(’xlb’,’var’),
xlb = 0*x0; end
if ~exist(’xub’,’var’),
xub = Inf*ones(size(x0)); end
if ~exist(’xstar’,’var’), xstar = []; end
if ~exist(’Jmin’,’var’), Jmin = Inf; end
% Now solve relaxed NL constrained problem
[x,j,exitflag] = fmincon(fun,x0,A,b,[],[],xlb,xub,NLCon,optns);
if exitflag <1; return, end
% disp(exitflag) if infeasible or not OK
if j>Jmin + tol; return, end
% for robustness; truncate
idx = find(abs(x-round(x)) > tol);
% which are non-integers?
if isempty(idx)
% All integer solutions
if j<Jmin, xstar=round(x); Jmin=j; end
% New integer optimum
return
% problem solved
end
% Depth first recursion
xnlb = xlb; xnlb(idx(1)) = ceil(x(idx(1)));
% new lower bound
[xstar,Jmin] = inlp_bb(fun,x,A,b,xnlb,xub,NLCon,xstar,Jmin);
xnub = xub; xnub(idx(1)) = floor(x(idx(1)));
[xstar,Jmin] = inlp_bb(fun,x,A,b,xlb,xnub,NLCon,xstar,Jmin);
end
Resolviendo
64/65
minimizar .x1
sujeta a
5/2 C .x2
0;4 .x1
7/2
x1 x2 .x2 C 1/
5;7/3 4;2 C x2 0
0;53x1 C x2 10
1;1x1 C x2 5:
Partiremos del punto x 0 D Œ2;5; 1;5T .
MINLP Optimization Contours + Constraints - Min: [3; 8] Fval: -211
-8
00
0
11
-80
10
0
00
200
9
0
-60
-4
-2
0
8
00
3
8
00
>> A=[0.53 1; -1.1 1]; b=[10,5]’;
>> lb=[0;0]; ub=[10;10];
>> x0=[2.5,1.5]’;
>> fun1=@(x) (x(1)-5).^2+(x(2)-7).^2-x(1)*x(2)*(x(2)+1);
>> [xstar,Jmin]=inlp_bb(fun1,x0,A,b,lb,ub,@circlecon)
xstar =
-200
12
-1
20
0
-1
00
0
-6
-40
0
13
x2
-40
7
0
Jmin =
6
-20
0
-211
5
0
>> plot(Opt)
4
3
-2
-1
0
1
2
3
x1
La rutina con la condición no lineal es
function [c,ceq] = circlecon(x)
c = 0.4*(x(1)-5.7)^3+x(2)-4.2;
ceq = [];
end
4
5
6
7
8
65/65
Con la herramienta OPTI Toolbox
MINLP Optimization Contours + Constraints - Min: [3; 8] Fval: -211
-8
0
20
0
00
0
-80
0
-60
-4
200
00
0
-2
0
00
8
x =
-40
0
7
3
8
-1
-1
00
0
00
-6
-40
fun=@(x) (x(1)-5).^2+(x(2)-7).^2-x(1)*x(2)*(x(2)+1);
A=[0.53 1; -1.1 1]; b=[10,5]’;
13
lb=[0;0]; ub=[10;10];
nlcon=@(x) 0.4*(x(1)-5.7)^3+x(2);
12
nlrhs=4.2;
11
x0=[2.5,1.5]’;
xtype = ’II’;
10
Opt=opti(’fun’,fun,’ineq’,A,b,’bounds’,lb,ub,’nlcon’,nlcon,’nlrhs’,nlrhs,’xtype’,xtype)
9
[x,fval] = solve(Opt,x0)
-200
>>
>>
>>
>>
>>
>>
>>
>>
>>
x2
6
-20
0
5
0
fval =
-211
>> plot(Opt)
4
3
-2
-1
0
1
2
3
x1
4
5
6
7
8