Relaxação Lagrangeana – Teoria

Fevereiro 5, 2012

Neste post vamos comentar brevemente sobre relaxações de formulações lineares inteiras e focar especialmente em relaxações lagrangeanas. Depois, vamos demonstrar a aplicação dessa técnica no problema da localização de instalações (facility location). A menos que dito o contrário, estamos falando de problemas de minimização.

Relaxações

Em programação linear inteira, usamos o termo relaxação para indicar que algumas das restrições do problema original foram relaxadas. Provavelmente o exemplo mais conhecido de relaxação é o da relaxação linear. Nesse caso em particular, removemos as restrições de integralidade do modelo.

Mais formalmente, se tivermos uma formulação F representada por

z = \min \{c(x): x \in P\}

sendo P o conjunto de pontos que satisfazem as restrições dessa formulação, c(x) a função objetivo e z o valor da solução ótima. Dizemos que uma formulação F_R, dada por

z_R = \min \{f(x): x \in R\}

é uma relaxação se todo ponto de P está em R, ou seja P \subseteq R e se f(x) \le c(x) para x \in P.

Note que devido a essas propriedades, a solução ótima da relaxação é um limitante dual (ou seja, no nosso caso é um limite inferior e superior no de maximização) para a formulação original ou seja, z_R \le z. Limitantes duais podem ser utilizados, por exemplo, para melhorar o desempenho do algoritmo de Branch and Bound.

A vantagem de se usar relaxações é que em certos casos elas são mais fáceis de serem resolvidas. Por exemplo, as relaxações lineares podem ser resolvidas em tempo polinomial, enquanto a versão restrita a inteiros, em geral, é NP-difícil.

Podemos também relaxar uma formulação removendo algumas de suas restrições. Note, entretanto, que quanto mais relaxamos a formulação pior tende a ser a qualidade do limitante dual (considere o caso extremo onde removemos todas as desigualdades do modelo: não deixa de ser uma relaxação, mas o limitante dual obtido não serve para nada).

O ideal é que encontremos uma relaxação fácil de resolver e que dê limitantes duais de qualidade, ou seja, bem próximos do valor ótimo da formulação original.

Relaxação Lagrangeana

Existe uma técnica denominada relaxação lagrangeana [1], que consiste em remover algumas das restrições da formulação original, mas tenta embutir essas desigualdades na função objetivo. A ideia é penalizar a função objetivo quando as restrições removidas forem violadas. O “peso” dessas penalidades é controlado por coeficientes chamados multiplicadores lagrangeanos.

Em termos de matemáticos, suponha que tenhamos uma formulação com m restrições e n variáveis:

z^* = \min cx

Sujeito a

\begin{array}{llcl}    & A_1x & \le & b_1 \\   & A_2x & \le & b_2 \\   & x & \ge & 0 \\   & x & \in & Z^n  \end{array}

Onde A_1x \le b_1 representa as desigualdades que queremos relaxar e A_2x \le b_2 as desigualdades que vamos manter. A relaxação lagrangeana é então dada por

z(u) = \min cx - u(b_1 - A_1x)

Sujeito a

\begin{array}{llcll}   & A_2x & \le & b_2 & \\  & x & \ge & 0 & \\  & x & \in & Z^n & \\  & u & \ge & 0, & u \in R^m_1  \end{array}

Onde u é um vetor representando os multiplicadores lagrangeanos (um por cada desigualdade removida). Por exemplo, considere a desigualdade j que será removida:

\sum_{i = 1}^n a_{ij} x_i \le b_j

Ao remover (11), o seguinte fator é subtraído da função objetivo

(1) u_j (b_j - \sum_{i = 1}^n a_{ij} x_i)

Note que quanto mais “violada” a desigualdade estiver, mais negativo o termo (1) vai ficar. Como a função é de minimização e estamos subtraindo, intuitivamente a solução ótima para o problema tenderá a não violar as desigualdades, melhorando o resultado do limitante dual.

Escolhe-se remover as desigualdades de forma que o problema resultante seja mais fácil resolver para um valor fixo de multiplicadores lagrangeanos. Porém, agora temos que resolver um outro problema: encontrar um conjunto de multiplicadores lagrangeanos u^* que minimize z(u).

Aplicação ao problema da localização de instalações

O problema da localização de instalações é bastante estudado na literatura. Uma de minhas iniciações científicas consistia em estudar algoritmos aproximados para esse problema.

Basicamente temos um grafo bipartido onde uma partição corresponde ao conjunto de clientes C e a outra ao conjunto de fábricas/instalações F. Cada cliente i pode ser potencialmente atendido pela fábrica j, a custo um c_{i,j}. Além disso, paga-se um custo c_j por abrir a fábrica j, dado por y_j. O objetivo é atender todos os clientes com o menor custo.

Vamos apresentar uma formulação linear inteira para esse problema. Seja x_{i,j} uma variável binária que vale 1 sse o cliente i é atendido pela fábrica j e a variável y_j uma variável binária que vale 1 sse a fábrica j foi aberta. Temos a seguinte formulação:

z = \min \sum_{i \in C, j \in F} x_{ij} c_{ij} + \sum_{j \in F} y_j c_j

Sujeito a

\begin{array}{llcll}   (2) & \sum_{j \in F} x_{ij} & \ge & 1 & \forall \, i \in C \\  (3) & x_{ij} & \le & y_j & \forall i \in C, \forall j \in F \\  & x_{i,j}, y_{j} & \in & \{0, 1\} & \forall i \in C, \forall j \in F  \end{array}

A desigualdade (2) diz que cada cliente deve ser atendido por pelo menos uma fábrica e a desigualdade (3), que uma fábrica deve estar aberta para que possa atender a algum cliente.

Essa é uma formulação bem simples e temos basicamente dois grupos de desigualdades. Vamos tentar relaxar o grupo de desigualdades indicado por (2). Teremos:

z_1(u) = \min \sum_{i \in C, j \in F} x_{ij} c_{ij} +
\sum_{j \in F} y_j c_j - \sum_{i \in C} u_i (\sum_{j \in F} x_{ij} - 1)

Sujeito a

\begin{array}{llcll}   & x_{ij} & \le & y_j & \forall i \in C, \forall j \in F \\  & u_i & \ge & 0 & \forall i \in C\\   & x_{i,j}, y_{j} & \in & \{0, 1\} & \forall i \in C, \forall j \in F  \end{array}

Com um pouco de álgebra, é possível “fatorar” o somatório da função objetivo:

z_1(u) = \min \sum_{j \in F} (\sum_{i \in C} (x_{ij} (c_{ij} - u_{i})) + c_j y_j)) + \sum_{i \in C} u_i

Como o termo \sum_{i \in C} u_i é constante, podemos resolver o problema sem ele e depois somá-lo ao valor da solução obtida. Além disso, note que cada fábrica contribui de maneira independente à função objetivo e em cada restrição só aparecem variáveis associadas a uma dada fábrica.

Desta forma, podemos decompor essa formulação por cada fábrica j:

\min \sum_{i \in C} (x_{ij} (c_{ij} - u_{i})) + c_j y_j

Sujeito a:

\begin{array}{llcll}   & x_{ij} & \le & y_j & \forall i \in C\\  & x_{i,j}, y_{j} & \in & \{0, 1\} & \forall i \in C  \end{array}

Cada um desses problemas pode ser resolvido por inspeção. Seja c'_{ij} = c_{ij} - u_{i}. Vamos considerar duas possibilidades. Se fizermos y_j = 0, não atenderemos nenhum cliente, levando a uma solução de custo 0. Se y_j = 1, só compensa atender aos clientes com c'_{ij} < 0. Assim, o custo ótimo dessa formulação é:

\min\{0, c_j + \sum_{c'_{ij} < 0} c'_{ij} \}

A outra alternativa consiste em relaxar a desigualdade (3). Nesse caso, depois de alguma álgebra, teremos:

z_2(u) = \min \sum_{i \in C, j \in F} x_{ij} (c_{ij} - u_{ij}) + \sum_{j \in F} y_j (c_j + (\sum_{i \in C} u_{ij}))

Sujeito a:

\begin{array}{llcll}   & \sum_{j \in F} x_{ij} & \ge & 1 & \forall \, i \in C \\  & u_{ij} & \ge & 0 & \forall i \in C, j \in F\\  & x_{i,j}, y_{j} & \in & \{0, 1\} & \forall i \in C, \forall j \in F  \end{array}

Note que como não há restrições misturando x e y, podemos resolver o modelo para cada tipo de variável separadamente. No caso da variável y, como c_j \ge 0 e u \ge 0, não compensa abrir fábrica nenhuma. No caso de x, podemos resolver de maneira independente para cada x_i, por inspeção. Seja c'_{ij} = c_{ij} - u_{ij}. Se c'_{ij} > 0, podemos fazer x_{ij} = 1. Caso contrário x_{ij} = 0, exceto no caso em que c'_{ij} \ge 0 para todo j \in F. Nesse caso, devido à restrição (25), temos que fazer x_{ij^*} = 1, onde j^* = \mbox{argmin}_{j \in F} \{c'_{ij} \}.

Temos duas formulações que podem ser facilmente resolvidas. Entretanto, devemos levar em conta qual possui a formulação mais forte. Em outras palavras, dado z_1^* = \max \{ z_1(u) : u \ge 0 \} e z_2^* = \max \{ z_2(u) : u \ge 0 \}, queremos o limitante dual que mais se aproxima de z. Entretanto, não sei dizer qual das duas relaxações apresentadas é a mais forte.

Conclusão

A relaxação Lagrangeana é uma técnica poderosa para se obter limitantes duais de problemas combinatórios que podem ser modelados como programas lineares inteiros. Esses limitantes podem ser usados na tentativa de melhorar o desempenho de algoritmos branch-and-bound.

Além do mais, para certos problemas, existe a possibilidade de ajustar a solução obtida via relaxação de modo a obter uma solução viável para o problema original, esperando que a qualidade da solução não seja muito pior.

Esse post apresentou a teoria por trás da relaxação lagrangeana. Num próximo post vou apresentar uma heurística que visa encontrar bons multiplicadores lagrangeanos.

Referências

[1] The Lagrangian relaxation method for solving integer programming problems – M.L. Fisher Management Science, 1981.

Anúncios

Visualização ótima de desenhos fisicamente realizáveis

Outubro 23, 2011

Hoje vou descrever o trabalho que apresentei no Sibgrapi 2011, entitulado “Determining an Optimal Visualization of Physically Realizable Symbol Maps“, escrito em conjunto com meus orientadores (professores Pedro e Cid) e o professor Tallys da Universidade de Miami. Há um

Trata-se de uma extensão do estudo que realizamos para o problema de mapas de símbolos proporcionais usando desenhos em pilha, sobre o qual comentei em um post anterior.

Desenho fisicamente realizável

Um desenho fisicamente realizável é uma generalização de um desenho em pilha, pois este permite ordens cíclicas, como no exemplo abaixo.


Desenho fisicamente realizável e desenho em pilha

Podemos pensar que desenhos fisicamente realizáveis são aqueles feitos a partir de discos de papel, desde que estes não sejam dobrados nem cortados.

O problema é encontrar um desenho fisicamente realizável que maximize a soma dos comprimentos das bordas visíveis. Esse problema é chamado Max-Total para desenhos fisicamente realizáveis.

Tínhamos desenvolvido um modelo de programação linear inteira para o problema Max-Total para desenhos em pilha e mostramos que um subconjunto das desigualdades desse modelo, resolvia o problema Max-Total para desenhos fisicamente realizáveis.

Discutimos aspectos teóricos sobre as desigualdades utilizadas no modelo, argumentando que elas tornam a formulação apertada. Na prática um modelo com formulação apertada tende a ser resolvido mais rapidamente pelo algoritmo de branch and bound.

Técnica de decomposição

Desenvolvemos uma técnica de decomposição nova, que só serve para esse tipo de desenho. A ideia básica das decomposições é quebrar o conjunto de discos de entrada em componentes menores e mostrar que podemos resolver a cada componente de maneira independente e então combiná-las em tempo polinomial para gerar a solução ótima da instância completa. Uma decomposição óbvia é considerar conjuntos de discos disjuntos no plano.

Para desenhos fisicamente realizáveis, podemos ignorar a interseção de discos em alguns casos. Isso por sua vez pode gerar novas componentes disjuntas. No exemplo abaixo, mostramos que é possível ignorar a interseção dos discos nas faces amarelas. Isso quebra a componente em duas que podem ser resolvidas isoladamente.

Exemplo de decomposição

Junto com as outras decomposições que já havíamos desenvolvido para desenhos em pilha, conseguimos diminuir o tamanho das instâncias de teste.

Resultados computacionais

Após resolvermos a formulação do problema para desenhos fisicamente realizáveis, constatamos que na grande maioria dos casos, o valor da solução era exatamente igual ao da solução para desenhos em pilha.

Mesmo para as instâncias onde houve diferença, esse valor foi muito pequeno e visualmente é muito difícil distinguir as soluções, conforme pode ser visto na figura a seguir.

Zoom de uma solução ótima usando desenho em pilha e desenho fisicamente realizável

Entretanto, também tivemos uma boa surpresa, que foi o tempo de execução desse novo modelo. Em média, ele foi aproximadamente 2 ordens de magnitude mais rápido do que modelo para desenhos fisicamente realizáveis.

Conclusão

Devido aos tempos obtidos com o novo modelo e o fato de as soluções obtidas serem muito próximas a desenhos em pilha, uma ideia é usarmos esse modelo para resolver o problema de desenhos em pilha, adicionando restrições para remover os eventuais ciclos que apareçam.

Podemos adicionar essas restrições de remoção de ciclo através de um algoritmo de planos de corte, de um modo similar ao usado para resolver o problema do caixeiro viajante.


Mapas de Símbolos Proporcionais

Abril 17, 2011

Sexta-feira passada apresentei uma palestra na série de seminários do LOCo. Falei sobre o problema que estou atacando no mestrado.

Um mapa de símbolos proporcionais emprega símbolos para exibir eventos associados a alguma localidade e intensidade. O símbolo é posicionado no local de ocorrência do evento e a sua área fica sendo proporcional à intensidade desse evento. Focamos em mapas que usam círculos opacos como símbolos. A figura abaixo é um exemplo representando as populações das maiores cidades do Brasil.

Note que há sobreposição entre os discos. Dependendo da ordem em que os empilhamos, pode ser que haja mais ou menos porções dos discos visíveis. Há casos ruins em que a borda do disco fica totalmente encoberta, como na figura abaixo.


Não podemos inferir o centro nem o raio do disco com bordo escondido.

Para tentar fazer com que o desenho tenha bastante borda visível, usamos uma métrica que consiste em maximizar o comprimento visível total dos discos. Com isso em mente, é possível definir um modelo de programação linear inteira, com um modelo que atribui cada disco a um nível.

Além do modelo básico, desenvolvemos algumas desigualdades adicionais para tornar o modelo mais forte, nos baseando em propriedades geométricas, que impedem certas configurações de acontecerem.

Também desenvolvemos duas técnicas de decomposição de instâncias. Um jeito trivial de decompor instâncias é considerar componentes de discos conexas de maneira independente.

Nossa primeira técnica de decomposição vem da observação que um disco que está contido no outro sempre pode ser desenhado por cima em uma solução ótima. Dessa forma, em uma instância como a da figura abaixo, a componente {a, b} pode ser desenhada de maneira independente de {c, d, e, f} e na hora de construir a solução ótima, basta desenhar a solução obtida para {a, b} em cima da solução {c, d, e, f}.

Componente {a,b} está contida em {f} e pode ser resolvida separadamente.

Podemos definir um grafo Gs sobre os discos de S, com conjunto de vértices correspondendo aos discos e há uma aresta (i, j) em Gs se os discos correspondentes se interceptam. Falei sobre esse tipo de grafo em um post anterior.

A outra técnica consiste em remover um ponto de articulação de Gs e replicá-lo nas componentes conexas resultantes, como na figura abaixo.


O disco f representa um ponto de articulação em Gs. As figuras (ii), (iii) e (iv) são as componentes resultantes de sua remoção, com f replicado.

Mostramos que é possível resolver cada componente com o ponto de articulação replicado de maneira independente. Depois, basta juntar os discos mantendo a ordem relativa de cada sub-solução. Isso pode ser feito através de uma ordenação topológica.

Implementamos todas essas ideias e reportamos os resultados. As instâncias testadas foram as mesmas de um trabalho anterior no qual nos baseamos. As técnicas de decomposição foram efetivas na redução do tamanho das instâncias e o resolvedor XPRESS conseguiu resolver nosso modelo para quase todas as instâncias. Porém, ficaram algumas componentes sem serem resolvidas, o que nos motiva a procurar novos modelos e novas desigualdades.

Esse nosso trabalho foi aceito na Computational Geometry and Applications 2011.


Lifted Cover Inequalities

Fevereiro 4, 2011

Nesse post falarei sobre as desigualdades de cobertura erguidas (em uma horrorosa tradução livre de lifted cover inequalities). A minha ideia inicial era comentar sobre cortes automáticos que resolvedores como o XPRESS inserem no programa para apertar o modelo. Há dois tipos de desigualdades usadas como cortes: as de cobertura erguidas e aquelas resultantes do procedimento de Chvátal-Gomory.

Como esse segundo método faz uso direto do tableau do simplex, antes eu gostaria de escrever um post revisando o primal simplex e o dual simplex. Por isso, ater-me-ei às coberturas erguidas.

Essas desigualdades surgiram a partir do problema da mochila binária. Vamos relembrar o problema: dado um conjunto de n itens — cada qual com um peso e valor — e uma mochila de capacidade b, queremos escolher um subconjunto de itens cuja soma dos pesos respeite à capacidade da mochila e cuja soma dos valores seja máxima.

Vamos supor que os valores dos itens sejam dados por e os pesos por (para facilitar, suponha que os itens estão ordenados em ordem não crescente de peso, ou seja ). Para cada item i, definimos a variável binária que vale 1 se vamos colocar o item na mochila e 0 caso contrário. O modelo de programação linear inteira é dado por:

Sujeito a:

(1)

(2)

Definições

Definimos qualquer subconjunto de itens que cabem na mochila como conjunto independente. Assim, se S é um conjunto independente, temos que

Se, por outro lado, o subconjunto não cabe na mochila, denominamo-no uma cobertura.

Adicionalmente, uma cobertura é dita minimal se a remoção de qualquer item do subconjunto torna-o independente. Dada uma cobertura C, uma desigualdade de cobertura é definida como:

Ela é claramente válida pois por definição incluir todos os itens de uma cobertura violará a capacidade da mochila.

Cobertura Estendida

Podemos estender essa desigualdade da seguinte forma: Seja j o item mais pesado dessa cobertura. Se colocarmos um item mais pesado que j, ainda assim não podemos ter mais do que |C|-1 itens. Aplicando esse argumento várias vezes, podemos inserir todos os itens j’ < j (que são mais pesados que j por causa da ordenação inicial) na cobertura, para deixar a desigualdade mais apertada. Definimos então a cobertura estendida E(C) como sendo a cobertura C mais os itens mais pesados que o item j. A desigualdade abaixo é então válida:

Lifting das Desigualdades de Cobertura

Uma outra maneira de melhorar uma desigualdade de cobertura é através de lifting, que consiste em aumentar os coeficientes das variáveis do lado esquerdo de uma desigualdade na forma (<=), desde que ela continue viável. Obviamente só conseguimos fazer lifting em desigualdades que não representam facetas, uma vez que a desigualdade que sofreu lifting domina a original.

As desigualdades de cobertura estendida é uma forma de lifting, pois estamos aumentando de 0 para 1 o coeficiente dos itens mais pesados que j. Porém, usando um método diferente (porém mais caro), podemos conseguir aumentar os coeficientes para valores maiores do que 1. Vamos considerar um exemplo de mochila, extraído de [1]:

Uma desigualdade de cobertura válida é

(3)

A cobertura estendida consiste em incluir os itens e , os itens mais pesados que , resultado na desigualdade de cobertura estendida: .

Ao invés de aumentar o valor de para 1 em (3), vamos usar uma variável , levando a

(4)

A pergunta é: qual o maior valor que pode assumir, desde que (4) se mantenha válida? Para isso, temos que encontrar o valor de considerando todas as possibilidades e escolher o menor. Primeiro, se então é ilimitado, o que não ajuda muito. Vamos considerar então que . Com isso, podemos reduzir o problema:

Agora queremos encontrar o maior valor de em (4) que respeite a desigualdade. Para isso, devemos considerar o pior caso, ou seja, em que a soma dos outros termos é máxima,

(5)

A soma máxima pode ser determinada resolvendo-se outro problema da mochila:

Sujeito a

A solução ótima desse problema é 1, o que significa que por (5), . Temos então a desigualdade

(6)

Repetindo esse procedimento para fazer um lifting em (6) aumentando o coeficiente de para , concluiremos que , levando a desigualdade

(7)

Que domina inclusive a desigualdade de cobertura estendida. Porém, isso tem seu preço: para estender uma desigualdade de cobertura só precisamos inserir os itens mais pesados, enquanto nesse último procedimento tivemos que resolver dois novos problemas da mochila. Ou seja, para resolver um problema NP-difícil precisamos resolver sub-problemas que também são NP-difíceis? Veremos…

Programação Dinâmica

Como se sabe, o problema da mochila pode ser resolvido através de um algoritmo clássico de programação dinâmica que possui complexidade pseudo-polinomial de O(nb). Basicamente define-se uma submochila m(i, j), que quer representa o maior valor que se consegue com uma mochila usando os i primeiros itens e com capacidade j. A recorrência é

Porém, é possível indexar o valor da mochila ao invés da capacidade. Nesse sentido definimos m'(i, v) como o menor peso que conseguimos usando os i primeiros itens e com valor v. A recorrência é semelhante:

Com a restrição de que os valores dos itens sejam inteiros. Para encontrar maior valor da mochila, basta percorrer a matriz m’ e determinar o maior v, tal que m'(n, v) <= b.

A complexidade desse segundo algoritmo é O(max_v n), onde max_v é o maior valor que a mochila pode ter (ou seja, um limite superior para o custo da função objetivo). No nosso caso, nossa função objetivo corresponde a um pedaço do lado esquerdo da desigualdade de cobertura. Como tal desigualdade é válida, sabemos que o lado esquerdo nunca ultrapassará |C|-1, logo max_v . Concluímos que os subproblemas da mochila podem ser resolvidos em .

Desigualdades de Cobertura para outros problemas

As desigualdades levantadas para coberturas podem ser usadas para outros problemas. Se temos um PLI na forma para i=1,…,m onde representa a i-ésima linha de A, então podemos particionar P em vários problemas da mochila na forma e encontrar desigualdades fortes para esses subproblemas.

Como temos que

A esperança é que as desigualdades fortes para os subproblemas (em especial as desigualdades de cobertura) sejam úteis para P na prática. Mais detalhes sobre essa técnica, ver em [1].

Bibliografia

[1] Crowder, H. and Johnson, E.L. and Padberg, M. – Solving large-scale zero-one linear programming problems


Desigualdades Definidoras de Facetas

Outubro 22, 2010

No último post da série sobre combinatória poliédrica, apresentei o algoritmo de Branch and Cut. Um aspecto essencial no desempenho desse algoritmo é o quanto uma desigualdade é capaz de aproximar uma formulação da envoltória convexa. Vimos que as desigualdades definidoras de facetas são as mais eficientes na prática.

Hoje vou discutir um pouco sobre métodos para determinar se uma dada desigualdade define faceta para a envoltória convexa. Para tal, precisamos revisar alguns conceitos de álgebra linear.

Álgebra Linear

Independência linear: o conjunto de vetores  é linearmente independente se

(1)

então, para todo i = 1, …, n

Para . Em outras palavras, não existe um vetor que possa ser escrito como combinação linear de outro.

Independência (linear) afim: o conjunto de vetores  é afim independente se

(2)

e (3)

então, para todo i = 1, …, n

Para . Aqui temos a restrição adicional de que a soma dos coeficientes deve ser 0. Observe que se o vetor nulo estiver no conjunto de vetores, esse nunca será linearmente independente, pois se a igualdade em (1) valesse, qualquer coeficiente não-nulo multiplicando o vetor nulo manteria a igualdade. Porém, a adição da desigualdade (3), garante nesse caso que o coeficiente do vetor nulo seja 0. Essa notação é interessante para não termos que tratar o vetor nulo como um caso especial ao fazer algumas contas.

Vamos apresentar a seguir dois possíveis métodos que podem ser utilizados para decidir se uma dada desigualdade define faceta de um poliedro.

Método Direto

Existe um teorema que diz que a dimensão de um poliedro é igual ao maior número de pontos afim independentes pertencentes a ele, menos um. Assim, se encontramos k soluções afim independentes do nosso problema, sabemos que ele tem dimensão maior ou igual a k-1. Se o modelo tem n variáveis e o poliedro tem dimensão n, dizemos que ele possui dimensão cheia. É muito mais difícil trabalhar com poliedros de dimensão menor que n, portanto, para fins de explicação, vamos assumir daqui pra frente que o poliedro tem dimensão cheia.

Uma face qualquer de um poliedro, é também um poliedro, só que com dimensão menor. Logo, para determinar a dimensão de uma face basta enumerar pontos afim independentes que pertençam a essa face.

Como vimos, a dimensão da faceta é uma a menos do que a dimensão do poliedro. Para determinar se uma desigualdade define faceta, precisamos calcular a dimensão da insterseção da mesma com o poliedro. Para isso, devemos encontrar pontos do poliedro, que satisfaçam essa desigualdade na igualdade e que ainda sejam afim independentes entre si. Se a dimensão dessa interseção for uma a menos do que a do poliedro, tal desigualdade define uma faceta do poliedro.

Método Indireto

Encontrar pontos afim independentes pode não ser uma tarefa simples. Uma alternativa consiste em usar um método indireto.

Dada uma desigualdade válida da forma , queremos decidir se a mesma define uma faceta do poliedro P. Considere a interseção dessa desigualdade com o poliedro P, representada por . Considere os pontos x em P que satisfazem e portanto satisfazem .

Cada uma dessas igualdades formará uma linha de um sistema linear onde as variáveis são os coeficientes . Se esse sistema tiver solução única, pode-se mostrar que a desigualdade define uma faceta e que é equivalente a , ou seja .

O método indireto procura resolver esse sistema linear por partes, escolhendo pontos de P que satisfazem (e por consequência ) com o objetivo de mostrar que os coeficientes são iguais a . Se formos capazes de definir todos os coeficientes a partir desses pontos, a desigualdade define uma faceta de P. A conveniência é que não precisamos determinar se os pontos são afim independentes, pois isso é considerado implicitamente ao resolvermos o sistema linear. Além disso, como sabemos o valor esperado para os coeficientes , fica mais fácil procurar pontos que levem a essa conclusão.

Os pontos p e q definem a reta . O ponto s sozinho não consegue definir unicamente a reta .

A intuição geométrica do método pode ser vista na figura acima. Os pontos da interseção de com P são ‘p’ e ‘r’. Como a única reta (solução única) que passa por esses dois pontos é , a equação da interseção é equivalente a . Essa interseção é uma faceta de P. Já para a reta , o único ponto de sua interseção com P é o ponto s. Observe que existem infinitas retas que podem passar por ‘s’, como por exemplo a própria , mas também . Ou seja, o sistema linear correspondente não tem solução única e como podemos ver, a interseção não é uma faceta de P.

Exemplo

Para dar um exemplo do método indireto, vou usar o problema do conjunto independente. Dado um grafo não direcionado, queremos encontrar o maior conjunto de vértices de forma que nenhum deles seja vizinho entre si. Para isso, definimos uma variável binária  que vale 1 se o vértice v está no conjunto independente e 0 caso contrário. O modelo é bem simples,

(4) 

(5) 


Grafo com conjunto independente máximo marcado em laranja.

Esse modelo é bem fraco, no sentido de que na prática o valor da relaxação linear está bem distante da solução inteira ótima (poucas podas no Branch and Bound). Uma desigualdade bastante eficaz é a da clique maximal. Seja M uma clique maximal, ou seja, um conjunto de vértices tal que todos são adjacentes entre si e tal que não exista outro vértice fora de M, vizinho a todos os vértices de M. Considere a seguinte desigualdade

(6)

Ela é válida para nosso problema pois nenhum conjunto independente pode conter mais de um vértice de qualquer clique. Vamos verificar se ela define uma faceta. Suponha então que a interseção de (6) com a envoltória convexa do problema do conjunto independente seja da forma

(7)

A ideia é encontrar pontos válidos para o problema e que satisfaçam (6) na igualdade, pois assim sabemos que esses pontos satisfazem (7). Vamos chamar de o coeficiente de associado ao vértice v.

Dado um v na clique M, podemos fazer e o resto 0, que ele será um conjunto independente válido e satisfaz (6) na igualdade. Substituindo esse x em (7), temos que

(8)

Agora, dado um vértice u fora de M, temos a propriedade de que existe algum vértice em M não adjacente a u (do contrário poderíamos incluir u em M e esta não seria maximal, uma contradição). Seja então v* tal vértice em M. O ponto x, com e o resto 0, é um conjunto independente válido uma vez que u e v* não são adjacentes. Como este ponto satisfaz (7) na igualdade, temos que . Mas, por (8), chegamos que

(9)

Substituindo por pela igualdade (8), temos que cada coeficiente da desigualdade (7) é igual ao coeficiente da desigualdade (6) multiplicado por , concluindo que (7) define uma faceta.

Dá para substituir as desigualdades (5) do modelo por desigualdade (6), desde que todas as arestas (u, v) estejam na desigualdade de pelo menos uma clique maximal. Felizmente, existe um algoritmo polinomial para encontrar uma cobertura de arestas por cliques maximais [2].

Exemplo de uma cobertura por cliques maximais (arestas podem ser cobertas por mais de uma clique).

Referências

[1] Disciplina de Programação linear inteira com o prof. Cid.
[2] Nemhauser, G.L. and Sigismondi, G.A – Strong cutting plane/branch-and-bound algorithm for node packing


Algoritmo de Branch and Cut

Outubro 15, 2010

Esse é o primeiro post de uma série sobre combinatória poliédrica. Ela consiste no estudo da estrutura de problemas combinatórios e representa uma ferramenta importante na construção de modelos eficientes de progração linear inteira. Aqui vou falar sobre o algoritmo de Branch and Cut. Estou assumindo conhecimento prévio sobre Branch and Bound.

Abordagem Geométrica

Podemos pensar em programas lineares de forma geométrica. Para exemplificar, considere um programa linear de minimização (daqui pra frente, estarei sempre me referindo a esse tipo de PL), com duas variáveis (). Uma desigualdade da forma

representa um semi-plano. Em um programa linear, temos um conjunto de desigualdades que devem ser satisfeitas ao mesmo tempo. Geometricamente, isso equivale à interseção dos semi-planos. Um semi-plano é um polígono convexo (ilimitado) e a interseção de polígonos convexos é sempre um polígono convexo. Assim, o espaço de soluções é representado por um polígono convexo, ou para dimensões maiores, poliedros convexos.

Além disso, a função objetivo é da forma

Variando o valor de z, vemos que essa equação representa uma família de retas paralelas, e nosso objetivo é encontrar a reta para a qual z que seja mínimo (seja z* o z mínimo). Para dimensões maiores, a função objetivo é um hiperplano. A solução ótima (seja ela ) é um ponto do poliedro que satisfaz . É possível mostrar que a solução ótima de um programa linear está em um vértice do poliedro.


Poliedro representando a formulação de um problema com duas variáveis. As retas paralelas representam a função objetivo.

Para um programa linear inteiro, temos a restrição de que as variáveis são inteiras. Considere o conjunto P de pontos de coordenadas inteiras que são viáveis para nosso problema. Definimos a envoltória convexa como o menor polígono que contém P.

Na figura abaixo, os pontos inteiros que representam soluções viáveis para nosso problema, estão marcados em vermelho. A envoltória convexa é o polígono de borda preta. Já os polígonos em verde e vermelho são formulações válidas para nosso problema. Uma formulação é válida se não deixa nenhum ponto viável de fora e não inclui nenhum ponto inteiro que não seja viável. De fato, as duas formulações abaixo são válidas, uma vez que não deixam nenhum ponto vermelho de fora e não incluem nenhum ponto branco. Observe que existem infinitas possíveis formulações válidas para um problema.

A envoltória convexa e formulações válidas.

Para problemas NP-difíceis, não podemos encontrar todas as desigualdades que definem a envoltória convexa em tempo polinomial a menos que P=NP, já que do contrário poderíamos resolver o modelo como um programa linear (sem as restrições de integralidade) em tempo polinomial.

Entretanto, só o fato de aproximar a formulação da envoltória convexa já é interessante para o algoritmo de Branch and Bound. Lembrando que em cada nó é resolvida a relaxação linear do modelo, ou seja, remover suas restrições de integralidade. O valor obtido, chamado limitante dual, é usado para realizar podas por limitantes na árvore. A poda por limitantes funciona assim: se a relaxação for maior ou igual à melhor solução conhecida até agora, então não precisamos explorar sua subárvore, uma vez que o limitante dual é sempre menor ou igual ao valor de qualquer solução obtida nessa subárvore.

Aproximar a formulação da envoltória, tende a aumentar o valor da relaxação linear nos nós, que por sua vez tende a aumentar o número de podas por limitantes, melhorando o desempenho do algoritmo. Para isso precisamos encontrar desigualdades válidas e que aproximem o modelo da envoltória convexa. Porém, o número  de dessas desigualdades pode ser muito grande, tornando a resolução da relaxação linear muito lenta, sem contar que algumas das desigualdades podem ser irrelevantes para diminuir o valor do limitante dual. A ideia então é adicionar tais desigualdades apenas por demanda, através de algoritmos de plano de corte.

Algoritmo de Planos de Corte

Como pode ser visto na figura abaixo, o ponto ótimo da relaxação linear (azul) do polígono verde, não possui coordenadas inteiras. Assim, não é uma solução válida para o problema. A ideia é utilizar retas que “cortem” o ponto fora. Essa reta (ou hiperplano em dimensões maiores), chamada de plano de corte, nada mais é do que uma desigualdade que é satisfeita por todos os pontos viáveis (vermelhos) e não é satisfeita pelo ponto que queremos remover (azul). Às vezes, decidir se uma desigualdade representa um plano de corte é NP-difícil! Na prática, usam-se heurísticas.

Plano de corte remove a solução ótima

O que se faz em geral é encontrar famílias de desigualdades válidas a priori e, ao resolver a relaxação linear em cada nó do Branch and Bound, verificamos se a solução ótima encontrada viola alguma dessas desigualdades. Em caso positivo, inserimos a nova restrição no modelo e resolvemos a relaxação linear novamente, repetindo o processo enquanto conseguirmos melhor o limitante dual. Ao mudar de nó no algoritmo, descartamos os planos de corte adicionados.

Um corte pode não ser muito eficaz, porém. Se a desigualdade remover apenas a solução ótima, não ajudou muito a aproximar a formulação da envoltória convexa. Dado um plano de corte, ele pode ser classificado de acordo com a dimensão da sua interseção com a envoltória convexa. Na figura abaixo, vemos os três tipos de desigualdades para duas dimensões: a interseção de l1 com a envoltória é nula e tem dimensão -1; l2 com a envoltória é um ponto e tem dimensão 0; e l3 com a envoltória forma uma reta, com dimensão 1.

Possíveis planos de corte

É intuitivo pensar que os cortes mais eficazes são aqueles de dimensão maior (no exemplo, a reta l3). A prática também mostra que tais cortes são melhores. Assim, sendo n a dimensão da envoltória convexa, estamos interessados em encontrar desigualdades cuja interseção com a envoltória tenham dimensão n-1, que são justamente aquelas desigualdades que definem facetas dessa envoltória.

O Branch and Bound onde aplicamos algoritmos de plano de corte é chamado de Branch and Cut. No próximo post da série falarei sobre como determinar se uma desigualdade representa uma faceta da envoltória convexa.

Leituras adicionais e referências

[1] Disciplina de Programação linear inteira com o prof. Cid.
[2] Wolsey, L. – Integer Programming [link]