mt19937 c++ example

La mejor manera de sembrar mt19937_64 para simulaciones de Monte Carlo

Question

Estoy trabajando en un programa que ejecuta la simulación Monte Carlo; específicamente, estoy usando un algoritmo de Metrópolis. El programa necesita generar posiblemente miles de millones de números "aleatorios". Sé que el tornado de Mersenne es muy popular para la simulación de Monte Carlo, pero me gustaría asegurarme de que estoy sembrando el generador de la mejor manera posible.

Actualmente estoy calculando una semilla de 32 bits usando el siguiente método:

mt19937_64 prng; //pseudo random number generator
unsigned long seed; //store seed so that every run can follow the same sequence
unsigned char seed_count; //to help keep seeds from repeating because of temporal proximity

unsigned long genSeed() {
    return (  static_cast<unsigned long>(time(NULL))      << 16 )
         | ( (static_cast<unsigned long>(clock()) & 0xFF) << 8  )
         | ( (static_cast<unsigned long>(seed_count++) & 0xFF) );
}

//...

seed = genSeed();
prng.seed(seed);

Tengo la sensación de que hay formas mucho mejores de asegurar semillas nuevas que no se repiten, y estoy bastante seguro de que mt19937_64 puede sembrarse con más de 32 bits. ¿Alguien tiene alguna sugerencia?

Answer 1

De los comentarios entiendo que quiere ejecutar varias instancias del algoritmo, una instancia por hilo. Y dado que la semilla para cada instancia se generará prácticamente al mismo tiempo, debe asegurarse de que estas semillas sean diferentes. Si eso es realmente lo que estás tratando de resolver, entonces tu función genSeed no necesariamente garantizará eso.

En mi opinión, lo que necesitas es un generador de números aleatorios paralelos (RNG). Lo que esto significa es que solo necesita un RNG que crea una instancia con una sola semilla (que puede generar con su genSeed) y luego la secuencia de números aleatorios que normalmente se gerena en un entorno secuencial se divide en X sin superposición secuencias; donde X es el número de hilos. Hay una biblioteca muy buena que proporciona este tipo de RNG en C ++, sigue el estándar de C ++ para RNG y se llama TRNG ( http://numbercrunch.de/trng ).

Aquí hay un poco más de información. Hay dos formas de lograr secuencias no superpuestas por hilo. Supongamos que la secuencia de números aleatorios de un solo RNG es r = {r (1), r (2), r (3), ...} y solo tiene dos hilos. Si sabe de antemano cuántos números aleatorios necesitará por hilo, digamos M, puede dar el primer M de la secuencia r al primer hilo, es decir, {r (1), r (2), ..., r (M)}, y el segundo M al segundo hilo, es decir, {r (M + 1), r (M + 2), ... r (2M)}. Esta técnica se llama división de bloques ya que divide la secuencia en dos bloques consecutivos.

La segunda forma es crear la secuencia para el primer hilo como {r (1), r (3), r (5), ...} y para el segundo hilo como {r (2), r (4), r (6), ...}, que tiene la ventaja de que no necesita saber de antemano cuántos números aleatorios necesitará por hilo. Esto se llama saltar en salto.

Tenga en cuenta que ambos métodos garantizan que las secuencias por hilo no se superponen. El enlace que publiqué arriba tiene muchos ejemplos y la biblioteca en sí es extremadamente fácil de usar. Espero que mi publicación te ayude

Answer 2

La función POSIX gettimeofday(2) proporciona el tiempo con una precisión de microsegundos.

La función de hilo POSIX gettid(2) devuelve el número de ID del hilo actual.

Debería poder combinar el tiempo en segundos desde la época (que ya está usando), el tiempo en microsegundos y la identificación del hilo para obtener una semilla que siempre es única en una máquina.

Si también necesita que sea único en varias máquinas, también podría considerar obtener el nombre de host, la dirección IP o la dirección MAC.

Supongo que 32 bits es probablemente suficiente, ya que hay más de 4 mil millones de semillas únicas disponibles. A menos que esté ejecutando miles de millones de procesos, lo que no parece probable, debería estar bien sin tener que recurrir a semillas de 64 bits.

Answer 3

Por lo que puedo decir de sus comentarios, parece que lo que le interesa es asegurarse de que si un proceso inicia varias de sus simulaciones al mismo tiempo, obtendrá diferentes semillas.

El único problema importante que puedo ver con su enfoque actual es una condición de carrera: si va a iniciar múltiples simulaciones simultáneamente, debe hacerlo desde hilos separados. Si se hace desde hilos separados, debe actualizar seed_count de forma segura para la ejecución de subprocesos, o múltiples simulaciones podrían terminar con el mismo seed_count . Simplemente puede hacer que sea std::atomic<int> para resolver eso.

Más allá de eso, parece más complicado de lo que debe ser. ¿Qué ganas usando dos temporizadores separados? Podrías hacer algo tan simple como esto:

1. al inicio del programa, tome la hora actual del sistema (usando un temporizador de alta resolución) una vez y almacénela.
2. asigne a cada simulación una ID única (esto podría ser un número entero inicializado a 0, (que debe generarse sin condiciones de carrera, como se mencionó anteriormente) que se incrementa cada vez que se inicia una simulación, efectivamente como su cuenta de seed_count .
3. cuando siembra una simulación, solo use la marca de tiempo generada inicialmente + la ID única. Si haces esto, cada simulación en el proceso tiene asegurada una semilla única.

Answer 4

Qué tal si...

Hay un código principal que inicia los hilos y hay copias de una función ejecutada en esos hilos, cada copia con su propio Marsenne Twister. ¿Estoy en lo correcto? Si es así, ¿por qué no usar otro generador aleatorio en el código principal? Sería sembrado con sello de tiempo, y enviaría sus números pseudoaleatorios consecutivos para que las instancias funcionen como sus semillas.

Answer 5

Recapitulemos (también los comentarios), queremos generar diferentes semillas para obtener secuencias independientes de números aleatorios en cada una de las siguientes ocurrencias:

1. El programa se relanza en la misma máquina más tarde,
2. Dos hilos se lanzan en la misma máquina al mismo tiempo,
3. El programa se lanza en dos máquinas diferentes al mismo tiempo.

1 se resuelve utilizando el tiempo desde epoch, 2 se resuelve con un contador atómico global, 3 se resuelve con un id de plataforma dependiente (ver Cómo obtener (casi) identificador de sistema único en una plataforma cruzada? )

Ahora el punto es ¿cuál es la mejor manera de combinarlos para obtener uint_fast64_t (el tipo de std::mt19937_64 de std::mt19937_64 )? Supongo aquí que no conocemos a priori el rango de cada parámetro o que son demasiado grandes, por lo que no podemos simplemente jugar con los cambios de bits obteniendo una semilla única de una manera trivial.

Un std::seed_seq sería el camino más fácil, sin embargo, su tipo de devolución uint_least32_t no es nuestra mejor opción.

Un buen hasher de 64 bits es una opción mucho mejor. El STL ofrece std::hash bajo el encabezado functional , una posibilidad es concatenar los tres números anteriores en una cadena y luego pasarlo al hasher. El tipo de devolución es un size_t que en 64 máquinas es muy probable que coincida con nuestros requisitos.

Las colisiones son poco probables pero, por supuesto, es posible. Si desea asegurarse de no generar estadísticas que incluyan una secuencia más de una vez, solo puede almacenar las semillas y descartar las ejecuciones duplicadas.

Un std::random_device también se podría usar para generar las semillas (las colisiones aún pueden ocurrir, es difícil decirlo más o menos a menudo), sin embargo, dado que la implementación depende de la biblioteca y puede descender a un generador pseudo aleatorio, es obligatorio verifique la entropía del dispositivo y evite utilizar un dispositivo de entropía cero para este fin, ya que probablemente romperá los puntos anteriores (especialmente el punto 3). Lamentablemente, puede descubrir la entropía solo cuando lleva el programa a la máquina específica y prueba con la biblioteca instalada.

Answer 6

Use std::random_device para generar la semilla. Proporcionará números aleatorios no deterministas, siempre que su implementación lo admita. De lo contrario, está permitido usar algún otro motor de números aleatorios.

std::mt19937_64 prng;
seed = std::random_device{}();
prng.seed(seed);

operator() de std::random_device devuelve un unsigned int , por lo tanto, si su plataforma tiene int 32 bits y desea una semilla de 64 bits, deberá llamarla dos veces.

std::mt19937_64 prng;
std::random_device device;
seed = (static_cast<uint64_t>(device()) << 32) | device();
prng.seed(seed);

Otra opción disponible es usar std::seed_seq para sembrar el PRNG. Esto permite que el PRNG llame a seed_seq::generate , que produce una secuencia no sesgada en el rango [0 ≤ i <2 ³² ) , con un rango de salida lo suficientemente grande como para ocupar todo su estado.

std::mt19937_64 prng;
std::random_device device;
std::seed_seq seq{device(), device(), device(), device()};
prng.seed(seq);

random_device al random_device 4 veces para crear una secuencia inicial de 4 elementos para seed_seq . Sin embargo, no estoy seguro de cuál es la mejor práctica para esto, en cuanto a la longitud o fuente de los elementos en la secuencia inicial.