💒 🖍️ 👨🏼‍🤝‍👨🏻 SHISHUA: o gerador de números pseudo-aleatórios mais rápido do mundo 🏞️ 🖕🏾 🍖

Seis meses atrás, eu queria criar o melhor gerador de números pseudo-aleatórios (PRNG) com alguma arquitetura incomum. Eu pensei que o começo seria fácil e, à medida que você trabalha, a tarefa se torna lentamente mais complexa. E eu pensei que se eu pudesse aprender tudo rápido o suficiente para lidar com as mais difíceis.

Para minha surpresa, a complexidade não aumentou linearmente. Teste de qui-quadrado byte provou ser muito difícil! Mais tarde, foi igualmente difícil passar em testes obstinados. Eu publicou os resultados atuais para entender o que outras dificuldades esperam por mim. No entanto , o teste do PractRand falhou naquele momento .

Então foi muito difícil passar no teste BigCrush .

Foi muito difícil transferir 32 terabytes de dados ao passar o PractRand. A velocidade se tornou um problema. Não foi suficiente criar um design que gerasse dez megabytes por segundo, porque a aprovação do PractRand levaria um mês. Mas devo admitir que foi muito difícil passar neste teste a uma velocidade de gigabytes por segundo .

Quando você chega a tal altura ... você quer saber se pode chegar à fronteira de Pareto. Você deseja criar o PRNG mais rápido do mundo, que passará nos testes estatísticos mais complexos.

Eu consegui.

Em um artigo anterior, falei sobre o que aprendi para alcançar meu objetivo. E aqui vou lhe dizer como a arquitetura final funciona.

objetivo

Vamos começar com o óbvio: a velocidade depende da plataforma . Eu me concentrei na otimização da arquitetura moderna x86-64 (processadores Intel e AMD).

A métrica cpb clássica é usada para comparar o desempenho : este é o número de ciclos do processador gastos na geração de um byte. Essa métrica é calculada e comparada em todos os trabalhos criptográficos. Uma cpb um pouco menor no mundo de software ou hardware pode garantir a vitória na competição ou o uso em sites em todo o mundo.

Para melhorar o cpb, você pode:

Gere mais bytes com a mesma quantidade de trabalho,
Ou faça menos trabalho para gerar o mesmo número de bytes,
Ou paralelizar o trabalho.

Faremos tudo o que precede.

De acordo com o primeiro ponto, precisamos produzir mais bits a cada iteração.

Eu tinha medo de que eles me dissessem: “Se não der números de 32 bits, esse não é o PRSP”, ou o mesmo com números de 64 bits. Ou: “O PRNG deve ser apenas para a arquitetura x86-64”, como se instruções como POPCNT ou registros como% xmm7 fossem proibidas.

No entanto, o PRNG é de engenharia: os geradores vêm tentando há várias décadas extrair todo o possível dos processadores! Quando o ROL apareceu, eles começaram a confiar nele. Com o advento dos processadores de 64 bits, eles começaram a confiar no% rax. Obviamente, no ARM, esses algoritmos podem funcionar mais lentamente (embora isso ainda deva ser visto), no entanto, os PRNs de 64 bits foram usados ativamente mesmo antes do Android começar a exigir suporte para 64 bits em 2019!

Ou seja, esta área está se desenvolvendo junto com o hardware. E hoje, os processadores Intel e AMD, devido ao AVX2, já suportam operações de 256 bits. O RC4 produziu 1 byte, o drand48 pode produzir 4 bytes por vez, pcg64 - 8 bytes, e agora podemos gerar imediatamente 32 bytes.

8 bytes podem ser um número de 64 bits, e a maioria das linguagens de programação possui tipos internos para isso. Mas poucos idiomas fornecem tipos para 16 bytes (uma exceção notável é __uint128_t em C). Menos idiomas ainda têm tipo para 32 bytes (exceto interno).

Portanto, podemos dizer adeus ao protótipo usual da função PRNG (exemplo do benchmark Vigny HWD ):

static uint64_t next(void);

Em vez disso, você pode criar um gerador que preencha o buffer (exemplo do meu benchmark ):

void prng_gen(prng_state *s, __uint64_t buf[], __uint64_t size);

Quais são as desvantagens desta solução?

Se o seu gerador produzir 32 bytes de cada vez, você precisará que o consumidor forneça uma matriz com um múltiplo de 32 (idealmente alinhado a 32 bytes). Embora você possa ficar sem ele, preencheremos o buffer. Removeremos os dados não utilizados e os preencheremos novamente conforme necessário. O atraso se torna imprevisível: algumas chamadas apenas lerão o buffer. Mas, em média, tudo será o mesmo.

Agora geramos mais bytes, fazendo a mesma quantidade de trabalho. Como nós a paralelizamos?

Paralelismo

Os processadores oferecem um incrível conjunto de ferramentas de paralelização em todos os níveis.

Em primeiro lugar, estas são as instruções SIMD (instrução única, dados múltiplos). Por exemplo, o AVX2 executa simultaneamente quatro adições de 64 bits ou oito de 32 bits, etc.

É usado em criptografia há cerca de quinze anos. A simultaneidade proporcionou o incrível desempenho do ChaCha20 . É usado pelas primitivas mais importantes que não usam AESNI. Por exemplo, NORX e Gimli são projetados com o paralelismo em mente.

Recentemente, o interesse neste tópico também aumentou na comunidade PRNG não criptográfica. Em particular, as primitivas existentes que não foram projetadas para SIMD podem ser a base para a criação de PRNs muito rápidos.

Quando Sebastiano Vigna promoveu sua arquitetura xoshiro256 ++ na biblioteca padrão Julia, ele descobriu que os resultados de oito instâncias PRNG competitivas e diferentemente inicializadas podem ser concatenadas muito rapidamente se cada operação for executada simultaneamente em todas as PRNRs.

SIMD é apenas um dos níveis de paralelização no processador. Eu recomendo a leitura do artigo anterior sobre este tópico para ter uma idéia melhor, mas darei algumas explicações.

Os pipelines do processador permitem que várias instruções sejam processadas em diferentes estágios. Se a ordem de execução for bem organizada para reduzir as dependências entre os estágios, o processamento das instruções poderá ser acelerado.

A execução superescalar permite processar simultaneamente as partes computacionais das instruções. Mas para isso, eles não devem ter dependências de leitura e gravação. Você pode adaptar a arquitetura para reduzir o risco de tempo de inatividade gravando muito antes da leitura.

A execução extraordinária permite que o processador execute instruções não na ordem da sequência, mas como estão prontas, mesmo que as instruções anteriores ainda não estejam prontas. Mas, para isso, não deve haver dependência de leitura e gravação.

E agora passamos à implementação!

Arquitetura

Considere um esquema chamado semi-SHISHUA. A origem do nome se tornará gradualmente aparente enquanto você lê.

O esquema é assim:

Considere sua linha por linha.

typedef struct prng_state {
  __m256i state[2];
  __m256i output;
  __m256i counter;
} prng_state;

O estado é dividido em duas partes, que são colocadas no registro AVX2 (256 bits). Para aumentar a velocidade, mantemos o resultado próximo ao próprio estado, mas ele não faz parte do estado.

Também temos um contador de 64 bits. Para simplificar o cálculo, também é um registro AVX2. O fato é que o AVX2 possui um pequeno recurso: os registros comuns (% rax e similares) não podem ser transferidos diretamente para o SIMD via MOV, eles devem passar pela RAM (geralmente pela pilha), o que aumenta o atraso e custa duas instruções do processador (MOV na pilha, VMOV da pilha).

Agora vamos olhar para a geração. Vamos começar carregando, depois percorrer o buffer e preenchê-lo com 32 bytes a cada iteração.

inline void prng_gen(prng_state *s, __uint64_t buf[], __uint64_t size) {
  __m256i s0 = s->state[0], counter = s->counter,
          s1 = s->state[1],       o = s->output;
  for (__uint64_t i = 0; i < size; i += 4) {
    _mm256_storeu_si256((__m256i*)&buf[i], o);
    // …
  }
  s->state[0] = s0; s->counter = counter;
  s->state[1] = s1; s->output  = o;
}

Como a função está em linha, o preenchimento imediato do buffer na inicialização permite que o processador execute imediatamente as instruções, dependendo disso, por meio de um mecanismo de execução extraordinário.

Dentro do loop, executamos rapidamente três operações de estado:

SHI ft
SHU ffle
A dd

Daí o nome SHISHUA!

Primeira marcha

u0 = _mm256_srli_epi64(s0, 1);              u1 = _mm256_srli_epi64(s1, 3);

Infelizmente, o AVX2 não suporta rotações. Mas eu quero misturar os bits de uma posição em um número de 64 bits com os bits de outra posição! E a mudança é a melhor maneira de perceber isso. Mudaremos para um número ímpar, para que cada bit visite todas as posições de 64 bits, e não metade delas.

Durante o turno, os bits são perdidos, o que leva à remoção de informações do nosso estado. Isso é ruim, você precisa minimizar as perdas. Os menores números ímpares são 1 e 3, usaremos diferentes valores de deslocamento para aumentar a discrepância entre as duas partes. Isso ajudará a reduzir a semelhança de sua auto-correlação.

Mudaremos para a direita, porque os bits mais à direita têm a menor difusão durante a adição: por exemplo, o bit menos significativo em A + B é apenas o XOR dos bits mais baixos A e B.

Agitação

t0 = _mm256_permutevar8x32_epi32(s0, shu0); t1 = _mm256_permutevar8x32_epi32(s1, shu1);

Usaremos a mistura de 32 bits, porque fornece uma granularidade diferente em comparação às operações de 64 bits que usamos em qualquer lugar (o alinhamento de 64 bits é violado). Também pode ser uma operação de faixa cruzada: outros shuffles podem mover bits dentro dos 128 bits esquerdos se começarem à esquerda, ou dentro dos 128 bits direitos se começarem à direita.

Constantes de mistura:

__m256i shu0 = _mm256_set_epi32(4, 3, 2, 1, 0, 7, 6, 5),
        shu1 = _mm256_set_epi32(2, 1, 0, 7, 6, 5, 4, 3);

Para que a mistura realmente melhore o resultado, moveremos as partes fracas (baixa dispersão) de 32 bits das adições de 64 bits para posições fortes, de modo que as próximas adições as enriquecem.

A parte inferior de 32 bits do bloco de 64 bits nunca se move para o mesmo bloco de 64 bits da parte de alta ordem. Assim, ambas as partes não permanecem no mesmo pedaço, o que melhora a mistura.

No final, cada parte de 32 bits passa por todas as posições em um círculo: de A a B, de B a C, ... de H a A.

Você deve ter notado que a mistura mais simples que leva em conta todos esses requisitos são dois de 256 bits rotatividade (rotações de 96 bits e 160 bits para a direita, respectivamente).

Adição

Vamos adicionar pedaços de 64 bits de duas variáveis temporárias - shift e mixagem.

s0 = _mm256_add_epi64(t0, u0);              s1 = _mm256_add_epi64(t1, u1);

A adição é a principal fonte de dispersão: nesta operação, os bits são combinados em combinações irredutíveis de expressões XOR e AND distribuídas em posições de 64 bits.

Armazenar o resultado da adição dentro de um estado preserva permanentemente essa dispersão.

Função de saída

De onde obtemos a saída?

É simples: a estrutura que criamos nos permite gerar duas partes independentes do estado s0 e s1, que não se afetam de forma alguma. Aplique o XOR a eles e obtenha um resultado completamente aleatório.

Para fortalecer a independência entre os dados aos quais aplicamos o XOR, obtemos um resultado parcial: a parte deslocada de um estado e a parte mista de outro.

o = _mm256_xor_si256(u0, t1);

Isso é semelhante à redução das dependências de leitura e gravação entre instruções em um processador superescalar, como se u0 e t1 estivessem prontos para ler para s0 e s1.

Agora discuta o balcão. Nós a processamos no início do ciclo. Primeiro, altere o estado e aumente o valor do contador:

s1 = _mm256_add_epi64(s1, counter);
counter = _mm256_add_epi64(counter, increment);

Por que primeiro alteramos o estado e atualizamos o contador? s1 se torna disponível anteriormente, isso reduz a probabilidade de que as instruções subsequentes que a leem parem no pipeline do processador. Além disso, essa sequência ajuda a evitar a dependência direta do contador de leitura e gravação.

Aplicamos o contador a s1, e não a s0, porque ambos afetam a saída de qualquer maneira; no entanto, s1 perde mais bits devido à mudança, de modo que ajuda a "se levantar" após a mudança.

O contador pode não gravar o teste do PractRand. Seu único objetivo é estabelecer um limite inferior de 2 ⁶⁹bytes = 512 exbibytes para o período PRNG: começamos a repetir o ciclo somente após um milênio de trabalho a uma velocidade de 10 gibytes por segundo. É improvável que seja muito lento para uso prático nos próximos séculos.

Incremento:

__m256i increment = _mm256_set_epi64x(1, 3, 5, 7);

Os números ímpares são escolhidos como incrementos, porque apenas os números básicos do coprime cobrem todo o ciclo do campo finito GF (2 ⁶⁴ ), e todos os números ímpares são coprime para 2. Em outras palavras, se você incrementar por um número inteiro par no intervalo de 0 a 4, retornando a 0 após 4, ocorre a sequência 0-2-0-2- ..., que nunca levará a 1 ou 3. E o incremento ímpar passa por todos os números inteiros.

Para todos os números de 64 bits no estado, usaremos números ímpares diferentes, que os separarão ainda mais e aumentarão ligeiramente a mistura. Eu escolhi os menores números ímpares para que eles não pareçam mágicos.

É assim que as funções de transição e saída de estado funcionam. Como inicializá-los?

Inicialização

Inicializamos o estado usando os dígitos hexadecimais Φ, o número irracional que é menos aproximado pela fração.

static __uint64_t phi[8] = {
  0x9E3779B97F4A7C15, 0xF39CC0605CEDC834, 0x1082276BF3A27251, 0xF86C6A11D0C18E95,
  0x2767F0B153D27B7F, 0x0347045B5BF1827F, 0x01886F0928403002, 0xC1D64BA40F335E36,
};

Pegue uma semente de 256 bits. Isso geralmente é feito em criptografia e não prejudica o trabalho de PRNGs não criptográficos:

prng_state prng_init(SEEDTYPE seed[4]) {
  prng_state s;
  // …
  return s;
}

Não queremos redefinir toda a parte do estado (s0 ou s1) com esse número inicial, precisamos apenas afetar a metade. Dessa forma, evitaremos o uso de números iniciais atenuantes, que acidental ou intencionalmente originam um estado inicial fraco conhecido.

Como não alteramos metade de cada estado, mantemos o controle sobre 128 bits de status. Essa entropia é suficiente para iniciar e manter uma posição forte.

s.state[0] = _mm256_set_epi64x(phi[3], phi[2] ^ seed[1], phi[1], phi[0] ^ seed[0]);
s.state[1] = _mm256_set_epi64x(phi[7], phi[6] ^ seed[3], phi[5], phi[4] ^ seed[2]);

Em seguida, repetimos ( ROUNDS) várias vezes a seguinte sequência:

Execute as etapas ( STEPS) das iterações SHISHUA.
Atribuímos uma parte do estado a outro estado e a outra parte à saída.

for (char i = 0; i < ROUNDS; i++) {
  prng_gen(&s, buf, 4 * STEPS);
  s.state[0] = s.state[1];
  s.state[1] = s.output;
}

Atribuir um resultado de saída aumenta a dispersão do estado. Durante a inicialização, o trabalho adicional e a correlação de estados não importam, porque essa série de operações é executada uma vez. Estamos interessados apenas na dispersão durante a inicialização.

Depois de avaliar o efeito na correlação dos valores iniciais, escolhi STEPS2 para o valor e ROUNDS10 para 10. Calculei a correlação calculando as anomalias “incomuns” e “suspeitas” decorrentes da ferramenta de controle de qualidade PRNG no PractRand.

atuação

É difícil medir a velocidade por vários motivos:

A medição do relógio pode não ser precisa o suficiente.
, , - , -, .
, . .
: , .

Eu uso a instrução do processador RDTSC, que calcula o número de ciclos.

Para que qualquer um possa reproduzir meus resultados, eu uso uma máquina virtual baseada em nuvem. Isso não altera o nível dos resultados de benchmark em comparação com os testes locais. Além disso, você não precisa comprar o mesmo computador que o meu. Por fim, há muitas situações em que o PRNG é lançado nas nuvens.

Eu escolhi o Google Cloud Platform N2 (processador Intel) e N2D (processador AMD). A vantagem do GCP é que eles oferecem servidores com processadores de ambos os fabricantes. Neste artigo, focaremos na Intel, mas para a AMD os resultados serão na mesma ordem.

Para aprofundar o assunto, primeiro vamos nos livrar do antigo gerador criptográfico RC4. Incapaz de paralelizar o trabalho, consegui7,5 cpb (ciclos por byte gerado).

Agora vamos executar um MCG muito popular e rápido: o Lehmer128 PRNG mais simples , que passa no teste BigCrush, mostrou 0,5 cpb . Uau, ótimo!

Em seguida, executaremos o desenvolvimento mais recente, usado para tabelas de hash rápidas - wyrand . 0,41 cpb , um pouco melhor!

Alguns PRSPs não passam no teste PractRand de 32 terabytes, mas funcionam muito rapidamente. O Xoshiro256 + atingiu apenas 512 mebibytes, mas mostrou uma velocidade muito alta: 0,34 cpb .

Outro desenvolvimento recente do RomuTrio . Ela afirma ser o PRNG mais rápido do mundo - 0,31 cpb.

Ok, é o suficiente. O que mostrou semi-SHISHUA?

0,14 cpb . Duas vezes mais rápido que o RomuTrio.

Legal. Agora teste o gerador criptográfico ChaCha8. Ele atingiu ... 0,12 cpb .

Oh.

SIMD é verdadeira magia!

Para a comunidade criptográfica, isso não foi uma surpresa especial . ChaCha8 é extremamente fácil de paralelizar. Este é apenas um contador bem-hash em um estado difuso.

E lembre-se de como a equipe de idiomas da Julia tentou combinar várias instâncias da arquitetura de Vigny para criar um PRNG rápido baseado em SIMD? Vejamos o resultado deles usando esta técnica ( 8 pedaços de Xoshiro256 + ). 0,09 cpb !

Tecnicamente, meu laptop pode afetar os resultados. Não sei por que o desenvolvimento da equipe de Julia é mais rápido que o ChaCha8 no GCP, mas mais lento quando testado localmente. Na minha máquina, o semi-SHISHUA roda mais rápido que o desenvolvimento da equipe Julia, mas mais lento que o ChaCha8.

É necessário derrotar todos os concorrentes.

Você provavelmente já está se perguntando por que chamamos a versão anterior do gerador semi-SHISHUA? Porque acabou sendo fácil dobrar a velocidade se você executar duas cópias do semi-SHISHUA.

Semelhante à idéia do comando Julia, inicializamos separadamente dois PRNGs (quatro blocos de um estado de 256 bits), fornecendo alternadamente a saída de seu trabalho.

Mas se criarmos mais estados, podemos produzir ainda mais dados, combinando quatro estados em pares:

o0 = _mm256_xor_si256(u0, t1);
o1 = _mm256_xor_si256(u2, t3);
o2 = _mm256_xor_si256(s0, s3);
o3 = _mm256_xor_si256(s2, s1);

Então nós adquirimos o SHISHUA, que mostrou uma velocidade de 0,06 cpb .

Isso é duas vezes mais rápido que o concorrente mais rápido do mundo, que passou no teste PractRand de 32 terabytes. O resultado está no gráfico.

Acredito que o desenvolvimento acabou sendo competitivo. Funciona ainda mais rápido no meu laptop - 0,03 cpb, mas seguirei meus princípios em relação ao benchmark.

Espero que por mais algumas semanas meu gerador fique no pódio dos mais rápidos do mundo (por favor, faça isso).

Qualidade

O gerador passa honestamente no BigCrush e no teste PractRand de 32 terabytes. E tudo graças a quatro fluxos de saída.

As desvantagens da arquitetura incluem sua irreversibilidade . Isso pode ser visto reduzindo-se para um estado de 4 bits com s0 = [a, b]e s1 = [c, d]. Com uma mudança, obtemos [0, a]e [0, d], e com agitação, [b, c]e [d, a]. Novos s0iguais [b, c] + [0, a] = [b⊕(a∧c), a⊕c], mas s1iguais [d, a] + [0, c] = [d⊕(a∧c), a⊕c]. Se a = ¬c, então, a⊕c = 1e a∧c = 0, portanto, s0 = [b, 1]e s1 = [d, 1]. Ou seja, temos duas combinações de ae c, que nos dão o mesmo estado final.

No nosso caso, isso não é um problema, porque o contador de 64 bits também faz parte do estado. Acontece que o ciclo mínimo de 2 ⁷¹bytes (128 bytes por transição de estado), a uma velocidade de 10 gibytes / s. durará sete mil anos. Isso equilibra os estados perdidos.

Além disso, apesar da irreversibilidade, o período médio de transição entre estados é de 2 ^ ((256 + 1) ÷ 2). Isso fornece um ciclo médio de ^2.135 bytes (a uma velocidade de 10 gibytes / s. Durará mais de um trilhão de vezes mais do que o universo existe). Embora eu acredite que os ciclos médios estejam superestimados, porque eles não nos dizem nada sobre a qualidade do gerador.

Aqui estão os resultados do benchmark:

Gerador	atuação	Qualidade	Correlação de sementes
SHISHUA	0,06	> 32 TiB	> 256 GiB
xoshiro256 + x8	0,09	1 KiB	0 KiB
ChaCha8	0,12	> 32 TiB?	> 32 TiB?
RomuTrio	0,31	> 32 TiB	1 KiB
xoshiro256 +	0,34	512 MiB	1 KiB
Wyrand	0,41	> 32 TiB	8 KiB
Lehmer128	0,44	> 32 TiB	1 KiB
RC4	7,48	1 TiB	1 KiB

Desempenho : o número de ciclos do processador gastos em um byte gerado. Recebido nas máquinas em nuvem N2 GCP e N2D (AMD), o pedido é o mesmo.
Qualidade : o nível em que o gerador falha no teste do PractRand. Se não falhar, há um sinal>. Se o resultado não for comprovado, há um ponto de interrogação.
Correlação de números de sementes : travessia PractRand com bytes alternados de oito fluxos com números de sementes 0, 1, 2, 4, 8, 16, 32, 64. Usamos o PractRand com convolução dupla e testes avançados.

Mais longe

Embora no nosso caso não haja problemas com a irreversibilidade, ainda podemos melhorar o SHISHUA.

Na minha opinião, o PRNG ideal tem as seguintes propriedades:

, 2¹⁰²⁴. 10 /. 10²⁸² , . «» ( ). , . , 128- NEON ARM? , , .
. , SHISHUA XOR . , .
, 2¹²⁸ ( ). SHISHUA, , . , ( ) (, , . 2).
A inicialização do estado tem dispersão perfeita : todos os bits do número inicial afetam todos os bits do estado com a mesma probabilidade. Eu quero descobrir em relação a SHISHUA.

Um dos problemas que atrasam o desenvolvimento de PRNGs e criptografia como um todo é a falta de melhores ferramentas de uso geral. Preciso de uma ferramenta que possa me fornecer imediatamente o resultado exato da medição, para que eu possa comparar diferentes arquiteturas em tempo real.

O PractRand é ótimo comparado ao que era antes, no entanto:

Não permite avaliar geradores de alta qualidade, tornando impossível compará-los entre si. Temos que dizer: "bem, depois de 32 terabytes eles não têm anomalias ..."
Demora semanas para executá-lo ...

Espero que a situação melhore muito em breve.

SHISHUA: o gerador de números pseudo-aleatórios mais rápido do mundo