GRASP para o PQA: um limite de aceitação para soluções iniciais
1. Introdução
O crescimento exponencial do número de soluções com o tamanho da instância do
PQA, torna inviável a procura direta de uma solução de valor ótimo. Nestes
casos, métodos heurísticos são empregados para encontrar soluções sub-ótimas de
qualidade aceitável. A literatura é rica em estudos dedicados ao PQA, dentre os
quais podemos citar dois livros, [PW94] e [Çe98], com extensas referências
bibliográficas. Neste trabalho é realizada uma modificação na heurística GRASP
aplicada ao Problema Quadrático de Alocação (PQA) desenvolvida por Li, Pardalos
e Resende [LPR94] com uma proposta para aceitar ou não a solução inicial gerada
na fase de construção evitando uma busca que, eventualmente, exigiria muito
esforço computacional. Tal proposta está baseada no cálculo dos custos
normalizados em um intervalo de limites das soluções do PQA. No desenvolvimento
do tema, o GRASP é apresentado na Seção 2. A Seção 3 discute uma formulação do
problema e a normalização dos custos. A Seção 4 consiste na descrição do GRASP
aplicado ao PQA, conforme Li et al. [LPR94]. A proposta do trabalho se
desenvolve na Seção 5, seguida pela discussão da implementação e resultados dos
experimentos aqui realizados. Por fim, apresentam-se as conclusões.
2. O Algoritmo GRASP
Em linhas gerais, o GRASP consiste em um método iterativo probabilístico, onde
a cada iteração é obtida uma solução do problema em estudo. Cada iteração GRASP
é composta de duas fases: a construtiva, que determina a solução que será
submetida à busca local, segunda fase do algoritmo, cujo objetivo é tentar
obter alguma melhoria na solução corrente. A seguir, o algoritmo GRASP em
pseudo-código é apresentado.
Procedimento GRASP( );
1 DadosEntrada( );
2 Enquanto"critério de parada não for satisfeito" faça
3 ConstSolInicGulosaAleatória(sol);
4 BuscaLocal(sol,Viz(sol));
5 AtualizSol(sol,melhorSolEnc);
6 FimEnquanto;
7 Retorna(melhorSolEnc);
FimGRASP;
Na maioria das aplicações, o critério de parada é baseado no número máximo de
iterações. Podem-se definir outros critérios, como por exemplo parar quando a
solução procurada for encontrada ou estabelecer um tempo máximo de execução.
Na fase de construção, uma solução viável é construída elemento a elemento. Os
melhores candidatos a compor a solução são ordenados em uma lista chamada de
lista restrita de candidatos(LRC). A escolha do próximo elemento é dita
adaptativa, pois é guiada por uma função gulosa que mede, de forma míope, o
benefício que o mais recente elemento adiconado à solução concede à parte já
construída. O GRASP possui uma componente probabilística, em vista da escolha
aleatória na lista de candidatos (em um guloso simples seria selecionado o
primeiro elemento da lista). Esta técnica de escolha permite que diferentes
soluções sejam geradas a cada iteração GRASP. Mostramos, em pseudo-código, a
fase de construção da heurística.
ProcedimentoConstSolInicGulosaAleatória(sol);
1 sol = { };
2 Enquanto"solução não estiver completa" faça
3 ConsLRC(LRC);
4 s = SelecAleatElem(LRC);
5 sol = sol È{s};
6 FuncAdapGul(s);
7 FimEnquanto;
Fim ConstSolInicGulosaAleatória;
Uma vez obtida uma solução s, consulta-se a estrutura de vizinhança Viz(s)
relativa à essa solução s. Uma solução é dita localmente ótima se não existir
nenhuma solução melhor em Viz(s). As soluções iniciais do GRASP não são
necessariamente ótimos locais. Como consequência, faz-se necessária a aplicação
de um procedimento de busca local para tentar melhorar as soluções advindas da
fase construtiva. Esta busca realiza sucessivas trocas da solução corrente,
sempre que uma melhor solução é encontrada na vizinhança. Este procedimento
termina quando nenhuma solução melhor é encontrada.
Acompanhando o pseudo-código a seguir, pode-se entender a idéia de uma busca
local genérica.
ProcedimentoBuscaLocal(sol,Viz(sol));
1 Enquanto"solução não é localmente ótima" faça
2 Encontrar uma melhor solução tÎViz(s);
3 s = t;
4 FimEnquanto;
Retorna(s como localmente ótima);
FimBuscaLocal;
O procedimento de otimização local pode exigir um tempo exponencial se a busca
partir de uma solução inicial qualquer, embora se possa constatar empiricamente
a melhoria de seu desempenho de acordo com a qualidade da solução inicial. O
tempo gasto pela busca local pode ser diminuído, portanto, através do uso de
uma fase de construção que gere uma boa solução inicial. É claro que uma
estrutura de dados eficiente e uma implementação cuidadosa são importantes.
Li, Pardalos e Resende [LPR94] submeteram o GRASP a 88 instâncias da biblioteca
QAPLIB [BKR97]. Esta metodologia encontrou as melhores soluções até então
conhecidas, superando-as em alguns casos. O algoritmo GRASP foi utilizado com
diversos problemas, tais como o problema de cobertura [FR95] e do planejamento
e escalonamento de produção [LG91], além de problemas de partição em grafos
[LFE93] e de localização [Kli92].
A técnica descrita neste trabalho, baseada em limites de aceitação das solucões
iniciais no GRASP, pode ser aplicada a qualquer problema de otimização
combinatória, desde que sejam conhecidos limites inferiores e superiores para
as soluções.
3. O Problema Quadrático de Alocação
Dados o conjunto N = {1,2,...,n} e as matrizes (n x n) simétricas de inteiros
não-negativos com diagonais principais nulas F = (fij) e D = (dkl), o Problema
Quadrático de Alocação pode ser estabelecido como
onde Pn é o conjunto de todas as permutações dos elementos de N.
Na Teoria da Localização encontramos uma das principais aplicações do PQA.
Trata-se de um problema de lay-out, onde são dadas n localidades, com a
respectiva matriz de distância D = (dkl) e n facilidades, com sua
correspondente matriz de fluxo F = (fij). O custo de alocar, simultaneamente,
as facilidades i,j nas localidades k,l é o produto fijdkl..O objetivo é
encontrar uma atribuição onde todas as facilidades sejam alocadas à todas as
localidades. Deste modo, desejamos determinar uma permutação j Î P n tal que j
(i)® k e j (j)® l de custo mínimo.
Dentre os problemas de otimização combinatória, o PQA é um dos mais difíceis no
que diz respeito à complexidade computacional. Ele pertence a classe de
problemas NP-Hard. Por isso, muitos métodos heurísticos tem sido desenvolvidos
na tentativa de resolver o PQA, com eficiência e rapidez [BCMP96, CSK93, FF94,
QAB97, Sk90, Wi87].
Consideremos a instância de Gavett e Plyter [GP66], dada pelas matrizes
A figura_3.1 mostra as cliques correspondentes KF e KD e a figura_3.2, a
atribuição ótima dada pela permutação com
custo C(j*) = 403. Denotaremos a permutação j pela imagem sua j (i) = j, assim,
j*= (2 4 3 1).
Sejam o conjunto E = {(i,j) / 1£i<j£n} e N = Cn,2 , a bijeção Yij = n(i-1)-i
(i+1)/2+1, onde Yij rotula as arestas das cliques, a partir dos pares de seus
vértices, de modo que os vetores definidos como Ft = (fz) e Dt = (dz), para z =
1,...,N tem suas coordenadas dispostas segundo a ordem lexicográfica das
arestas de KF e KD, quando z = Y-1(i,j), para (i,j)ÎE.
As coordenadas da matriz quadrada de ordem N,Q = FDt, contém todos os
coeficientes da função objetivo (3.1), para qualquer j Î P n.
O Problema de Alocação Linear (PAL) associado a Q contém todas as soluções
viáveis para a instância correspondente PQA. Uma solução r Î P N do PAL (uma
permutação das colunas sob as linhas de Q), corresponde a uma atribuição das
arestas de KF sobre KD, nem sempre é compatível com uma atribuição de vértices
de uma clique a vértices da outra. Lawler [La63], Querido, Abreu e Boaventura
[QAB97] entendem tal formulação linear como uma relaxação do PQA.
Tomando-se para os vetores (F-)t e (D+)t as respectivas coordenadas de Ft em
ordem não-crescente e das Dt em ordem não-decrescente, definimos a matriz Q* =
(F-)(D+)t e a família QA(Q*), como o conjunto de todas as instâncias P do PQA,
tendo Q* em comum. O valor de tr(Q*), traço de Q*, é um limite inferior
universal (LimInf) e a soma das coordenadas da diagonal secundária determina o
limite superior universal (LimSup) para os custos de todas as instâncias em QA
(Q*), [Wi58] e [HLP52].
A instância [GP66] tem Q representada na tabela_3.1 e Q* na tabela_3.2.
O custo ótimo da instância PAL de matriz Q é igual a tr(Q*) que representa
LimInf = 389 para a instância PQA correspondente, enquanto LimSup= 587.
Define-se o custo normalizado êCusto(j ) êde uma solução j do PQA por:
Como exemplo, consideremos a solução ótima j* = (2 4 3 1) de Custo(j*) = 403 e
êCusto(j*) ê = 0.0707. Na figura_3.3, podemos verificar que êCusto(j*) êestá
bem próximo do LimInf. Para comparação, seja uma solução viável j de Custo(j) =
479 e custo normalizado êCusto(j) ê= 0.4545, o que nos parece estar
"relativamente" longe do LimInf.
4. Aplicação do GRASP ao PQA
O algoritmo GRASP é formado essencialmente por quatro componentes básicos,
todos utilizados em cada iteração do algoritmo, que constrói uma solução e a
submete, como solução inicial, a busca local. Os componentes são: (a) uma
função gulosa; (b) uma estratégia de busca adaptativa; (c) um processo
probabilístico de seleção de elementos; (d) uma técnica de busca local.
Quando aplicada ao PQA, a fase de construção do GRASP pode ser dividida em duas
etapas: a primeira cria uma correspondência entre dois elementos da permutação
a ser gerada e a segunda constrói, passo a passo, a correspondência para os (n
- 2) elementos restantes.
O primeiro par de associações é escolhido a partir de uma lista ordenada de
custos resultantes da correspondência entre os vértices do grafo de
localidades, com os do grafo de facilidades. A ordenação é guiada por uma
função gulosa, que associa facilidades com alto fluxo a localidades próximas.
Para determinação dos demais pares, a escolha do próximo elemento de menor
custo é feita a partir do último (pertencente à solução) já escolhido. Isto
caracteriza uma estratégia adaptativa.
4.1. Fase de Construção
Uma solução viável é construída iterativamente em 2 (dois) estágios. No estágio
1,são determinados os primeiros dois pares localidade-facilidade da solução
inicial. A lista restrita de candidatos (LRC) é construída, considerando-se
para tal, um parâmetro 0<b<1. São ordenadas, em ordem não-decrescente, as n2 -
n entradas da matriz de distâncias D e escolhidas as ëb (n2-n)û menores. Daí
Da mesma forma, as n2 - n entradas da matriz de fluxos F são ordenadas em ordem
não crescente e são escolhidas as ëb(n2-n)û maiores,
Finalmente, os produtos das distâncias pelos fluxos são calculados e ordenados
em ordem crescente. São considerados os ëab(n2-n)û menores. Observe-se que a,
0<a<1, é o segundo parâmetro utilizado na construção da LRC.
Sendo constantes os dados de entrada, a ordem dos custos na LRC permanece
inalterada e todo o processo de ordenação é realizado apenas uma vez.
Em cada iteração do GRASP é selecionado, de forma aleatória, um par localidade-
facilidade dentre os ëab (n2-n)û com menores custos, pertencentes à lista e
acessar os índices do custo dijfkl escolhido, correspondentes aos pares de
associações {(i,k),(j,l)}. Esses pares constituem parâmetros de entrada para o
segundo estágio construção da solução inicial.
Chamaremos de W o conjunto de pares localidade-facilidade correspondentes à
solução inicial do problema.
O estágio 2 começa com |W| = 2. Nesta etapa, para cada par de associações (i,k)
é determinado o valor de cik dado por:
Observe-se que cik corresponde ao custo da facilidade i com respeito à
localidade k, considerando-se as últimas associações feitas. Dentre todas as
possibilidades, escolhe-se aquela que minimize o valor do somatório.
Considerando m o número de possíveis pares (i,k) a serem ainda escolhidos e o
parâmetro a citado anteriormente, a lista restrita de candidatos limita a busca
de (i,k) aos ëamû pares localidade-facilidade com custos cik mínimos.
Determinado o par (i,k), o conjunto W é então atualizado para W ¬ W È {(i,k)}.
O algoritmo determina os (n-2) pares localidade-facilidade a partir dos
definidos no estágio 1. Para a determinação de cada par, são armazenados em um
vetor am candidatos e uma posição é escolhida aleatoriamente, correspondente ao
par que será colocado em W, até que toda a solução inicial do problema tenha
sido construída. Esta solução é então submetida à segunda fase do GRASP,
correspondente à busca local.
O algoritmo que descreve toda a etapa 2 da fase de construção é o seguinte:
4.2. Busca Local
A idéia central de um algoritmo de busca local é procurar iterativamente uma
solução de melhor custo dentre todas pertencentes à vizinhança da solução
corrente. No que diz respeito ao PQA, são propostas na literatura diversas
estratégias de busca local.
Dadas duas permutações quaisquer j1 e j2, a diferença entre elas é definida
comosendod(j1,j2) = {i | j1(i)¹ j2(I)}e a distância entre j1 e j2 é definida
por d(j1,j2) = | d(j1,j2) |. A estrutura de vizinhança mais comumente usada
chama-se k-troca. Uma vizinhança k-troca para uma permutação j1 Î Pné definida
por:
Vizk(j1) = {j2 | d(j1,j2) £ k }, onde 2 £ k £ n
Neste trabalho foi utilizada a vizinhança 2-troca. No exemplo, para n = 4 se
tem |Viz(j0)| = C4,2 = 6. Seja j0 = (3 1 4 2), temos que Viz(j0) = { (1 3 4 2),
(4 1 3 2), (2 1 4 3), (3 4 1 2), (3 2 4 1), (3 1 2 4)}.
4.2.1. Estratégias de Busca Local
O procedimento de busca local no GRASP desenvolvido em [LPR94] a busca local é
iniciado a partir de uma permutação j0 gerada aleatoriamente na fase de
construção, tal como descrito acima. A cada iteração do GRASP, uma permutação
j1 de melhor custo é pesquisada na vizinhança da permutação corrente. Se esta
existir, j0 (corrente) é substituída por j1. Vendo essa busca como a construção
de uma árvore, podemos percorrê-la explorando sua largura e a profundidade.
Na estratégia de Busca em Largura e Profundidade, constrói-se uma árvore por
níveis. A partir de j0 é construída uma vizinhança com estrutura 2-troca.
Dentre as permutações que pertencem ao nível 1, escolhe-se a que fornece menor
custo. Se o custo for menor que j0, gera-se o nível 2 (permutações geradas por
2-troca sobre a melhor solução do nível 1) e novamente se escolhe a permutação
de menor custo. Este processo se repete até que não haja melhora no custo com
relação ao nível anterior. A profundidade da árvore é definida pelo custo da
permutação (Fig._4.1). Esta estratégia fornece resultados satisfatórios e tem
sido bastante utilizada.
5. O GRASP Restrito
Conforme foi dito, um algoritmo de busca local depende da escolha apropriada de
uma estrutura de vizinhança, uma técnica eficiente de busca e uma solução
inicial de boa qualidade. O procedimento de construção do GRASP descrito na
seção 4.1 se preocupa com este último ponto citado. Na primeira fase, escolhe-
se aleatoriamente um elemento da LRC (de dimensão ëab(n2-n)û ) de maneira a
fornecer as duas primeiras componentes da permutação. Para as (n - 2)
componentes restantes, a segunda fase minimiza o custo acumulado das
componentes que estão sendo incluídas. Os parâmetros 0 < a,b< 1 são importantes
neste processo pois "filtram" as melhores opções de escolha.
O ajuste destes parâmetros não é uma tarefa fácil. A a tendência imediata é de
se fixar os valores de a e b bem pequenos para restringir a escolha aos
melhores candidatos fazendo com que as soluções geradas sejam de boa qualidade.
No entanto, o espaço que tais soluções varrem é pequeno. Relaxando-se os
valores desses ,os parâmetros aumenta-se o espaço de busca, mas deteriora-se a
média das soluções iniciais. Na tabela_5.1 podemos analisar este ajuste de
parâmetros com algumas instâncias da biblioteca QAPLIB. Para melhor comparação,
a média dos custos das soluções iniciais (3000 iterações) foi normalizada
(seção 3) e os cáculos feitos para a = b = 0.1, 0.25, 0.5, 0.75, 1.0.
Para melhor aproveitamento das soluções geradas por parâmentros mais relaxados,
podemos limitar a aplicação do procedimento de busca às soluções mais
promissoras advindas da fase de construção. Como exemplo, tomemos a instância
Chr12b [BKR97] com a = b = 0.75, resultando média normalizada de 0.46 e o custo
normalizado da solução ótima igual a 0.048. Propomos aceitar as soluções
geradas com custo abaixo de um determinado limite, por exemplo Lim = 0.45. A
figura_5.1 ilustra as soluções aceitas cujo custo normalizado é inferior ou
igual a Lim = 0.45.
Com esta poda no espaço das soluções iniciais tem-se um menor tempo de execução
do algoritmo pois se descartam várias tentativas de busca que não seriam
interessantes, em vista da possibilidade de um longo caminho até o ótimo.
6. Resultados Computacionais e Conclusões
Este trabalho foi implementado em linguagem C e executado em estações de
trabalho Digital DEC-3000, com OSF1, 125 Mhz e 32 MRAM. Para comparação dos
resultados, foram implementados o GRASP desenvolvido por Li et al [LPR94] e o
GRASP Restrito, ambos utilizando a mesma linguagem, os mesmos dados de entrada
e a mesma estrutura de dados. Esta restrição foi testada com limites Lim = 0.3,
0.45 e 0.5. Ambos os algoritmos terminam quando a solução ótima é encontrada ou
quando se atinge o número máximo de iterações, aqui fixado em 3000. Para
algumas instâncias, 3000 iterações não foram suficiente para atingir o ótimo
(QAPLIB). Para análise de desempenho tais execuções foram desconsideradas.
Os testes foram realizados com um conjunto de 18 instâncias (enumeradas na
tabela_5.1), com os parâmetros a = b = 0.1, 0.25, 0.5, 0.75, 1.0. e Lim = 0.3,
0.45 e 0.5 totalizando 360 execuções. O GRASP Restrito apresentou melhor
desempemho com a = b = 0.5 e 0.75 e Lim = 0.45 e 0.5. Quando Lim = 0.3 o
algoritmo descartou muitas soluções, o que prejudicou o seu desempenho, pois
houve necessidade de gerar mais soluções iniciais.
Com base nesse estudo utilizando instâncias de dimensão até 20, o algoritmo foi
submetido à instâncias de dimensões maiores. Foram executados 32 testes no
conjunto de 8 instâncias variando os parâmetros a = b = 0.5 e 0.75, Lim = 0.45
e 0.5 e com número máximo de iterações igual a 500
A tabela_6.1 fornece a média normalizada das soluções iniciais geradas pelos
parâmetros a = b = 0.5 e o número de soluções descartadas pelos limites de
aceitação Lim = 0.45 e 0.5. Os resultados alcançados considerando a = b = 0.75
foram semelhantes aos anteriores e estão na tabela_6.2.
Com relação ao tempo computacional, o algoritmo proposto neste trabalho
apresentou um desempenho melhor. Os gráficos_6.1 e 6.2 mostram a economia do
tempo computacional. As instâncias Ste36a e Ste36b, cujas médias normalizadas
das soluções iniciais são baixas (tabela_6.1 e 6.2) não apresentaram bons
resultados.
A qualidade da solução encontrada pelo GRASP Restrito é tão boa quanto a
encontrada pelo GRASP (gráfico_6.3 e 6.4). Em alguns casos houve pequeno
aumento da diferença para o valor da melhor solucão viável conhecida, o que é
compensado pela significativa economia no tempo computacional.
Os resultados obtidos nos testes computacionais utilizando o Lim = 0.5 não
foram vantajosos para o GRASP Restrito devido às médias normalizadas das
soluções iniciais . Pode-se observar nas tabelas_6.1 e 6.2 que o número de
soluções descartadas é pequeno.
Para melhor escolher o valor do limite, deve-se levar em conta a média
normalizada das soluções iniciais geradas pela fase de construção.
Utlizando o algoritmo proposto em instâncias de dimensões n = 100, os
resultdados obtidos foram bastante satisfatórios. As instâncias escolhidas
foram sko100a-f, disponíveis na QAPLIB. Foram executadas 100 iterações para
cada instância. As médias das soluções iniciais estiveram um pouco acima de
0.46 (a = b = 0.5) e o limite escolhido foi Lim = 0.47. Ver tabela_6.3.
O uso de Lim = 0.47 permitiu economizar o tempo computacional gasto sem
qualquer prejuízo à qualidade das soluções alcançadas pelo GRASP Restrito, isto
é, a melhor solução viável alcançada em cada instância foi a mesma que o GRASP
atingiu. A proximidade com a melhor solução viável conhecida (QAPLIB) foi em
média 99%.
O gráfico_6.5 mostra que a proximidade da média normalizada das soluções
iniciais com Lim = 0.47fornece uma economia diretamente proporcional no tempo
computacional
Nota-se que um bom ajuste do limite com a média normalizada das soluções
iniciais é um ponto importante para a eficiência do algoritmo proposto. Se
fosse escolhido Lim = 0.45 não haveria economia alguma e se a escolha fosse Lim
= 0.5, não se aceitaria nenhuma solução para aplicar a busca local.
Na equação 3.2, os limites usados poderiam ser substituídos por outros mais
adequados (ver Li et al [LPRR94], Carraresi e Malucelli [CM94], Gilmore [Gi62],
Lawler [La63] e Karisch et al [KÇCE98]). O custo da determinação de melhores
limites (todos de ordem igual ou superior a O(n3)), poderia a reduzir a
vantagem obtida pelo algoritmo proposto, dada a baixa complexidade da
determinação do limite aqui utilizado.
A técnica aqui apresentada pode ser aplicada em outros problemas de otimização
combinatória, desde que seja possível a determinação de limites superiores e
inferiores para os custos das soluções viáveis.
