Calidad de generadores lineales congruentes para números aleatorios

Estoy haciendo algunas simulaciones de la ecuación de Langevin, para varias fuerzas externas. Se le dijo que el C es rand()de stdlib.hpuede introducir un sesgo en los resultados de mi, estoy usando un Mersenne Twister.

Sin embargo, me gustaría saber (y ver) exactamente qué tipo de errores puede introducir un generador congruencial lineal en mi simulación. Estas son cosas que he probado:

• Generando tuplas 3D de randoms para tratar de ver hiperplanos. No puedo ver nada
• Haciendo la FFT de un gran vector de números aleatorios. Es casi lo mismo tanto para Mersenne Twister como para rand().
• Comprobación del principio de equipartición de una partícula en movimiento browniano. Ambos están de acuerdo integradores en el valor esperado de con el mismo número de dígitos significativos.$\langle \text{KE}\rangle=\frac{1}{2}k_BT$
• Al ver qué tan bien se agrupan en una serie de contenedores que no son dos poderosos. Ambos dan los mismos resultados cualitativos, nadie es mejor.
• En cuanto a los caminos browniano para ver claras divergencias de . De nuevo, sin suerte.$\langle x\rangle = 0$
• Distribución de puntos en un círculo. Lleno, y solo en el perímetro. Entre todos ellos y entre vecinos más cercanos (respuesta de Shor, más abajo en los comentarios). Disponible en este resumen , solo ejecútelo con Julia 0.5.0 después de instalar las bibliotecas necesarias (consulte el resumen para obtener instrucciones).

Me gustaría enfatizar que estoy buscando sesgos introducidos en el contexto de las simulaciones físicas. Por ejemplo, he visto cómo rand()falla miserablemente las pruebas de Dieharder mientras que el Mersenne Twister no lo hace, pero por el momento eso no significa demasiado para mí.

¿Tiene algún ejemplo físico, concreto, de cómo un mal generador de números aleatorios destruye una simulación de Montecarlo?

Nota: He visto cómo los PRNG RANDUpueden ser horribles. Me interesan los ejemplos no obvios, de generadores que parecen inocentes pero que en última instancia introducen sesgos.

Una referencia interesante que describe una falla de una simulación Monte Carlo de un sistema físico debido a un RNG inadecuado (aunque no usaron un LCG) es:

A. Ferrenberg y DP Landau. Simulaciones de Monte Carlo: errores ocultos de generadores de números aleatorios "buenos". Physical Review Letters 63 (23): 3382-3384, 1992.

Los modelos de Ising que estudiaron Ferrenberg y Landua son buenas pruebas de RNG porque se pueden comparar con una solución exacta (para el problema 2-D) y encontrar errores en los dígitos. Estos modelos deben mostrar las fallas en un antiguo PMMLCG aritmético de 32 bits sin demasiada dificultad.

Otra referencia interesante es:

H. Bauke y Stephan Mertens. Las monedas pseudoaleatorias muestran más caras que colas. arXiv: cond-mat / 0307138 [cond-mat.stat-mech]

Bauke y Mertens hacen un fuerte argumento en contra de generadores de números aleatorios de registro de desplazamiento de retroalimentación lineal binaria. Bauke y Mertens tienen algunos otros documentos relacionados con esto.

Puede ser difícil encontrar los planos Marsaglia en un diagrama de dispersión 3D. Puede intentar rotar la trama para obtener una mejor vista y, a veces, simplemente saldrán hacia usted. También puede hacer pruebas 3D de uniformidad estadística: para los LCG de 32 bits más antiguos, estos fallarán en un número bastante pequeño de contenedores. Por ejemplo, una prueba de uniformidad con una cuadrícula de contenedores de 20x20x20 en 3 dimensiones es suficiente para detectar la falta de uniformidad para el PMMLCG ampliamente utilizado (y previamente bien considerado) con m = 2 ^ 31-1, a = 7 ^ 5.

Es posible utilizar el conjunto de pruebas PRNG TestU01 para averiguar cuál de esas pruebas randfalla. (Consulte TestU01: Biblioteca de CA para pruebas empíricas de generadores de números aleatorios para obtener una descripción general del conjunto de pruebas). Eso es más fácil que crear simulaciones propias de Monte Carlo. En cierto modo, también se trata de la capacidad de compilación del software (y la corrección del software): dado un PRNG que parece funcionar bien en pruebas pequeñas y simples, ¿cómo sabe que sus comportamientos patológicos no serán activados por un programa más grande?

Aquí está el código:

#include "TestU01.h"

int main() {
  // Same as rand() on my machine
  unif01_Gen* gen = ulcg_CreateLCG(2147483647, 16807, 0, 12345);

  bbattery_SmallCrush(gen);
  bbattery_Crush(gen);

  return 0;
}

Para la suite SmallCrush , hay 3 pruebas que fallan de 15 (consulte la guía longtestu01.pdf en TestU01 para obtener descripciones largas y todas las referencias; estas son 15 estadísticas de 10 pruebas).

• $n$ $t$$d^t$$d^t$$I_1,\ldots$$\{I_{j+1}-I_j\}$

• $n$ $t$$[0,1)^t$$d^t$

• $n$$t$$[0,1)$$X$$n$$\mathbb{P}(X<x) = x^t$$n=2\times 10^6$$t=6$$\chi^2$$<10^{-300}$

Asumiendo que todas estas son simulaciones "típicas" de Monte Carlo (aunque podrían no ser como los problemas que tenía en mente), la conclusión es que randfalla algún subconjunto desconocido de ellas. Sin embargo, no sé por qué es específicamente ese subconjunto, por lo que es imposible para mí decir si funcionará con su propio problema o no.

MaxOft parece particularmente sospechoso, dada la sencillez de la descripción.

Entre las pruebas en la suite Crush , randfalla 51 de 140 (140 estadísticas en 96 pruebas). Algunas de las pruebas fallidas (como Fourier3 ) se realizan en cadenas de bits, por lo que tal vez sea posible que no sean relevantes para usted. Otra prueba curiosa que falla es GCD , que prueba la distribución del GCD de dos enteros aleatorios. (Nuevamente, no sé por qué falla esta prueba en particular o si su simulación sufrirá esto).

PD : Otra cosa a tener en cuenta es que en rand()realidad es más lento que algunos PRNG que superan con éxito todas las pruebas SmallCrush , Crush , BigCrush , como MRG32k3a (ver el documento L'Ecuyer y Simard más arriba).