O padrão de projeto Visitor e Vtables em C++

Agosto 26, 2012

Outro dia estava desenvolvendo um código Java em que precisei usar instanceof, mas como já havia lido que o uso dessa palavra-chave é sinal de design ruim, busquei ajuda e me recomendaram o padrão de projeto Visitor.

Neste post vamos falar sobre esse padrão de projeto em Java, que foi onde minha dúvida se originou e também vamos discuti-lo em C++, porque sua implementação tem relação com a estrutura Vtables que eu queria estudar a algum tempo.

Visitor em Java

Consideremos o seguinte exemplo, onde temos as classes Cachorro e Gato que são especializações de Animal e que possuem métodos para emissão de som, respectivamente latir() e miar(), respectivamente.

Queremos implementar um método emitirSom(), que recebe uma instância de Animal e emite o som de acordo com a classe instanciada. Uma alternativa usando instanceof é a seguinte:

// Listagem 1
interface Animal {    
}
public class Cachorro implements Animal {
	public void latir(){
		System.out.println("Au!");
	}
}
public class Gato implements Animal {
	public void miar(){
		System.out.println("Miau");
	}
}
...
public static void emitirSom(Animal animal){
    if(animal instanceof Cachorro){
        Cachorro cachorro = (Cachorro) animal;
        cachorro.latir();
    }
    else if(animal instanceof Gato){
        Gato gato = (Gato) animal;
        gato.miar();
    }
}

Aqui cabe uma citação a Scott Meyers, no seu livro Effective C++:

Anytime you find yourself writing code of the form “if the object is of type T1, then do something, but if it’s of type T2, then do something else,” slap yourself.

Para nos livrarmos do instanceof, basta adicionarmos um método emitirSom() à interface Animal e substituir latir()/miar(). No método emitirSom() usamos polimorfismo para escolher a implementação correta.

// Listagem 2
interface Animal {
    void emitirSom();
}
public class Cachorro implements Animal {
    @Override
    public void emitirSom(){
        System.out.println("Au!");
    }
}
public class Gato implements Animal {
    @Override
    public void emitirSom(){
        System.out.println("Miau");
    }
}
...
public static void emitirSom(Animal animal){
    animal.emitirSom();
}

Outra possibilidade é delegar a implementação de emitirSom() para uma outra classe através do padrão de projeto Visitor. Esse padrão considera dois tipos de classes: um elemento e um visitante. O elemento implementa um método geralmente chamado de accept() e o visitante o método visit().

O método accept() recebe um visitante e chama o método visit() desse visitante passando a si mesmo como parâmetro. Desta forma, o visitante pode implementar o método visit() para cada tipo.

// Listagem 3
interface Animal {
    void emitirSom();
    void accept(AnimalVisitor visitor);
}
public class Cachorro implements Animal {
    @Override
	public void accept(AnimalVisitor visitor){
        visitor.visit(this);
    }
}
public class Gato implements Animal {
    @Override
	public void accept(AnimalVisitor visitor){
        visitor.visit(this);
    }
}
interface AnimalVisitor {
    void visit(Cachorro cachorro);
    void visit(Gato gato);
}
public class EmissorDeSom implements AnimalVisitor {
    @Override
	public void visit(Gato gato){
        System.out.println("Miau");
    }
    @Override
	public void visit(Cachorro cachorro){
        System.out.println("Au!");
    }	
}

A vantagem dessa abordagem é que a classe EmissorDeSom pode conter membros e métodos comuns ao modo como os animais emitem som, que de outra forma acabaria sendo implementado na classe Animal.

Outra vantagem é que a implementação do visitante pode ser escolhida em tempo de compilação e, sendo uma injeção de dependência, facilita testes.

Visitor em C++

Em C++ não temos a palavra-chave instanceof, mas podemos usar a função typeid().

Segundo [1], se o argumento dessa função for uma referência ou um ponteiro de-referenciado para uma classe polimórfica, ela retornará um type_info correspondendo ao tipo do objeto em tempo de execução.

// Listagem 4
struct Animal {
    virtual void foo(){};
};
struct Cachorro : public Animal {
    void latir(){
        cout << "Au!\n";
    } 
};
struct Gato : public Animal {
    void miar(){
        cout << "Miau\n";
    }
};
void emitirSom(Animal *animal){
    if(typeid(*animal) == typeid(Cachorro)){
        Cachorro *cachorro = dynamic_cast<Cachorro *>(animal);
        cachorro->latir();
    }
    else if(typeid(*animal) == typeid(Gato)){
        Gato *gato = dynamic_cast<Gato *>(animal);
        gato->miar();
    }
}

Novamente, podemos criar um método emitirSom() em uma interface e usar o polimorfismo. Em C++ não temos interface, mas podemos definir um método puramente virtual como no código abaixo.

// Listagem 5
struct Animal {
    virtual void emitirSom() = 0;
};
struct Cachorro : public Animal {
    void emitirSom(){
        cout << "Au!\n";
    }
};
struct Gato : public Animal {
    void emitirSom(){
        cout << "Miau\n";
    }
};
void emitirSom(Animal *animal){
    animal->emitirSom();
}

A implementação do padrão Visitor pode ser feita da seguinte maneira:

// Listagem 6
struct Gato;
struct Cachorro;

struct AnimalVisitor {
    virtual void visit(Gato *gato) = 0;
    virtual void visit(Cachorro *cachorro) = 0;
};
struct Animal {
    virtual void accept(AnimalVisitor *visitor) = 0;
};
struct EmissorDeSom : public AnimalVisitor {
    void visit(Gato *gato){
        cout << "Miau\n";
    }
    void visit(Cachorro *cachorro){
        cout << "Au!\n";
    }    
};
struct Cachorro : public Animal {
    void accept(AnimalVisitor *visitor){
        visitor->visit(this);
    }
};
struct Gato : public Animal {
    void accept(AnimalVisitor *visitor){
        visitor->visit(this);
    }
};
void emitirSom(Animal *animal){
    AnimalVisitor *visitor = new EmissorDeSom();
    animal->accept(visitor);
}

Despacho único e despacho múltiplo

Vamos agora analisar como C++ implementa funções virtuais. Antes disso, precisamos entender o conceito de despacho dinâmico.

Despacho dinâmico é uma técnica utilizada no caso em que diferentes classes implementam diferentes métodos, mas não se sabe qual o tipo do objeto que chama o método em tempo de compilação, como nos exemplos acima.

Em C++ e Java, esse tipo de despacho é também chamado de despacho único (single dispatch) porque ele só leva em conta o tipo do objeto que chama o método. Em Lisp e algumas outras poucas linguagens, há o chamado despacho múltiplo que também leva em conta o tipo de parâmetro.

Por exemplo, considere o código abaixo:

// Listagem 7
void visitor(Cachorro *cachorro){
    cout << "Au!\n";
}
void visitor(Gato *gato){
    cout << "Miau!\n";
}
void emitirSom(Animal *animal){
    visitor(animal);
}

Teremos um erro de compilação, porque não temos ligação dinâmica a partir dos parâmetros.

Agora veremos como o despacho dinâmico é implementado em C++.

C++ Vtables

Embora a especificação do C++ não diga qual a maneira de implementar o despacho dinâmico, a maioria dos compiladores fazem uso de uma estrutura chamada virtual method table, também conhecida por outros nomes como vtable. Essencialmente, essa tabela é um array de ponteiros para os métodos virtuais.

Cada classe que define ou herda pelo menos um método virtual possui sua vtable. Ela aponta para os métodos virtuais definidos pelo seu ancestral (incluindo a si mesmo) mais próximo que podem ser utilizados por aquela classe. Além disso, cada instância de uma classe com vtable possui um ponteiro para essa tabela. Esse ponteiro é setado quando instanciamos esse objeto.

Considere o trecho da Listagem 5:

Animal *animal = new Cachorro();

Aqui o ponteiro de animal aponta para a vtable de Cachorro, e não de Animal. Então, quando fazemos

animal->emitirSom();

é feita uma consulta para saber exatamente qual implementação emitirSom() chamar. Note que para métodos virtuais existe uma indireção que pode afetar o desempenho. Por isso, por padrão, o despacho dinâmico não é efetuado em C++, a menos que alguma classe adicione a palavra-chave virtual. Em Java o despacho é habilitado por padrão.

Vamos a um exemplo,

struct A {
    virtual void foo(){}
    virtual void bar(){}
    void baz(){}
};
struct B : public A {
    void bar(){}
    virtual void baz(){}
};

Podemos analisar as vtables das classes A e B compilando o código acima com gcc usando a opção -fdump-class-hierarchy. Um arquivo de texto com a extensão .class será gerado.

Vtable for A
A::_ZTV1A: 4u entries
0     (int (*)(...))0
4     (int (*)(...))(& _ZTI1A)
8     (int (*)(...))A::foo
12    (int (*)(...))A::bar

Class A
   size=4 align=4
   base size=4 base align=4
A (0xb71f6038) 0 nearly-empty
    vptr=((& A::_ZTV1A) + 8u)

Vtable for B
B::_ZTV1B: 5u entries
0     (int (*)(...))0
4     (int (*)(...))(& _ZTI1B)
8     (int (*)(...))A::foo
12    (int (*)(...))B::bar
16    (int (*)(...))B::baz

Class B
   size=4 align=4
   base size=4 base align=4
B (0xb7141618) 0 nearly-empty
    vptr=((& B::_ZTV1B) + 8u)
  A (0xb71f61f8) 0 nearly-empty
      primary-for B (0xb7141618)

Aqui podemos ver que a função A::foo e A::bar aparecem na vtable de A, mas a função A::baz não, porque ela não foi declarada virtual. Na vtable de B temos A::foo, porque ela não foi sobrescrita em B. Temos também B::bar, embora ela não tenha sido declarada virtual em B, essa propriedade foi herdada de A::bar. Finalmente, B::baz aparece na vtable porque foi declarada como virtual em B.

Podemos ver também o valor que os ponteiros que as instâncias dessas classes terão: vptr=((& A::_ZTV1A) + 8u) e vptr=((& B::_ZTV1B) + 8u), que é respectivamente onde os ponteiros para as funções começam nas respectivas tabelas.

Uma referência interessante para compreender melhor vtables pode ser vista em [2].

Referências

[1] The typeid operator (C++ only)
[2] LearnCpp.com – The Virtual Table


Introdução a Constraint Programming usando Gecode

Agosto 19, 2012

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 &amp;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)


Relaxação Lagrangeana – Prática

Março 11, 2012

Em um post anterior, expliquei brevemente o conceito de relaxações, focando em relaxação lagrangeana. Hoje vamos apresentar uma heurística que encontra multiplicadores lagrangeanos.

Aplicaremos o algoritmo no problema da localização de instalações.

Algoritmo

O algoritmo mais conhecido para encontrar bons multiplicadores lagrangeanos é dado por Beasley [1], e é conhecido por método do subgradiente.

Geometricamente, se tivéssemos dois multiplicadores lagrangeanos, estes representariam pontos em um plano e os valores das soluções ótimas do dual lagrangeano para um desses pontos (multiplicadores) formariam curvas de níveis. O que o algoritmo faz é procurar um ponto com o valor mais próximo do da solução ótima primal e para isso anda em uma direção (gradiente) em que há melhoria dessa função.

Se não houver melhoras por um longo tempo, pode ser que o “passo” que estamos dando a cada iteração esteja muito grande e pode ser que haja uma solução melhor do meio desse passo. Por isso, vamos diminuímos gradativamente o passo, até que eventualmente este fique menor que um tamanho mínimo.

Onde a função Atualiza_multiplicadores é dada por:

Neste algoritmo é usado o valor de uma solução primal, que é obtida a partir da relaxação através de uma heurística.

Aplicação ao problema da Localização de Instalações

Vamos usar a relaxação das desigualdades que forçam os clientes a serem cobertos por pelo menos uma instalação (como descrito em um post anterior).

Heurística lagrangeana. Como estamos trabalhando com a versão sem restrição de capacidade, dado um conjunto de fábricas abertas é trivial encontrar a solução de melhor custo, bastando ligar cada cliente à instalação aberta mais barata. É essa heurística que usaremos para obter uma solução primal a partir de uma solução relaxada.

Multiplicadores lagrangeanos. Note que o algoritmo requer um conjunto de multiplicadores lagrangeanos inicial. No nosso caso, fizemos u_i = \frac{1}{n}, onde n é o número de desigualdades relaxadas.

Para os outros parâmetros do algoritmo usamos: \pi = 2, T_\pi = 0.2 e N_{stuck} = 30.

Fiz uma implementação em C++, tentado separar o problema (classe FacilityLocation implementando a classe virtual LagrangeanModel) do algoritmo (classe GradientMethod). Como de costume, o código está no github.

Resultados computacionais

Para testar o algoritmo, usei as instâncias do problema Capacitated Facility Location da OR-Library [3]. Como estamos lidando com a versão não restrita (uncapacitated), ignorei os dados de capacidade.

Há também este arquivo com os valores ótimos das soluções de forma que podemos observar a qualidade das soluções obtidas via relaxação empiricamente. A tabela a seguir contém o nome da instância, o número de instalações, o número de clientes, a distância percentual entre o ótimo e o melhor valor primal obtido (Gap), a distância percentual entre os melhores valores valores primal e dual (Gap DP) e finalmente o tempo.

Observamos que a solução ótima foi encontrada para boa parte das instâncias. O pior resultado foi para a instância capb, devido à dificuldade do algoritmo em encontrar bons multiplicadores lagrangeanos (note a distância entre os valores primal e dual).

Note que mesmo para instâncias grandes, provavelmente inviáveis de serem resolvidas de maneira exata, são resolvidas em menos de 6 segundos (processador i7 2.2GHz, sem paralelismo).

Conclusão

Embora tenhamos que lidar com uma formulação do problema para decidir quais restrições relaxar e obter o problema relaxado, geralmente para resolvê-lo não precisamos trabalhar com a formulação explicitamente (em geral resolver PL’s, principalmente de formulações grandes, não é muito rápido). Assim, dependendo da complexidade do problema resultante, a relaxação lagrangeana pode ser bastante eficiente.

Vimos que na prática as soluções obtidas via heurística lagrangeana podem ser muito boas. O principal problema é que a escolha dos parâmetros podem influenciar bastante a qualidade dessas soluções.

Eu já havia feito um trabalho sobre relaxação lagrangeana durante uma disciplina na faculdade, sobre um problema chamado Axis Parallel Minimum Stabbing Number Spanning Tree [2].

Referências

[1] J. E. Beasley – Modern Heuristics for Combinatorial Optmization Problem, Chapter 6 (Lagrangean Relaxation). Blackwell Scientific Publications, 1993.
[2] Relatório Programação Linear Inteira (MO420) – Guilherme Kunigami
[3] OR-Library


Currying de funções em Haskell

Novembro 6, 2011

Currying de função é uma técnica de avaliação parcial dos argumentos de uma função. O nome da técnica é em homenagem a Haskell Curry, que segundo a Wikipedia re-inventou o método (Moses Schönfinkel já havia desenvolvido o método anteriormente). O curioso é que esse matemático também deu nome à linguagem Haskell.

Nesse post vou falar um pouco sobre currying de funções em Haskell, baseando-me principalmente no capítulo 4 do livro “Real World Haskell“.

Funções com múltiplos argumentos em Haskell

A sintaxe para especificarmos o tipo de uma função em Haskell com múltiplos argumentos parece meio confusa no começo. Por exemplo, se quisermos definir uma função que soma dois inteiros, devemos especificá-la da seguinte forma.

 
soma:: Int -> Int -> Int
soma a b = a + b

A partir do cabeçalho não dá a impressão de que a função recebe dois inteiros como argumento e retorna um inteiro. Isso porque as funções em Haskell recebem na verdade apenas um argumento. Quando há múltiplos argumentos, chamamos a função ‘f’ com o primeiro argumento e retornamos uma função ‘g’ com as ocorrências desse primeiro argumento substituídas pelo valor passado como argumento.

Por exemplo, considere uma função f(a, b, c) sobre três inteiros a, b e c. A primeira chamada da função f(1,2,3) retorna uma nova função g(b, c) = f(1, b, c), que por sua vez retorna uma outra função h(c) = g(2, c), que recebe apenas um argumento.

Para testar isso na prática, podemos fazer no terminal (ghci):

> :type soma
soma :: Int -> Int -> Int
> let somaTres = soma 3
> :type somaTres
soma3 :: Int -> Int
> somaTres 7
10 

No exemplo acima, o parâmetro a de soma, foi substituído por 3 e a função resultante atribuída a somaTres, que agora só recebe um argumento.

Em uma função com vários parâmetros, podemos trabalhar com versões intermediárias do processo de currying provendo apenas alguns dos parâmetros. Por exemplo,

-- Função que retorna uma tripla a partir de 3 argumentos
foo a b c = (a, b, c)
foo1 = foo 3
foo2 = foo 4 5
foo3 = foo 6 7 8

Nesse caso, foo1 recebe dois argumentos, foo2 apenas um e foo3 já é o resultado da função foo.

Omissão de parâmetros usando currying

Suponha que queiramos selecionar os números primos de uma lista de entrada. Uma versão bem simples, onde tentamos encontrar um divisor de p indo de 2 até \sqrt{p}, pode ser dada por:

ehPrimo p | p < 2 = False
          | otherwise = ehPrimoAux p 2 
                      
ehPrimoAux p q | q*q > p        = True
               | p `mod` q == 0 = False 
               | otherwise      = ehPrimoAux p (q+1)

Na biblioteca padrão do Haskell, existe a função filter, que recebe um array e uma função que recebe um elemento desse array e retorna verdadeiro ou falso. Podemos passar nossa função ehPrimo para gerar uma lista com os elementos primos de 1 a 50 por exemplo.

> :type filter
filter :: (a -> Bool) -> [a] -> [a]
> filter ehPrimo [1..50]
[2,3,5,7,11,13,17,19,23,29,31,37,41,43,47]

Para facilitar, podemos criar uma função que recebe um array e retorna a lista de primos:

listaPrimos xs = filter ehPrimo xs

Utilizando currying, é possível definir a função acima sem o parâmetro

listaPrimos = filter ehPrimo

No primeiro caso, estamos definindo uma função que chama filter ehPrimo xs. No segundo, estamos apenas associando o resultado da avaliação parcial da função filter com seu primeiro argumento, no caso o ehPrimo.

Currying sobre o segundo parâmetro

Por padrão, a substituição dos argumentos é feita a partir do primeiro argumento, ou seja, da esquerda para a direita. Entretanto, para funções com dois parâmetros (binárias), podemos usar a notação infixa para fazer a substituição do segundo argumento.

A notação infixa consiste em chamar a função entre acentos graves (backtick). Por exemplo, a função soma poderia ser chamada assim:

> 1 `soma` 2
3

Essa sintaxe é interessante para fins de legibilidade de código. Ao chamar a função na forma infixa passando somente o segundo argumento, a substituição é feita somente para ele. Como um exemplo, vamos definir a função potência que recebe dois inteiros a e b e eleva a a b (a função soma não tem graça porque ela é comutativa).

> let potencia a b = a ^ b

Agora, podemos fazer currying e aplicar a função a um argumento apenas

-- Faz a substituição do parâmetro b
> let quadrado = (`potencia` 2)
> quadrado 9
81

Currying em C++

A STL do C++ provê algumas ferramentas para trabalhar com o paradigma funcional, através da biblioteca functional. Através dela podemos reproduzir o nosso primeiro exemplo em Haskell:

#include <functional>
#include <iostream>
using namespace std;

struct adder : public binary_function<int, int, int> {
    int operator() (int a, int b) const {
        return a + b;
    }
};
// currying sobre o primeiro argumento da funcao 
int main (){

    adder soma;
    binder1st<adder> somaTres = bind1st(soma, 3);  
    cout << somaTres(7) << endl; // 10
    cout << somaTres(-1) << endl; // 2
}

Nesse caso trabalhamos com um functor, que é um objeto representando uma função. A classe adder deve derivar da classe binary_function. A função bind1st faz um bind do primeiro argumento da função, e retorna um objeto do tipo binder1st, que deriva da classe unary_function e que recebe um argumento.

É possível fazer o bind do segundo argumento através da função bind2nd.


Haskell – Tipos de Dados Algébricos

Setembro 25, 2011

Haskell usa uma sintaxe elegante para representar um tipo de estrutura de dados, análoga a diversas outras estruturas presentes em linguagens como C e C++. Essas estruturas são chamadas tipos de dados algébricos. Minha principal referência para esse post foi o capítulo 3, do livro Real World Haskell.

Este é o segundo post sobre meus estudos em Haskell. O primeiro se encontra aqui.

Struct

Considere a seguinte struct em C++ representando um livro através de seu código, título e autores:

#include <string>
#include <vector>
struct {
    int isbn;
    std::string title;
    std::vector < std::string > authors;
};

Em haskell podemos escrever:

data BookInfo = Book Int String [String]
                deriving (Show)

Para definir essa estrutura, deve-se iniciar com a keyword data e dar um nome para esse novo tipo, no caso BookInfo (chamado de type constructor). Além disso, definimos um construtor desse tipo, que no exemplo tem o nome Book (chamado de value constructor ou data constructor), mas poderia ter o mesmo nome do tipo. Tanto o nome do construtor quanto o nome do tipo devem ser iniciados por letra maiúscula.

A linha deriving (Show) está sendo usada nesse e nos outros exemplos deste post para que uma instância desse estrutura possa ser impressa no terminal. Podemos instanciar a estrutura BookInfo da seguinte forma:

myBook = Book 9780135072455 "Algebra of Programming"
         ["Richard Bird", "Oege de Moor"]

Note que os parâmetros de entrada devem estar na mesma ordem em que eles são definidos na estrutura.

Note também que os “membros” dessa estrutura são anônimos, sendo representados apenas por seu tipo. Para contornar a possível perda de legibilidade, é possível usar sinônimos de tipos, através da keyword type (parecido com o typedef de C++).

type ISBN = Int
type Title = String
type Author = String

data BookInfo = Book ISBN Title [Author]
                deriving (Show)

Enum

Outra estrutura em C++ que pode ser representada sucintamente em Haskell é o enum. Considere um enum representando as cores RGB:

enum Color { Red, Green, Blue };

Um tipo de dado algébrico admite mais de um construtor, desde que tenham nomes diferentes. Além disso, os construtores não precisam ter parâmetros, o que permite escrever o enum acima da seguinte forma:

data Color = Red | Green | Blue
             deriving (Show)

Union

Ainda aproveitando que podemos ter vários construtores com diferentes parâmetros, podemos também representar uma union. Considere o seguinte exemplo em C++ de uma union podendo representar duas formas geométricas – um círculo ou um polígono:

struct Point {
    double x, y;
};
struct Circle {
    Point center;
    double radius;
};
struct Polygon {
    int nvertices;
    Point *vertices;
};
union Shape {
    Circle circle;
    Polygon polygon;
};

Em Haskell, é possível fazer o seguinte:

type Point = (Double, Double)
type Center = Point
type Radius = Double

data Shape = Circle Center Radius
           | Polygon [ Point ]
             deriving (Show)

Curiosidade: Na especificação atual do C++ (2003), membros de uma union não podem ter construtores não-triviais. Parece porém, que na nova especificação (C++11) será possível. Portanto, o exemplo abaixo não compila com gcc-4.2.1:

#include <utility>
#include <vector>

typedef std::pair <double, double> Point;
typedef std::vector<Point> Polygon;

struct Circle {
    Point center;
    double radius;
};

union Shape {
    Circle circle;
    Polygon polygon;
};

Getters e Setters

Quando trabalhamos com classes, eventualmente queremos proteger o acesso de leitura/escrita a seus membros usando getters e setters, como por exemplo:

#include <string>
class Person {
    int id;
    std::string name;

public:

    int get_id(){ return id; }
    void set_id(int _id){ id = _id; };
    std::string get_name(){ return name; }
    void set_name(std::string &_name) { name = _name; }
};

Há um mecanismo parecido em Haskell, implementado pelo chamado record syntax. Num exemplo análogo ao anterior:

type ID = Int
type Name = String

data Person = Person {
      personID      :: ID,
      personName    :: Name
      } deriving (Show)

Essa sintaxe mais verbosa permite instanciar a estrutura atribuindo os parâmetros aos membros explicitamente como:

you = Person {personName = "YourName", personID = 2}

O getter de um membro é dado por:

personName you

Enquanto o setter seria:

you {personID 3}

Tipos Parametrizados

Tipos de dados algébricos representam ainda outra funcionalidade presente em C++: templates. Tomemos como exemplo uma classe representando um ponteiro genérico para qualquer tipo:

template <typename T>
struct GenericPointer {
    T* ptr;
};

Ao instanciar essa classe podemos definir os tipos que queremos representar:

int i = 2;
GenericPointer<int> p1;
p1.ptr = &i;
GenericPointer<double> p2;
p2.ptr = NULL;

Aqui temos ponteiros para os tipos int e double, sendo que o NULL é um valor genérico para dizer que o ponteiro não aponta pra lugar algum.

Uma analogia (forçada) em Haskell seria a seguinte:

data GenericPointer a = Ptr a |
                        Null
                        deriving(Show)

Embora em Haskell não tenhamos o valor NULL, podemos simulá-lo através de um construtor como no exemplo acima. Esse encapsulamento é particularmente interessante para prover uma condição de exceção. Por exemplo, considere uma função que retorna o primeiro elemento de uma lista:

first_elem (x:xs) = x

O que retornar caso a lista seja vazia? Uma solução seria usar o encapsulamento:

first_elem (x:xs) = Ptr x
first_elem _      = Null

No exemplo, o caractere '_' funciona como um coringa, ou seja, qualquer parâmetro casa com ele. Como foi posto depois do outro padrão, ele se comporta como um “else” e só será chamado quando a lista não tiver elementos.

Tipos Recursivos

Outra possibilidade de um tipo de dado algébrico é usar tipos recursivos, isto é, membros da estrutura são do tipo dessa estrutura. Em C e C++ é possível fazer isso utilizando ponteiros. Por exemplo, podemos implementar uma lista ligada de elementos genéricos da seguinte maneira:

template <typename T>
struct LinkedList {
    LinkedList<T>* next;
};

Em Haskell temos algo como:

data LinkedList a = Node a (LinkedList a)
            | Null 
            deriving(Show)

Para representar uma lista ligada contendo 3, 1 e 2, nessa ordem, usamos a seguinte sintaxe:

Node 3 (Node 1 (Node 2 Nil))

Um árvore binária também é bem simples de se representar:

data Tree a = Node a (Tree a) (Tree a)
            | Null
              deriving (Show)

Conclusão

Haskell tem uma sintaxe sucinta para descrever tipos de dados algébricos, que por sua vez são bastante versáteis, podendo representar tanto estruturas quanto uniões, enumerações, classes, etc.

Ainda sei muito pouco sobre a linguagem, mas por enquanto estou achando-a bastante interessante.


Shader OSL Multithread

Julho 24, 2011

Essas últimas semanas eu estava tentando implementar suporte a múltiplas threads do shader OSL. O BRL-CAD usa pthreads para paralelizar o trabalho.

A primeira vez que eu pus o shader OSL pra rodar, tive problemas, provavelmente devido à leitura e escritas simultâneas a dados na memória. Eu tinha adiado esse problema e estava trabalhando apenas com uma thread.

O que a minha função osl_render faz é preencher uma estrutura chamada RenderInfo com dados como ponto de interseção, a normal, um ponteiro para o shader dessa superfície, etc. Preenchida essa estrutura, ela faz uma chamada ao método QueryColor da classe OSLRender. O objeto dessa classe é global, o que significa que está em uma posição de memória compartilhada pelas threads. Como a chamada do método QueryColor modifica o objeto, temos o problema do acesso simultâneo a esse método.

Primeira tentativa

A primeira coisa que eu tentei, foi proteger a chamada dessa função com semáforos. O BRL-CAD tem uma biblioteca interna de utilitários, entre eles a implementação de um semáforo.

Basicamente, o que eu fiz foi o seguinte:

OSLRender oslr;
...
bu_semaphore_acquire(BU_SEM_SYSCALL);
oslr.QueryColor();
bu_semaphore_release(BU_SEM_SYSCALL);

BU_SEM_SYSCALL na verdade é uma macro representando um inteiro. Quando uma thread chama a função bu semaphore acquire(BU_SEM_SYSCALL), ela bloqueia o acesso a outras threads para o mesmo inteiro.

Com isso o código passou a funcionar para mais de uma thread! Porém, fui conservador demais ao bloquear a chamada do QueryColor. Creio que 90% do trabalho da função osl_render seja feito lá. Logo, colocando um semáforo nessa função, a paralelização vai por água abaixo. De fato, fiz testes com 1 e 2 processadores e os tempos para renderizar uma imagem de teste foram exatamente os mesmos.

Segunda tentativa

Fui estudar um pouco mais o código do OSL e descobri que ele tem suporte a aplicativos multithread. Basicamente, cada thread deve manter uma estrutura local chamada ThreadInfo. Ao chamar o sistema de shaders do OSL, passamos essa informação, que será modificada.

Gambiarra!

Para manter uma ThreadInfo para cada thread, precisamos identificá-las de alguma forma. O problema é que não temos acsso a essa informação diretamente. O único jeito que encontrei sem ter que mudar a estrutura interna do raytracer do BRL-CAD foi a seguinte: toda vez que a função osl_render é chamada, uma estrutura chamada struct application é passada. Essa estrutura tem um campo chamado resource, que é um ponteiro para um pedaço de memória exclusivo para cada thread.

Esse pedaço é alocado na criação das threads e não muda ao longo de suas vidas. Assim, usar o endereço desse pedaço de memória é um jeito (feio) de identificá-las. Outro problema é que, embora haja uma função de inicialização do shader (chamado osl_setup), as threads são criadas após isso, o que me obriga a inicializá-las na função osl_render mesmo.

A solução na qual cheguei foi a seguinte:

/* Every time a thread reaches osl_render for the 
   first time, we save the address of their own
   buffers, which is an ugly way to identify them */
std::vector<struct resource *> visited_addrs;
/* Holds information about the context necessary to
   correctly execute a shader */
std::vector<void *> thread_infos;

int osl_render(struct application *ap, ...){
    ...
    bu_semaphore_acquire(BU_SEM_SYSCALL);
    
    /* Check if it is the first time this thread is
       calling this function */
    bool visited = false;
    for(size_t i = 0; i < visited_addrs.size(); i++){
        if(ap->a_resource == visited_addrs[i]){
            visited = true;
            thread_info = thread_infos[i];
            break;
        }
    } 
    if(!visited){
        visited_addrs.push_back(ap->a_resource);
        /* Get thread specific information from 
           OSLRender system */
        thread_info = oslr->CreateThreadInfo();
        thread_infos.push_back(thread_info);
    }
    bu_semaphore_release(BU_SEM_SYSCALL);
    ...
}

Classe thread-safe

Com as ThreadInfo's, podemos nos despreocupar com o acesso concorrente ao objeto da classe OSLRender. Para ter uma garantia de que não será feita escrita nesse objeto, podemos declarar os métodos como const, como no exemplo abaixo:

struct A {
    int x;
    // OK
    void set_x(int _x) { x = _x; };
    // ERRO, tentando modificar membro do objeto 
    // num método const
    void const_set_x(int _x) const { x = _x; };
};

Testes

Com essa nova implementação, os tempos melhoraram bastante com mais processadores. A tabela abaixo mostra os tempos de execução para renderizar uma cena de exemplo, usando até 4 processadores.

Gráfico 1: Tempo x no. de processadores

Os ganhos foram praticamente lineares, sendo que para 4 processadores a medida de paralelização foi de 1.25 (se fosse perfeitamente linear seria 1.0).

Próximos passsos

Falei com meu mentor sobre o desenvolvimento daquele modo de visualização incremental, mas parece que ele não gostou muito da ideia. Isso exigiria uma mudança bastante grande no sistema, pois o renderizador do BRL-CAD foi desenvolvido para ray-tracing e não path-tracing.

Por enquanto estou estudando maneiras de adaptar o código OSL para suportar ray-tracing, mas não sei se isso é viável. A propósito, eu ainda confundo bastante os termos ray-tracing, path-tracing, photon mapping e todos esses algoritmos de iluminação e pretendo em breve escrever um post cobrindo esses tópicos bem por cima.


Usando bibliotecas C++ em código C

Junho 12, 2011

O BRL-CAD é escrito majoritariamente em C, enquanto o OSL é escrito em C++. Para que o código C possa usar código C++ é preciso usar uns truques. No meu caso, eu queria disponibilizar uma classe C++ para ser acessada pelo C.

Para indicar que um código C irá usar uma função sua, o código C++ deve usar a diretiva extern "C". Por exemplo, se temos a função void foo(int x), ela ficará assim:

extern "C" void foo(int x);

ou alternativamente

extern "C" {
    void foo(int x);
}

Entretanto, o código C não entende essa diretiva. Por isso, o cabeçalho que ele deverá ver é

void foo(int x);

Um truque para o código C e C++ compartilharem o mesmo header, é usar a macro __cpluscplus que só é setada se o código for C++. Podemos usar o seguinte cabeçalho para ambos os códigos:

#ifdef __cplusplus
extern "C" {
#endif
    void foo(int x);
#ifdef __cplusplus
}
#endif

Exportando uma classe

Agora suponha que queremos exportar uma classe, como por exemplo a classe A, definida como:

class A {
    int x;
 public:
    A(int _x) : x(_x){};
    int getx() { return x; };
};

Precisamos definir funções estilo C para encapsular os métodos da class A, inclusive o construtor e o destrutor.

/* wrapper para construtor */
A* init(int x){
    A* objA = new A(x);
    return objA;
}
/* wrapper destrutor */
void destroy(A** objA){
    delete *objA;
    *objA = NULL;
}
/* wrapper para A::getx() */
int getx(A* objA){
    return objA->getx();
}

Para que o código C mantenha uma referência para o objeto da classe A, podemos usar um ponteiro opaco, já que a estrutura da classe A é desconhecida para o código C. Além disso, para poder compartilhar as funções acima com o código C++, devemos criar um apelido para nosso ponteiro opaco, de A.

typedef void A;

Agora A* representa um ponteiro opaco que pode apontar para um objeto da classe A, independente de sua estrutura. Aliás, ao invés de usarmos void, recomenda-se usar struct A para que haja uma diferenciação dos tipos dos pointeiros opacos [2]. Isso é interessante no caso em que há duas classes A e B, e erroneamente fazemos um ponteiro da classe B apontar para um objeto da classe A, como no código abaixo:

typedef void A;
typedef void B;
A* objA = init(10);
B* objB = objA; 

Se ao invés de void usarmos structs,

typedef struct A A;
typedef struct B B;
A* objA = init(10);
B* objB = objA; /* Gera um warning */

Temos uma garantia adicional contra erros não-intencionais.

Exemplo completo

lib.h:

 
#include <stdio.h>

#ifdef __cplusplus
class A {
    int x;
 public:
    A(int _x);
    int getx();
};

#else
typedef struct A A;
#endif

#ifdef __cplusplus
extern "C" {
#endif

    A* init(int x);
    void destroy(A** objA);
    int getx(A* objA);

#ifdef __cplusplus
}
#endif

lib.cpp:

#include "lib.h"
A::A(int _x) : x(_x) {}
int A::getx(){ return x; }

A* init(int x){
    A* objA = new A(x);
    return objA;
}
void destroy(A** objA){
    delete *objA;
    *objA = NULL;
}
int getx(A* objA){
    return objA->getx();
}

main.c:

#include <stdio.h>
#include "lib.h"

int main (){

    A *objA = init(3); 
    printf("%d\n", getx(objA));
    destroy(&objA);
    return 0;
}

Compilando

Podemos usar o cmake para simplificar o processo de compilação:

add_library(lib lib.cpp)
add_executable(main main.c)
target_link_libraries(main lib)

Referências

[1] http://www.parashift.com/c++-faq-lite/mixing-c-and-cpp.html
[2] Stack overflow: A question about exporting a C++ class to C