Introdução a Constraint Programming usando Gecode

Neste post vamos falar sobre constraint programming (CP) ou programação com restrições.

Constraint Programming é uma técnica em que se define um espaço de busca, um conjunto de restrições e deseja-se encontrar soluções viáveis, isto é, pontos que pertencem ao espaço de busca e que satisfazem as restrições.

Segundo a Wikipedia, surgiu a partir da programação por restrições lógicas (constraint logic programming). Essa variante de programação lógica foi criada por Jaffar e Lassez, que introduziram restrições ao Prolog II, sendo que Prolog III possui suporte nativo a programação por restrições lógicas.

É possível trabalhar com CP em linguagens imperativas através de frameworks. Algumas das bibliotecas mais famosas são Gecode (C++), IBM ILOG CP (C++ proprietário), Choco (Java), JaCoP (Java) e Google OR-Tools (C++, com bindings para Java, Python e C#). A biblioteca que utilizaremos para implementar as soluções para os exemplos é a Gecode (v 3.7.3).

Em geral os resolvedores de CP fazem uso de backtracking para encontrar essas soluções e portanto além da modelagem em si, outro fator que influencia o desempenho é a maneira como se faz as ramificações (branching) na árvore de busca.

Um breve tutorial de Gecode

Gecode é uma biblioteca open-source escrita em C++. É disponibilizada com licença MIT, é bastante eficiente, tem suporte a paralelismo e é bem documentada [1].

Modelos no Gecode devem definir um “espaço”. Dentro desse espaço ficam as variáveis, os propagadores (implementações de restrições) e branchers (implementações de branchings).

Um espaço é representado pela classe Space, que o modelo deve estender.

Vamos agora detalhar alguns conceitos utilizados pelo Gecode. Após isso, apresentaremos alguns exemplos clássicos que pode ser atacados com programação por restrição.

Variáveis

O Gecode define novos tipos para especificar o domínio das variáveis, como por exemplo IntVar, para variáveis inteiras e BoolVar para variáveis binárias. Também há o tipo representando um vetor desses tipos, IntVarArray e BoolVarArray.

Podemos construir um vetor de inteiros da seguinte maneira:

IntVarArray (Space &home, int n, int min, int max)

Sempre precisamos passar uma referência de um Space como parâmetro para informar em qual espaço criaremos a variável. O parâmetro n define o número de elementos no vetor, min é o limite inferior e max o liminte superior dos elementos nesse vetor.

Restrições

As restrições delimitam o domínio de busca de uma solução. No Gecode existem funções que criam propagadores implementando essas restrições. Algumas das principais são:

Restrições de relação:

void rel (Home home, IntVar x0, IntRelType r, IntVar x1);

Aqui, Home é equivalente a &Space, pois, novamente, precisamos especificar em que espaço a restrição será gerada. Essa restrição define uma relação entre a variável x0 e x1 usando uma relação r.

Um exemplo de uso é

rel(*this, x, IRT_NQ, 0)

Que modela a restrição x!= 0. Esse exemplo pode ser escrito de maneira mais legível com o uso de sobreposição de operadores:

rel(*this, x != 0)

Restrições de distinção:

Uma restrição bastante útil é a distinct, que exige que todos os valores de um vetor tenham valores diferentes:

void distinct (Home home, IntVarArray x0);

Restrições com expressões regulares:

É possível especificar que o vetor de solução satisfaça uma expressão regular através da restrição extensional.

extensional(Home home, const IntVarArgs &x, DFA d)

Essa restrição espera um vetor de inteiro e um DFA (Deterministic Finite Automata ou Autômato Determinístico Finito), que pode ser criado a partir de uma expressão regular. Essa função restringe os valores de x de forma que seja aceitável pelo DFA.

O tipo REG é uma maneira conveniente de se construir uma expressão regular a partir de expressões regulares menores.

Podemos criar a “base” de uma expressão regular com um único valor ou um conjunto possível de valores:

REG(int i);
REG (const IntArgs &x);

No exemplo abaixo, r corresponde a “1” e s corresponde a “(0|1)”.

REG r(1);

int v[] = {0, 1};
REG s(IntArray(vector<int>(v, v+2)));

Os operadores * e + foram convenientemente sobrescritos para modelar respectivamente “0 ou mais repetições” e “1 ou mais repetições”, por exemplo

REG r(1);
REG s = *r;
REG t = +r;

A variável s corresponde à expressão regular “1*” e t à “1+”.

O operador () também foi sobrescrito com duas versões: recebendo um inteiro n como parâmetro, que indica que a expressão da variável deve aparecer pelo menos n vezes; e a outra versão recebendo n e m, indicando que a expressão deve aparecer entre n e m (inclusive) vezes. Assim, no exemplo

REG r(1);
REG s = r(10);
REG t = r(3, 4);

A variável s corresponde à “1{10,}*” e t corresponde a “1{3,4}”.

Note que uma expressão regular seguida por “{n}” deve aparecer exatamente n vezes (e não pelo menos n vezes). Portanto, para modelar esse caso usando REG, devemos fazer

REG s = r(n, n);

Veremos a utilidade desse tipo de restrição no exemplo do Nonograma!

Cópia e clone

Como o algoritmo usa backtrack, ele precisa manter cópias para eventualmente voltar a uma solução anterior. Para isso precisamos implementar o método puramente virtual da classe Space, o copy(), retornando um ponteiro para a cópia criada.

virtual Space* copy(bool share) = 0;

A flag share está relacionada com o uso de múltiplas threads, mas não vamos nos preocupar com ela nesse post introdutório. Na cópia, o seguinte construtor da classe Space deve ser chamado:

Space::Space(boolean clone, Space &space);

Para fazer o clone de uma variável, podemos usar o método update(), presente na classe base dos tipos IntVar e IntVarArray:

template<class VarImp>
void Gecode::VarImpVar< VarImp >::update(Space &home, bool share, VarImpVar< VarImp > &y) 

Ele equivale a um clone. Existe um construtor de IntVar que recebe uma referência para outra IntVar, mas ele não cria um cópia e sim uma referência. Por exemplo,

IntVar x(y);

Aqui x e y referem-se à mesma variável.

Branching

Precisamos especificar o modo como será feito o backtracking através do método branch().

branch (Home home, const IntVarArgs &x, IntVarBranch vars, IntValBranch vals)

O primeiro parâmetro especifica o espaço onde a regra de branching será criada, x é um vetor de variáveis sobre o qual faremos branch, vars indica define como será feita a escolha das variáveis. O parâmetro vals indica como será feita a escolha dos valores.

Nos exemplos usaremos sempre a seguinte configuração:

branch(*this, sol, INT_VAR_SIZE_MIN, INT_VAL_MIN);

Aqui, sol representa a solução, INT_VAR_SIZE_MIN indica que o backtracking sempre fará branching sobre as variáveis com menor domínio primeiro e o parâmetro INT_VAL_MIN indica que os valores “chutados” serão feitos em ordem crescente.

Resolvedores

No Gecode há três tipos de motores de busca (search engines), um que usa busca em profundidade (chamado de DFS), outro que usa Branch-and-Bound (BAB) e um chamado Depth First Restart, que não sei exatamente o que é, mas talvez seja uma busca em profundidade iterativa (Iterative deepening depth search).

Os dois últimos motores são usados para quando se deseja otimizar uma função objetivo. Nos nosso exemplos estamos preocupados apenas em encontrar uma solução para o problema, então vamos usar apenas o DFS.

O resolvedores são um template para um tipo genérico e podem ser instanciados usando uma instância desse tipo. No nosso primeiro exemplo vamos implementar a classe SendMoreMoney e portanto podemos instanciar o motor de busca fazendo:

SendMoreMoney *m = new SendMoreMoney();
DFS<SendMoreMoney> e(m);

Os resolvedores definem o método next(), que retorna uma solução do problema (que é uma instância do modelo) e avança para a próxima. Se nosso modelo definir a função de imprimir uma solução, print(), podemos imprimir todas as soluções encontradas da seguinte maneira:

while(SendMoreMoney* s = e.next()){
    s->print();
    delete s;
}

A seguir, vamos comentar sobre 4 problemas clássicos e apresentar soluções usando CP. Os códigos-fontes estão disponíveis no Github e são adaptações dos códigos providos pela própria Gecode.

1. Puzzles Alpha-Numéricos

Um puzzle que pode ser resolvido com programação por restrições é aquele em que são dadas equações de palavras e cada letra da palavra deve ser substituída por um algarismo de modo que a equação seja satisfeita.

Letras iguais devem ser substituídas pelo mesmo algarismo, letras diferentes devem ter algarismos diferentes e a primeira letra de cada palavra não pode ser substituída por 0.

Um exemplo clássico é:

SEND + MORE = MONEY

Podemos resolver esse exemplo definindo 8 variáveis D, E, M, N, O, R, S, Y sujeitas às seguintes restrições:

1) As variáveis devem pertencer ao conjunto de inteiros [0-9].
2) M e S devem ser diferentes de 0.
3) As variáveis devem ter valores diferentes
4) A seguinte igualdade deve ser satisfeita

S*1000 + E*100 + N*10 + D +
M*1000 + O*100 + R*10 + E =
M*10000 + O*1000 + N*100 + E*10 + Y

Podemos implementar essas restrições da seguinte forma:

// Referencias para melhor legibilidade
IntVar 
s(sol[0]), e(sol[1]), n(sol[2]), d(sol[3]),
    m(sol[4]), o(sol[5]), r(sol[6]), y(sol[7]);

// Números nao podem comecar com 0
rel(*this, s != 0);
rel(*this, m != 0);

// Todas as variaveis devem ter valores distintos
distinct(*this, sol);

// Vetor de inteiros
IntArgs c(4+4+5);
// Vector de variaveis inteiras
IntVarArgs x(4+4+5);

// SEND
c[0]=1000; c[1]=100; c[2]=10; c[3]=1;
x[0]=s; x[1]=e; x[2]=n; x[3]=d;
// MORE
c[4]=1000; c[5]=100; c[6]=10; c[7]=1;
x[4]=m; x[5]=o; x[6]=r; x[7]=e;
// MONEY
c[8]=-10000; c[9]=-1000; c[10]=-100; c[11]=-10; c[12]=-1;
x[8]=m; x[9]=o; x[10]=n; x[11]=e; x[12]=y;

linear(*this, c, x, IRT_EQ, 0);

2. O problema das N-damas

O problema das N-damas (do inglês N-queens) é também bastante conhecido. Dado um tabuleiro de xadrez NxN, encontrar uma solução em que N damas sejam dispostas neste tabuleiro de forma que nenhuma esteja em posição de ataque com a outra. Segundo Hoffman, Loessi e Moore [2], sempre existe solução para N maior ou igual a 4.

Podemos definir a solução em outras palavras: queremos uma solução em que nenhum par de damas esteja na mesma coluna, na mesma linha ou na mesma diagonal. Claramente cada coluna deverá conter exatamente uma dama, de modo que podemos modelar a solução através de um vetor linhas de N posições, onde a linhas[i] representa a linha da rainha na i-ésima coluna.

Se colocarmos uma restrição de que os valores de linhas devam ser diferentes, garantimos que as damas estão em linhas diferentes. Entretanto, isso ainda não impede o posicionamento na mesma diagonal.

Se i e j estão na mesma diagonal principal, então dada a célula da dama i é dada por (linha[i], i) e a de j é (linha[i]+k, i+k). Se estão na diagonal invertida, então a posição de j é (linha[i]+k, i-k).

É fácil ver que i e j compartilham alguma diagonal se e somente se linha[i]-i = linha[j]-j ou se linha[i]+i = linha[j]+j.

Podemos implementar essas restrições da seguinte forma:

for (int i = 0; i<n; i++)
    for (int j = i+1; j<n; j++) {
        rel(*this, sol[i] != sol[j]);
        rel(*this, sol[i]+i != sol[j]+j);
        rel(*this, sol[i]-i != sol[j]-j);
    } 

3. Sudoku

No Sudoku clássico tem-se um tabuleiro 9×9, formado por 9 sub-tabuleiros 3×3. O objetivo é preencher o tabuleiro com algarismos de 1 a 9 de forma que nenhuma coluna, nenhuma linha e nenhum sub-tabuleiro contenha números repetidos. Além disso, na grande maioria das vezes os Sudokus já vêm com alguma de suas células preenchidas.

O Gecode oferece uma estrutura bastante simples de se trabalhar com matrizes, que é a classe Matrix. Temos os métodos row(i) e col(j) que retorna respectivamente a i-ésima linha e a j-ésima coluna.

Além disso, temos o método slice(i1, i2, j1, j2) que retorna a submatriz correspondendo a [i1, i2-1] e [j1, j2-1], o que é útil para modelar a restrição de sub-tabuleiros.

As restrições podem ser implementadas praticamente a partir da definição:

// Linhas e colunas devem ter valores diferentes
for (int i=0; i<N; i++) {
    distinct(*this, m.row(i));
    distinct(*this, m.col(i));
 }

// Restricoes dos sub-quadrados
for (int i=0; i<N; i+=n) {
    for (int j=0; j<N; j+=n) {
        distinct(*this, m.slice(i, i+n, j, j+n));
    }
 }

// Celulas ja preenchidas
for (int i=0; i<N; i++)
    for (int j=0; j<N; j++)
        if (mat[i][j] != '.')
            rel(*this, m(i,j), IRT_EQ, mat[i][j] - '0');

4. Nonograma

No nonograma é dado um tabuleiro com n linhas e m colunas. Para cada linha e coluna é dado um conjunto de números. Isso indica que naquela linha ou coluna deve blocos com exatamente aqueles comprimentos. Considere o exemplo abaixo:

Já tinha falado sobre esse puzzle anteriormente e à época não tinha conseguido pensar em uma solução para ele.

Vimos que o Gecode oferece uma maneira de especificar as restrições usando DFA’s. Se usarmos uma matriz de variáveis binárias como nosso vetor de soluções, então dada uma especificação de linha ou coluna dada por

3 4 5

Podemos criar a seguinte expressão regular:

0* 1{3} 0+ 1{4} 0+ 1{5} 0*

Basicamente, nossa linha ou coluna deve começar com 0 ou mais ‘0‘s, seguido por um bloco de exatamente 3 ‘1‘s. Entre dois blocos consecutivos de ‘1‘s deve haver pelo menos um ‘0‘.

Podemos transformar uma especificação em uma expressão regular da seguinte maneira:

// Expressao regular a partir de uma linha
DFA regexp(std::vector<int> &v) {
    REG r0(0), r1(1);
    REG border = *r0;
    REG between = +r0;
    REG r = border;
    for(size_t i = 0; i < v.size(); i++){
        if(i > 0) r += between;
        r += r1(v[i], v[i]);
    }
    return r + border;
}

Agora basta adicionar uma expressão regular dessas para cada coluna e linha. Nossa entrada é dada por um vetor de vetores com as especificações das linhas (_rows) e das colunas (_cols):

Matrix<BoolVarArray> m(sol, _cols.size(), _rows.size());
    
// Restricoes das colunas
for (size_t w=0; w < cols.size(); w++)
    extensional(*this, m.col(w), regexp(cols[w]));
// Restricao das linhas
for (size_t h=0; h < rows.size(); h++)
    extensional(*this, m.row(h), regexp(rows[h]));
// Branch   
for (int w = 0; w < cols.size(); w++)
    branch(*this, m.col(w), INT_VAR_NONE, INT_VAL_MIN);

Constraint Programming vs. Programação linear Inteira

Muitos dos problemas apresentados aqui também podem ser modelados como programação linear inteira e resolvidos através de um algoritmo de branch-and-bound.

Lembrando, no algoritmo de branch-and-bound resolve-se a relaxação linear do modelo em cada nó da árvore de busca. O valor dessa relaxação serve como limitante inferior (superior) para um problema de minimização (maximização). Com isso o algoritmo pode descartar soluções viáveis do problema sem analisá-las se o valor da relaxação em um nó for pior do que o valor da melhor solução encontrada até o momento.

Em CP não temos o valor desse limitante e portanto em um problema com função objetivo, em geral não se pode afirmar que uma solução é ótima a menos que se tenha analisado todas as soluções viáveis do problema.

Por outro lado, existem algumas modelagem que são complicadas de se fazer em programação linear inteira. A mais comum delas é a restrição de que duas variáveis devam ter valores diferentes.

Por exemplo no caso das N-damas, uma das restrições é de que cada elemento de linha tenha um valor diferente. Uma abordagem que pode ser usada para que variáveis x e y inteiras tenham valores diferentes é

(1) x >= y + 1

ou

(2) x <= y - 1

Lembrando que restrições disjuntas podem ser modeladas através da inclusão de uma variável binária, z que vale 0 se (1) for utilizada e 1 se (2) for utilizada. Dado uma constante M suficientemente grande temos:

(3) x \ge y + 1 - M \cdot z
(4) x \le y - 1 + M \cdot (1 - z)

Se z é 0, (3) vira (1) e (4) torna-se redundante. Se z é 1, (3) torna-se redundante e (4) vira (2).

Uma alternativa mais interessante é modelar os valores das variáveis inteiras como variáveis binárias. No caso do problema das N-damas, ao invés de trabalhar com linhas, poderíamos definir a variável c_{i,j} para cada célula do tabuleiro, sendo c_{i,j} = 1 se a dama está na posição (i,j) do tabuleiro e 0 caso contrário.

As restrições lineares a ser adicionadas são tais que a soma em cada linha e cada coluna seja exatamente 1 e que a soma em cada diagonal seja menor ou igual a 1.

Conclusão

Neste post fizemos uma pequena introdução a programação por restrições usando Geocode. Há muitos detalhes que faltam ser explorados nessa biblioteca como a customização de propagadores e branchers, bem como a comparação do desempenho de diferentes modelos.

Além disso, os exemplos empregados aqui são relativamente simples, com poucas restrições. A técnica de CP vem sendo aplicada a problemas do mundo real onde existem bastante detalhes que são modelados na forma de restrições. Pretendo estudar alguns papers e escrever sobre eles.

Ao decidir estudar CP, fiquei na dúvida sobre qual biblioteca utilizar. Eu pretendia usar a OR-Tools usando Python porque deve ser mais simples escrever os modelos. Entretanto, a documentação dessa bilioteca ainda está muito escassa e por isso acabei escolhendo a Gecode. De qualquer maneira Gecode também tem binding para uma linguagem simples de se escrever, que é o Ruby. Talvez seja uma boa oportunidade para aprender um pouco dessa linguagem.

Referências

[1] Gecode – GEneric COnstraint Development Environment
[2] Constructions for the Solution of the m Queens problem – E. Hoffmann , J. Loessi e R. Moore (Mathematics Magazine, Vol. 4, No. 2)

Os comentários estão fechados.

%d bloggers like this: