Modelo numérico de Ising: algoritmo de Wolff y correlaciones

Estoy haciendo un análisis numérico de Monte Carlo en el modelo de Ising bidimensional en el punto crítico. Estaba usando la evolución de 'single flip' de Metropolis al principio con éxito, aunque sufre una ralentización crítica y hace que sea poco probable que sea posible estudiar redes grandes. Ahora estoy mirando los algoritmos de cambio de clúster, específicamente el algoritmo de Wolff.

Me las arreglé para implementarlo, y parece estar funcionando como debería (volteando un giro único en T = + , toda la red en T = 0 , coincide con la densidad de energía correcta en el límite termodinámico...) pero no obtengo el comportamiento correcto para los dos puntos < σ i σ j > función de correlación.

Según CFT, debería comportarse como:

< σ i σ j > 1 | i j | 1 4

Cada vez estoy más convencido de que tiene que ver con las condiciones de contorno, uso límites libres no periódicos. La literatura sobre el tema no dice mucho sobre este punto.

¿Me estoy perdiendo una sutileza (o una evidencia) en este procedimiento, o en el uso de este algoritmo?

Probablemente debería migrar a uno de los sitios de código, pero sin las implementaciones de su algoritmo es imposible decir dónde (o si) se equivocó.
No, esto es apropiado para un sitio de física. Sólo un físico puede responder a esto.
@CarlWitthoft: Realmente no creo que sea la implementación: primero lo hice por mi cuenta y lo verifiqué doble, triple, etc. luego copié una versión de Internet ( physics.buffalo.edu/phy411-506-2004/Topic3/topic3-lec7.pdf ) y todavía no funciona (al nivel de dar buenas correlaciones).
Eso es también lo que me lleva a pensar que tiene que ver con el propio algoritmo y la información que utiliza al componer clústeres realmente grandes como los que aparecen en la criticidad.
Para cualquier tipo de sutileza en la simulación del modelo de Ising, puedo recomendar "Simulación de Monte Carlo en física estadística" de Binder y Heermann. Aparte de eso, ¿has probado con condiciones periódicas?
@Alexander: Gracias por la referencia, lo investigaré. De lo contrario, sí, desde que abrí esta pregunta, puse la red en un toro (condiciones de contorno periódicas en ambas direcciones) y da resultados totalmente diferentes pero aún incorrectos. fui de Δ σ = 0,4 a 0,05, siendo el resultado correcto 0,125. También estoy investigando eso ahora.
@suresh no, si solo un físico puede responderlo es irrelevante si se trata de un tema aquí. El criterio relevante es si se trata de física .
@DavidZ Sí, se trata de física.
@suresh no le hace ningún bien a nadie simplemente decir algo sin ninguna justificación.
¿Obtienes un comportamiento sensato lejos del punto crítico? Las correlaciones espín-espín también se conocen para estos límites y convergerían más rápido y serían menos susceptibles a problemas de tamaño finito.
@DavidZ es difícil imaginar cómo se podría decir que una pregunta sobre las correlaciones en el modelo de Ising no se trata de física.
@Aprendiendo qué tan grande es su cuadrícula? Dado que el rango espacial de las correlaciones también diverge en el punto crítico, las condiciones de contorno tendrán un fuerte efecto en todo el sistema. Supongo que esto se puede mitigar un poco en la práctica haciéndolo muy grande, por lo que tengo curiosidad sobre el tamaño de su sistema.
@Nathaniel Sí, pero esta pregunta puede ser sobre el efecto de las condiciones de contorno en el algoritmo de Wolff, y las correlaciones en el modelo de Ising son simplemente el contexto en el que surgió la pregunta.
@DavidZ como experto en mecánica estadística, simular numéricamente el modelo de Ising es una de las primeras cosas en las que pienso cuando alguien dice "física". El centro de ayuda dice que las preguntas sobre "diseños y resultados experimentales" están en el tema, y ​​aunque uno puede entrar en argumentos filosóficos sobre si los resultados computacionales son "empíricos", me parecería absurdo decir que un algoritmo para probar las implicaciones de la El modelo de Ising es menos relevante que un experimento para medir las correlaciones en un sistema de espín real. O un enfoque más analítico para el caso.
Por supuesto, acordamos como comunidad que las preguntas sobre software están fuera de tema. Sin embargo, una pregunta sobre los resultados obtenidos usando un algoritmo en particular es bastante diferente de una pregunta sobre una pieza de software en particular .
@Nathaniel, tengo entendido que los algoritmos numéricos se incluyen específicamente en la categoría de preguntas computacionales que acordamos que estaban fuera de tema.
@DavidZ, ¿es esta la discusión relevante? meta.physics.stackexchange.com/questions/43/… Me parece, por la respuesta aceptada y la pregunta en sí, que están específicamente excluidos .
@DavidZ Del resumen en la meta-pregunta: " parece que estamos tomando la posición de que las preguntas sobre la interpretación o justificación de un algoritmo o sus resultados, o sobre el diseño de algoritmos motivado físicamente, están bien ". Si alguna pregunta alguna vez cayó en esa categoría, esta lo hace.
@Nathaniel Bueno, en realidad creo que este es relevante, ya que es más nuevo, pero la conclusión es la misma. Así que tienes razón. Probablemente deberíamos actualizar el Centro de ayuda en consecuencia.
@DavidZ Como recién llegado a este foro, no estoy seguro de cómo se prueba que algo es física. ¿Necesito mostrar mis credenciales o dedicar tiempo a explicar cómo el modelo de Ising es un problema clásico de mecánica estadística y su simulación, especialmente en el caso 3D, es un enfoque importante para resolver esta clase de problemas de universalidad?
@suresh Nadie discute que el modelo de Ising es física, pero esta pregunta no se trata realmente del modelo de Ising, se trata del comportamiento de un algoritmo numérico particular. Para "probar" que la pregunta es sobre física, recurre a las políticas relevantes como se describe en las publicaciones meta y en el centro de ayuda, como lo hizo Nathaniel. Si desea obtener más detalles, comuníquese conmigo en Physics Chat , estaré encantado de aclarar esto más.

Respuestas (2)

Intenté reproducir su problema usando mi propio código, pero no pude: obtuve el valor correcto para el exponente. No puedo decirte qué salió mal en tu cálculo, ¡pero puedo decirte exactamente lo que hice!

Puedes inspeccionar mi código en GitHub . Fue lo primero que escribí en C++ y hoy haría muchas cosas diferentes, así que no juzgues el código con demasiada severidad. Hay un archivo Léame en el repositorio que explica casi todo al respecto.

Cargué mi script de ejecución, la parte relevante de los resultados correspondientes y el pequeño script de Pyxplot que encaja en este GitHub Gist . Mi valor para el exponente es:

β = 0.243 ± 0.001
Puede ver el ajuste trazado junto con los datos de Monte Carlo en el archivo sscorr_fit.pdf .

Aquí hay una lista de cosas que creo que podrían ser relevantes.

  • Estoy seguro de que lo sabe, pero solo para que conste: en la ecuación de la pregunta, no se supone que uno tome el valor absoluto de la diferencia de los índices. i y j de los dos sitios. Los índices que asignas a los sitios son completamente arbitrarios. Lo que importa es la distancia entre los dos sitios, por lo que preferiría escribirlo así:

    s i s j 1 | r j r i | 1 4

  • Solo medí la correlación a lo largo de la dirección de los vectores de red y no a lo largo de ninguna diagonal. Sin embargo, eso no debería importar, ya que existe una prueba de que la correlación espín-espín es rotacionalmente simétrica en el punto crítico.

  • Solo encajé tu formulario en el rango 1 | r j r i | 50 para un modelo de 512x512. Creo que si vas más allá, obtienes demasiada influencia en las condiciones de contorno periódicas. Definitivamente no puedes encajar todo 1 | r j r i | 256 rango, ya que la correlación debe tener un mínimo en | r j r i | = 256 .

  • La función de correlación cambia bastante rápido con la temperatura alrededor del punto crítico, así que asegúrese de tener la temperatura correcta y que todo esté equilibrado. Según mi experiencia, es mejor comenzar el algoritmo de Wolff con un modelo en el que todos los giros apunten en la misma dirección. Si realiza giros aleatorios y luego se enfría desde T = a T C el algoritmo de Wolff es extremadamente ineficiente al principio, ya que no puede construir grandes grupos en este ruido de giros aleatorios.

¡Espero que ayude!

Muchas gracias por la entrada, voy a mirar más a fondo su código. ¿Puedo preguntarle cuántos sistemas simuló (es decir, está promediando)? Logré derivar un buen coeficiente en mi caso también (muy cerca de 0.249) pero mis datos eran más ruidosos, probablemente la suma utilizada para ajustar (linealmente en escala logarítmica) canceló la mayor parte del ruido. Además, ¿cuántos pasos de evolución de Wolff realiza para termalizar a partir de configuraciones T = 0 ordenadas? Muchas gracias de nuevo ;)
@Learningisamess: Permítanme reiterar lo que dije como comentario a la otra respuesta. Puede evitar por completo el problema de determinar la cantidad de pasos necesarios para termalizar: simplemente use un algoritmo de simulación perfecto ; cuando este último se detiene, ¡tiene la garantía de que el sistema se ha termalizado! En este artículo se describe uno de estos algoritmos (de clúster): M. Huber, abounding chain for Swendsen-Wang , Random Structures & Algorithms, 22(1):43–59, 2003. Es más complicado de implementar que el algoritmo estándar, pero es muy eficiente y obtienes todos los beneficios de una simulación perfecta.
Todo lo que publiqué se obtuvo de una sola simulación, cuyo script de ejecución es parte del repositorio de git que vinculé . Así que no hice ningún promedio más allá del promedio dentro de la propia cadena de Markov. Realicé algunas simulaciones a temperaturas ligeramente diferentes, pero no incluí sus resultados en lo que publiqué anteriormente.
Para el equilibrio utilicé 100 pasos de Monte Carlo (MCS). Se supone que un MCS debe voltear aproximadamente una cantidad de giros igual a la cantidad de sitios de celosía. La cantidad de clústeres de Wolff que se necesitan para eso se determina en tiempo de ejecución midiendo el tamaño promedio de los últimos clústeres. El objetivo de todo este negocio es hacer que el MCS sea independiente del número de sitios de celosía, de modo que el tiempo de autocorrelación medido en MCS sea solo un número constante que no dependa del tamaño del sistema. Verifique la configuración en el script de ejecución con el archivo Léame de SSMC para ver qué hice exactamente.
Interesante, tengo una forma diferente de alcanzar el equilibrio. Primero monitorearía el valor de la energía media por giro (para cada configuración) y leería a partir de ahí cuántos pasos necesito para cada tamaño de sistema (generalmente me atengo a 3,4 tamaños también ). Comprobé si partir de una configuración totalmente ordenada es más rápido que una totalmente desordenada y, curiosamente, depende de las condiciones de contorno: con límites fijos, desordenada es muy, muy ventajosa (la razón por la que el algoritmo de Wolff debe modificarse con una configuración fija) . límites con un factor de aceptación 'à la Metropolis'
dependiendo del cambio de energía desde la ruptura de enlaces hasta el límite y con límites fijos se produce mucho rechazo), pero con límites libres, un punto de inicio de T bajo es aproximadamente 4 veces más rápido.
  1. Es solo exactamente a la temperatura crítica que este resultado CFT funciona. No has mencionado si has usado la temperatura crítica cuando hiciste el monte-carlo.

  2. En/cerca del punto crítico, el tiempo de autocorrelación se vuelve enorme. (Si no me equivoco, el tiempo de autocorrelación debe explotar exactamente a la temperatura crítica, sin embargo, se corta debido a la finitud del sistema). Por lo tanto, es mejor registrar las mediciones una vez cada 3 a 5 tiempos de autocorrelación. Tenga en cuenta que, para una semilla aleatoria dada, los efectos de autocorrelación pueden generar un error sistemático, mientras que el error estadístico debido al efecto de autocorrelación se puede estimar utilizando diferentes semillas aleatorias. Busque en este libro: http://www.amazon.com/Monte-Carlo-Methods-Statistical-Physics/dp/0198517971

  3. Mi conjetura es que, siempre que esté a 10-20 espines del límite cuando realice la medición de la corrección de espín-espín, no verá los efectos del límite. Esto es solo una suposición. ¿Por qué no pones condiciones de contorno periódicas y ves? El vecino de i,j es (i+1)%L, (j+1)%L y 3 otros como este.

¿Es fácil ajustar la temperatura para que esté en T C ?
Para ising, conocemos el valor exacto (2.27...). En general, para las transiciones de segundo orden, la temperatura crítica se estima mejor trazando el cumulante magnético de cuarto orden y notando el punto de cruce.
Primero, gracias por la entrada Srivatsan. Para responder a sus preguntas: 1. Sí, la temperatura se establece en el valor crítico: empiezo con una red totalmente aleatoria (por lo que T = ) y luego lo 'enfrio' usando la evolución del algoritmo de Wolff con probabilidad de agregar el estado vecino de 1 mi 2 k C donde Kc es un valor lo suficientemente cercano a 1 T C para que la longitud de correlación en la criticidad sea mayor que el tamaño de mi sistema (estoy usando redes entre 512x512 y 4096x4096, por lo que la precisión de 6-7 dígitos en el acoplamiento).
2. Gracias por la referencia, en realidad ya había estudiado este libro. De hecho, muestran cómo el tiempo de autocorrelación aumenta drásticamente en el punto crítico, especialmente con una "dinámica" de giro único (es decir, Metrópolis). Es por eso que me gustaría hacer que funcione con el algoritmo de Wolff, logré calcular el peso correcto con el algoritmo de Metropolis (y las condiciones de contorno libres), pero tomó mucho tiempo de cálculo y un método más efectivo se volverá realmente útil en los próximos pasos de mi proyecto.
3. Desde mi primera publicación, puse el sistema en condiciones de contorno periódicas y encontré un nuevo valor para el peso conforme del operador de espín: 0.05. Todavía no es el correcto (0.125) pero está en el "otro lado" ya que las condiciones de contorno fijas me dieron 0.4. Al menos me dice que las condiciones de contorno influyen en esta derivación, por lo que podemos ir por el camino correcto, pero tampoco estoy tan seguro de este cálculo, ya que solo se realizó sobre 5000 muestras. Verificaré eso y tal vez te responda más tarde. De nuevo, agradezco tu aporte ;)
Recién ahora miro sus comentarios. Debo decir que los datos que das (beta=0.05 con pbc y beta=0.4 con open bc) serían tomados como artefactos numéricos. Supongo que esto se debe a los tamaños de sistema muy grandes que usa y, por lo tanto, es posible que el equilibrio no haya ocurrido en absoluto. Le sugiero que use tamaños como 64x64.
Un comentario sobre los dos últimos puntos: el algoritmo de Wolff no sufre de ralentización crítica en absoluto. Si define un paso de Monte Carlo como una cantidad de volteos de clúster que, en promedio, voltea cada giro una vez, entonces Wolff se las arregla para descorrelar en unos pocos pasos de Monte Carlo, independientemente del tamaño del sistema. Para mantener pequeña la influencia de los límites, recomendaría optar por el sistema más grande que sea numéricamente factible.
¿Por qué no implementar un algoritmo de simulación perfecto ? Entonces tendrá la garantía de que cada muestra se extrae sin ningún error estadístico de la distribución de Gibbs. Es posible hacer eso también con el algoritmo de clúster. El que implementé (hace años, para hacer ilustraciones para un curso que estaba dando) estaba basado en este papel: M. Huber. Una cadena delimitadora para Swendsen-Wang. Random Structures Algorithms, 22(1):43–59, 2003. Es un poco difícil de implementar correctamente, pero es muy eficiente y obtiene todos los beneficios de una simulación perfecta.