Geradores de números pseudo-aleatórios

Geradores de números pseudo-aleatórios GNPA são uma caixa preta, em geral usamos as indicações dos estudantes mais antigos e dos professores. Provavelmente, graças a esse mecanismo, os GNPA mais usados são o ran2 ou o rand, e a implementação destes é encontrada no clássico Numerical Recipes.

O cerne das simulações de Monte-Carlo são justamente os GNPA e dificilmente nossos programas passarão fazendo mais tempo algo que não seja gerar números aleatórios. Por isso é importante conhecer e se atualizar quanto este assunto.

Embora seja importante entender um pouco da minúcias dos GNPA eu não vou nem tocar nesse assunto, para mim eles permaneceram uma caixa preta por muito tempo. Bom..., as atualidades ? Ah, essas são bem mais fáceis. 

O GNPA padrão nos dias de hoje o Mersenne-Twister, introduzido em 1997 http://en.wikipedia.org/wiki/Mersenne_twister (alguem tem que escrever o artigo na wikipedia portuguesa, nem que seja tradução). Este é provavelmente o algoritmo mais popular, ele possui um período de 2^19937 − 1, passou em vários testes e é superior aos ran2,3,4 .

Outro ponto de fundamental importância é a sua performance.

Como implementar o algoritmo:

Eu sugiro fortemente a utilização da biblioteca GSL ( Gnu Scientific Library ) http://www.gnu.org/software/gsl/ e o manual http://www.gnu.org/software/gsl/manual , essa biblioteca já esta instalada no abax e você pode encontra-la nos repositórias da sua distribuição Linux.

E muito fácil usar, segue um exemplo

#include < stdio.h >
#include < gsl/gsl_rng.h >
int main ()
  {
  int i;
  unsigned long a;
/////////////////////////////////////////////////////
  const gsl_rng *r;
  r = gsl_rng_alloc(gsl_rng_mt19937);
////////////////////////////////////////////////////
  gsl_rng_set(r, 123); // semente 123
   for (i = 0; i < 100000000; i++)
    {
    printf ("%u\n",  gsl_rng_get(r));
    }
  printf ("\n");
  printf ("%lu %lu\n", gsl_rng_max(r), gsl_rng_min(r) );
  return 0;
}

Compilar com: gcc rng.c -lgsl -lgslcblas -lm

Comentários:

"const gsl_rng *r;
r = gsl_rng_alloc(gsl_rng_mt19937); "
- Define o nome do gerador, no caso é "r",  e o GNPA usado ( mt19937 = Mersenne Twisted ).

"gsl_rng_set (r, 123);" - Semente (123) do gerador r.

"gsl_rng_get (r);" - Retorna um número inteiro (unsigned long) pseudo-aleatório

"gsl_rng_max(r) e gsl_rng_min(r)" - Valor máximo e mínimo do gerador.


Outras dicas:

"gsl_rng_uniform(r);" - Retorna um real (double) pseudo-aleatório o intervalo [0,1)
"gsl_rng_uniform_int(r,j);" Retorna um número inteiro (unsigned long) pseudo-aleatório entre [0, j-1]

Existe uma lista de extensa de implementações de GNPA nesta biblioteca, inclusive o ran2, basta substituir o "mt19937" pelo algoritmo desejado.

Veja a lista de performance ( http://www.gnu.org/software/gsl/manual/html_node/Random-Number-Generator-Performance.html ):

     1754 k ints/sec,    870 k doubles/sec, taus
     1613 k ints/sec,    855 k doubles/sec, gfsr4
     1370 k ints/sec,    769 k doubles/sec, mt19937
      565 k ints/sec,    571 k doubles/sec, ranlxs0
      400 k ints/sec,    405 k doubles/sec, ranlxs1
      490 k ints/sec,    389 k doubles/sec, mrg
      407 k ints/sec,    297 k doubles/sec, ranlux
      243 k ints/sec,    254 k doubles/sec, ranlxd1
      251 k ints/sec,    253 k doubles/sec, ranlxs2
      238 k ints/sec,    215 k doubles/sec, cmrg
      247 k ints/sec,    198 k doubles/sec, ranlux389
      141 k ints/sec,    140 k doubles/sec, ranlxd2
    
     1852 k ints/sec,    935 k doubles/sec, ran3
      813 k ints/sec,    575 k doubles/sec, ran0
      787 k ints/sec,    476 k doubles/sec, ran1
      379 k ints/sec,    292 k doubles/sec, ran2


Uma pequena explicação sobre as características de cada algoritmo é encontrada no manual http://www.gnu.org/software/gsl/manual/html_node/Random-number-generator-algorithms.html