¿Cómo puedo sembrar un generador de números pseudoaleatorio congruente lineal paralelo para un período máximo?

Normalmente cuando siembro un generador secuencial de números aleatorios en C, uso la llamada

srand(time(NULL))

luego usa

rand() mod N

para obtener un entero aleatorio entre 0 y N-1. Sin embargo, cuando hago esto en paralelo, las llamadas al tiempo (NULL) están tan cerca una de la otra que terminan siendo exactamente el mismo número.

He intentado usar un generador lineal de números aleatorios congruenciales:

$x_{n+1} := a x_n + c\;\;\; (\bmod m)$ para algunos enteros $a,c,$ y $m$ .

Sé que elegir $m=2^k$ para algún entero grande $k$ produce resultados rápidos porque el operador del módulo se puede calcular truncando dígitos. Sin embargo, me resulta difícil establecer semillas que produzcan secuencias aleatorias con un gran período en paralelo. Sé que la duración de un período es máxima si

1. c y m son números primos entre sí
2. a-1 es divisible por todos los factores primos de m
3. si m es un múltiplo de 4, a-1 también debe ser un múltiplo de 4.

(fuente: wikipedia )

Pero, ¿cómo me aseguro de que todas las secuencias de números aleatorios tengan esta propiedad máxima? En términos de MPI, ¿cómo incorporo ranky sizepara producir períodos máximos utilizando el método lineal congruencial? ¿Sería más fácil usar un Fibonacci rezagado o un Mersenne Twister para producir flujos aleatorios paralelos más largos?

El truco es intercalar la secuencia LCG de cada proceso: para procesos, modificamos el LCG$p$

ser - estar

donde y efectivamente avanzan pasos. Podemos derivarlos rápidamente expandiendo el paso LCG original:$A_p$$C_p$$p$

y el patrón es que y es el resultado de, comenzando con el número , multiplicando sucesivamente por y sumando veces, luego multiplicando por , todo .$A_p \equiv a^p \mod m,$$C_p$$0$$a$$1$ $p$$c$$\bmod m$

El último paso es asegurarse de que el LCG -stride de cada proceso no se solape: simplemente inicialice el proceso con rango con y un LCG paralelo con períodos individuales de esté listo, donde es el período original y es el Número de procesos. Si el LCG modificado de cada proceso se usa por igual, el período completo se recupera en paralelo.$p$$r$$x_r$$N/p$$N$$p$

Implementé esto hace unos seis meses (quizás ingenuamente), y puedes encontrar el código aquí .

Hay un muy buen documento de resumen de Katzgrabber, Random Numbers In Scientific Computing: An Introduction , al que apunto a las personas que desean ser usuarios de PRNG para el cómputo científico. Los generadores lineales congruentes son rápidos, pero eso es todo lo que tienen para ellos; tienen períodos cortos y pueden equivocarse muy fácilmente; combinaciones perfectamente razonables de a, c y m pueden terminar con resultados horriblemente correlacionados, incluso si cumple con los requisitos habituales entre a, c y m.

Peor aún, en un caso común donde m es una potencia de dos (por lo que la operación de modificación es rápida), los bits de orden inferior tienen un período mucho más corto que la secuencia en su conjunto, por lo que si está haciendo rand ()% N, tienes un período aún más corto de lo que esperas.

Como regla general, los generadores de fibbonacci rezagado, MT y WELL tienen propiedades mucho mejores y siguen siendo bastante rápidos.

En términos de siembra en paralelo, el método de Jack Poulson es bueno porque proporciona una secuencia bien definida de números divididos equitativamente entre los procesadores. Si eso no importa, puede hacer cualquier cosa razonable para sembrar los diferentes PRNG; el mismo documento sugiere algo que mucha gente ha inventado de forma independiente, cambiando el número de tarea PID o MPI con el tiempo. La fórmula particular sugerida allí es

long seedgen(void)  {
    long s, seed, pid;

    pid = getpid();
    s = time ( &seconds );

    seed = abs(((s*181)*((pid-83)*359))%104729); 
    return seed;
}

No tengo opiniones particulares sobre esa implementación específica, pero el enfoque general es ciertamente razonable.

Una idea simple para difundir un RNG secuencial típico sobre un número decente de subprocesos es hacer que un solo subproceso avance la semilla lo más rápido posible y envíe solo cada milésima semilla a la memoria. Luego, haga que cada uno de sus otros subprocesos recoja una de estas semillas de referencia espaciadas y procese los 1000 valores en ese bloque, es decir, regenere nuevamente las 1000 semillas en el bloque, genere sus sorteos psuedorandom y luego haga cualquier otro procesamiento que su tarea implique.

Esto funciona porque para los RNG que no computan tanto (LCG es ciertamente uno, pero muchos otros deberían estar en esta categoría) el verdadero cuello de botella está enviando las semillas a la memoria (y probablemente también el procesamiento posterior). Si ejecuta un LCG sin enviar nada a la memoria, todo debería permanecer en los registros de la CPU y ser extremadamente rápido. Incluso para un RNG más complicado, debe permanecer en el caché L1 y ser muy rápido.

He usado este enfoque muy simple con un LCG que por razones heredadas tenemos que mantener. Básicamente, obtenemos una aceleración lineal de hasta 4-8 hilos en una estación de trabajo multinúcleo típica. Pero ahora probaré el método de la respuesta de Jack Poulson y espero ser aún más rápido :).

OTOH, creo que este simple truco debería funcionar para otros RNG inherentemente secuenciales.

Esta respuesta llega tarde, pero deberías echar un vistazo a SPRNG . Está específicamente diseñado para la escalabilidad en paralelo y admite un puñado de tipos de PRNG.