Números aleatorios, radiaciones y métodos de Monte Carlo

1.- Introducción

Como bien nos explicó Salvador García en El Casino de Monte Carlo y las Radiaciones, el método de Monte Carlo es una herramienta numérica de gran utilidad para resolver problemas complejos de naturaleza probabilística en muchas ramas de la ciencia y la ingeniería. De su análisis quedaba claro que para su uso práctico eran necesarios recursos computacionales potentes por un lado, y generadores de números aleatorios rápidos y eficientes por otro.

En la entrada de hoy voy a analizar más en profundidad el asunto de los generadores de números aleatorios, su naturaleza y sus limitaciones.

2.- Generadores de Números Aleatorios Verdaderos (GNAV) y de Números PeudoAleatorios (GNPA)

Los generadores de números aleatorios ideales o verdaderos serían aquellos sistemas o algoritmos informáticos capaces de producir secuencias infinitas de números uniformemente distribuidos, no correlacionados, que superaran cualquier prueba estadística de aleatoriedad y con periodos de repetición infinitos.

¿Existen realmente esos generadores o son una leyenda urbana más?

Ciertamente existen los generadores de números aleatorios verdaderos (GNAV). Son aquellos que utilizan un fenómeno físico aleatorio conocido (por ejemplo el ruido electrónico o el decaimiento radiactivo) y lo incorporan a un ordenador que es el que en definitiva produce el número aleatorio. Los GNAV se emplean en simulaciones, seguridad informática, criptografía, loterías y juegos de azar, por ejemplo.

Curiosamente, el Sistema Nacional de Loterías del Reino Unido utiliza un equipo llamado ERNIE (Electronic Random Number Indicator Equipment) para producir los números de la lotería nacional. La versión actual (ERNIE 4) utiliza el ruido térmico de los transistores para generar números aleatorios verdaderos. La versión original del año 1957 utilizaba el ruido generado por tubos de gas neón. Esta notable diferencia entre nuestro sistema de loterías de bombo con bolas y el británico de fluctuaciones de la naturaleza, y su relación con la idiosincrasia de cada país, seguro que daría para un debate más que interesante en el café de la mañana.

De acuerdo con lo anterior, podríamos construir nuestro propio ERNIE con una fuente radiactiva de Cs-137 para generar números aleatorios y atacar rigurosamente cualquier problema con cálculos de Monte Carlo. Sin embargo, la baja eficiencia y la naturaleza no determinista de los GNAV imponen serias limitaciones que hacen que no sea factible su uso en simulaciones de Monte Carlo. La baja eficiencia hace que sean demasiado lentos para realizar cálculos computacionales en tiempos razonables y la naturaleza no determinista hace que las secuencias generadas no sean reproducibles y, por lo tanto, hacen muy complejo un análisis de su influencia en el cálculo.

Con la llegada de los primeros ordenadores los programadores reconocieron la necesidad de introducir los números aleatorios en sus programas y crearon los primeros algoritmos generadores de números pseudoaleatorios (GNPA). Dado que un ordenador sigue ciegamente las instrucciones que se le dan de una manera secuencial y totalmente predecible nunca podrán producir por sí solo números aleatorios verdaderos. Lo que generan realmente son números casi aleatorios (pseudoaleatorios) que no están uniformemente distribuidos, están correlacionados, son perfectamente reproducibles a partir de un número inicial que se toma como semilla y tienen periodos de repetición finitos. A pesar de estas limitaciones, que hasta cierto punto se pueden superar, los GNPA son mucho más eficientes (rápidos) estando mejor adaptados que GNAV para su implantación en programas que realizan cálculos de Monte Carlo.

3.- Algoritmos para generar números pseudoaleatorios

Un algoritmo ideal de números pseudoaleatorios sería uno que fuera fácil de implementar, rápido, que consumiera pocos recursos computacionales y que produjera secuencias de números uniformemente distribuidos, con una correlación baja y con periodos de repetición muy largos.

La mayoría de los generadores de números pseudoaleatorios están basados en Generadores de Congruencia Lineal (GCL), Generadores de Fibonacci (GF) o combinaciones y variaciones de ambos. Todos ellos emplean relaciones de recurrencia del tipo Ni= F(Ni-1,…, Ni-k) que parten de un número inicial llamado semilla de tal forma que para generar un nuevo número aleatorio Ni se emplean uno o varios de los números generados con anterioridad.

i) Generadores de Congruencia Lineal (GCL)

Son los más sencillos y comunes. Se basan en la siguiente relación de recurrencia:
Ni+1 = (aNi + c) mod m, donde a, c y m son constantes y N0 sería la semilla.

El resultado de una operación mod m es el resto del cociente (aNi + c)/m, por lo tanto hay m valores posibles de Ni. La secuencia es periódica y el período de repetición depende de las constantes a, c y m, así como de la semilla N0. Una buena elección de estas constantes no es trivial y es lo que determina la validez del GNPA.

Si alguien quiere pasar el rato puede entretenerse con este GCL que ha demostrado ser muy eficiente:

Ni+1 =65539 Ni mod 231

Las ventajas más importantes de los GCL son su sencillez, su bajo consumo de recursos computacionales y facilidad de implementación. Además, son fáciles de paralelizar lo que permite dividir el proceso de generación de números aleatorios en subprocesos que se pueden calcular de manera simultánea (paralela) en varios procesadores ahorrando así tiempo de cálculo.

Las limitaciones más importantes de los GCL son su periodicidad (normalmente no suele ser más grande de 232), que depende de la semilla inicial, y la presencia de correlaciones entre términos consecutivos de la secuencia de los números generados. Existen técnicas como desordenar la secuencia generada, por ejemplo, para minimizar estos problemas. Y, por supuesto, con las nuevas arquitecturas de 64-bits los periodos se pueden ampliar hasta 264.

ii) Generadores de Fibonacci (GF)

Se basan en relaciones de recurrencia del tipo:
Ni+1 = (Ni-r * Ni-s) mod m, donde r<s son enteros dados, * representa alguna de las operaciones binarias (+, -, x, ^) y son necesarias s semillas.

Este tipo de generador requiere el uso de otro generador para iniciarse y se tiene que mantener una lista de los últimos s números generados.

Aunque los GF son más complejos de implantar y requieren más recursos computacionales que los GCL, son mucho más rápidos y tienen periodos de repetición más largos.

Me gustaría recalcar que la creación de algoritmos para obtener buenos GNPA tiene una parte empírica fundamental que reside en la elección de las constantes a, c, m, r, s y de las operaciones binarias en las relaciones de recurrencia. Pero, ¿cómo podemos saber si unas son más adecuadas que otras? La única manera que hay es poniendo a prueba los algoritmos y para ello existen diversas pruebas de aleatoriedad que es necesario realizar antes de aceptar un nuevo generador.

4.- Pruebas de aleatoriedad

Existen pruebas de uniformidad (Figura 1), correlación y análisis espectral para analizar las propiedades de las series de números pseudoaleatorios. En principio, si un GNPA supera estas pruebas parece que podríamos asegurar que estamos ante un generador de números aleatorios bastante competente.

Figura 1: Representación bidimensional de la uniformidad para un GNAV (arriba) y para un mal GNPA (abajo). Fuente: www.random.org

Figura 1: Representación bidimensional de la uniformidad para un GNAV (arriba) y para un mal GNPA (abajo). Fuente: www.random.org

Sin embargo, a veces esto no es cierto y para ponerlo de manifiesto se pueden emplear precisamente los mismos métodos de Monte Carlo que utilizan los GNPA. La idea es muy sencilla. Supongamos que, por un lado, conocemos la solución analítica exacta de un problema estadístico determinado y que, por otro, resolvemos ese mismo problema con cálculos de Monte Carlo. Entonces, por simple comparación podríamos determinar fácilmente la influencia de los distintos GNPA en el cálculo.

¿Existe un problema estadístico real con una solución analítica exacta? Pues sí y se conoce como Modelo de Ising 2D (Figura 2). Éste es un modelo propuesto para estudiar el ferromagnetismo de la materia y consiste en una red bidimensional cuadrada de espines que sólo pueden tener dos valores + y – (arriba y abajo) y que interaccionan con sus vecinos próximos según el hamiltoniano cuántico correspondiente. A diferencia del Modelo de Ising 1D (unidimensional), éste sí presenta propiedades físicas útiles y reales como la existencia de una transición de fase ferromagnética a una determinada temperatura crítica Tc. Onsager fue el primero en determinar la función de partición Z correspondiente a este modelo 2D (lo que le sirvió, entre otras aportaciones, para recibir el premio Nobel) y a partir de Z el resto de propiedades termodinámicas: Energía, calor específico, entropía, etc. Para un modelo 3D general aún no se conoce una solución analítica exacta y ni siquiera está del todo claro si existe.

3

Figura 2: Red cuadrada de espines (+,-) del modelo bidimensional de Ising (arriba) y representación de la magnetización, susceptibilidad magnética y capacidad calorífica obtenidos a partir del modelo (abajo)

Figura 2: Red cuadrada de espines (+,-) del modelo bidimensional de Ising (arriba) y representación de la magnetización, susceptibilidad magnética y capacidad calorífica obtenidos a partir del modelo (abajo)

Así pues, estudiando las diferencias entre el valor calculado y teórico de la energía o el calor específico, por ejemplo, para un número dado de iteraciones de Monte Carlo uno puede evaluar fácilmente los diferentes GNPA y encontrarse con sorpresas desagradables que no quedaban al descubierto con otras pruebas: Generadores que superan las pruebas de uniformidad y correlación presentan problemas de convergencia, por ejemplo, en un cálculo de Ising 2D.

5.- Conclusiones

La naturaleza aleatoria de las desintegraciones radiactivas, de las que tantas veces contamos cosas en Desayuno con Fotones, sirve para generar números aleatorios verdaderos. A su vez estos números aleatorios son los que nos permiten realizar cálculos de Monte Carlo para simular el transporte de la radiación en la materia y conocer, por ejemplo, la dosis de radiación depositada en un tejido biológico.

Sin embargo, en la práctica resulta mucho más eficiente trabajar con números que a pesar de no ser realmente aleatorios, nos funcionan casi como si lo fueran. No es una tarea fácil conocer cómo y cuándo pueden fallarnos y por ello esto es un campo de investigación abierto y en pleno desarrollo. Aquí he tratado de ilustrar algunas de las herramientas que existen para ayudarnos a entender un poco mejor la naturaleza del número pseudoaleatorio.

Para saber más

En esta entrada sólo he mostrado algunas de las ideas y conceptos más sencillos sin profundizar demasiado en ellos. Por ello, al que tenga curiosidad le recomiendo las siguientes lecturas:

Random.org: https://www.random.org/analysis/

Random Numer Generation and Monte Carlo Methos, Gentle, James E., Ed. Springer-Verlag New York, 2003

Random Number Generators: http://random.mat.sbg.ac.at/links/rando.html

4 Respuestas a “Números aleatorios, radiaciones y métodos de Monte Carlo

  1. Buen post, Agustín, creo que ni siquiera los más familiarizados con Monte Carlo tenemos siempre una visión correcta de este asunto de fondo tan relevante para nuestros cálculos. Al menos en mi caso, reconozco no haberle prestado toda la atención que debía.

  2. Gracias Manolo. Transmitir asuntos tan técnicos como éste no es nada fácil, pero muchos de nuestros cálculos se apoyan en ellos y creo que es importante al menos tener un conocimiento básico de su existencia.

  3. Este post parece un paper!!

    Por qué se les llama generadores de Fibonacci? Tras una lectura rápida no veo la relación con la famosa serie…

    • Ja ja ja, un paper, está oxidado ;-)

      Buena pregunta Igor. El verdadero nombre es Lagged Fibonacci Generator, pero por simplicidad he omitido el término Lagged (retardado).

      La relación de recurrencia correspondiente es una generalización de la serie de Fibonacci en la que se utilizan dos términos anteriores para generar uno nuevo: Xn=Xn-1 + Xn-2

      Saludos.

Responder

Introduce tus datos o haz clic en un icono para iniciar sesión:

Logo de WordPress.com

Estás comentando usando tu cuenta de WordPress.com. Cerrar sesión / Cambiar )

Imagen de Twitter

Estás comentando usando tu cuenta de Twitter. Cerrar sesión / Cambiar )

Foto de Facebook

Estás comentando usando tu cuenta de Facebook. Cerrar sesión / Cambiar )

Google+ photo

Estás comentando usando tu cuenta de Google+. Cerrar sesión / Cambiar )

Conectando a %s