Génération de configurations d'état stable du modèle Ising

Question

Quelle est la manière la plus efficace de simuler des configurations d'état stationnaire du modèle Ising? Je suis juste intéressé par un grand ensemble de configurations aléatoires d'état stationnaire du modèle 1D Ising (avec des constantes de couplage homogènes). Quelques idées me sont venues à l'esprit:

1. Échantillonnage par force brute. Puisque le modèle d'Ising est exactement résoluble en 1D et 2D, on a des expressions exactes pour les probabilités de chaque état. Cependant, l'échantillonnage aléatoire sur un ensemble de$2^N$ causera probablement des problèmes de mémoire pour les petits $N$ déjà.
2. Dynamique de Monte Carlo. On pourrait exécuter les algorithmes de Monte Carlo habituels (par exemple la dynamique de Glauber) sur des états initiaux aléatoires et attendre que le système converge vers l'équilibre thermique. Cependant, cela semble inefficace lorsque vous n'êtes pas intéressé par la dynamique et que vous ne voulez que des configurations en régime permanent.
3. Utilisation de la densité d'états. On pourrait aussi d'abord échantillonner au hasard l'énergie du système, selon$P(E) \sim N(E) \exp(-\beta E)$, où $N(E)$est la densité d'états, qui est calculable (au moins numériquement). Ensuite, on génère une configuration aléatoire avec cette énergie, par exemple en utilisant un algorithme de spin flip où l'on retourne des spins simples pour augmenter / diminuer l'énergie jusqu'à ce qu'elle corresponde à l'énergie cible. Mais je ne suis pas sûr que les configurations obtenues de cette manière suivent statistiquement la distribution de Boltzmann.

Remarque: en 1D, il existe également une expression exacte de la densité d'états d'Ising, $g(E(k)) = 2 \binom{N-1}{k}$ avec $E(k) = -N + 2k + 1$. Voir cette autre question: Ising modèle de densité d'états .

Avez-vous des idées sur la meilleure façon d'aborder cela?

Answer 1

Pour le modèle unidimensionnel, le moyen le plus efficace, de loin, de simuler le modèle d'Ising est d'utiliser une chaîne de Markov sur $\{-1,1\}$, générant un tour à la fois, conditionnellement aux valeurs prises par les tours précédents. Notez également que de cette manière, vous échantillonnez exactement à partir de la distribution de Gibbs, sans approximation (contrairement à une approche de Monte Carlo).

Pour simplifier, permettez-moi de considérer le modèle avec condition aux limites libre, c'est-à-dire le modèle avec hamiltonien $$ \beta\mathcal{H} = - \beta\sum_{i=2}^N \sigma_{i-1}\sigma_i . $$(Vous pouvez également ajouter un champ magnétique, mais je ne le ferai pas ici pour simplifier l'exposition). Ensuite,$\sigma_1$ est égal à $+1$ ou $-1$ avec probabilité $\tfrac12$par symétrie. De plus, pour tout$k\geq 2$, $$ \mathrm{Prob}(\sigma_k=\sigma_{k-1} \,|\, \sigma_1, \dots, \sigma_{k-1}) = \mathrm{Prob}(\sigma_k=\sigma_{k-1}) = \frac{e^{\beta}}{e^{\beta} + e^{-\beta}} = \frac{1}{1+e^{-2\beta}}. $$ Appelons cette probabilité $p$.

Résumer:

• Vous échantillonnez $\sigma_1$: il est $+1$ avec probabilité $\tfrac12$ et $-1$ avec probabilité $\tfrac12$.
• Donné $\sigma_1$, vous définissez $\sigma_2 = \sigma_1$ avec probabilité $p$ et $\sigma_2 = -\sigma_1$ avec probabilité $1-p$.
• Donné $\sigma_2$, vous définissez $\sigma_3 = \sigma_2$ avec probabilité $p$ et $\sigma_3 = -\sigma_2$ avec probabilité $1-p$.
• etc...

C'est très facile à mettre en œuvre et extrêmement rapide (bien sûr, calculer $p=1/(1+e^{-2\beta})$juste une fois). Ensuite, la plupart du temps est pris par la génération de nombres pseudo-aléatoires. De cette façon, vous pouvez simuler des chaînes de longueur arbitraire sans aucun problème.

(Voir aussi cette réponse pour un autre point de vue sur la relation entre les modèles unidimensionnels et les chaînes de Markov.)

---

Explication de la formule pour $p$.

La façon la plus simple de voir pourquoi la formule pour $p$donnée ci-dessus est en utilisant soit le cluster aléatoire soit les représentations à haute température du modèle d'Ising, si vous les connaissez (elles sont décrites, par exemple, dans les sections 3.7.3 et 3.10.6 de ce livre ) .

Si vous n'êtes pas familier avec ces représentations, permettez-moi d'essayer de fournir un argument direct.

Laisser $s_1,\dots,s_N \in \{-1,1\}$ et écris $s=(s_1,\dots,s_{k-1},s_k,\dots,s_N)$ et $s'=(s_1,\dots,s_{k-1},-s_k,\dots,-s_N)$ (c'est-à-dire la configuration $s'$ est obtenu à partir de la configuration $s$ en retournant les tours à $k, k+1, \dots N$).

Maintenant, $$ \frac{{\rm Prob}(\sigma = s)}{{\rm Prob}(\sigma = s')} = \frac{\exp\bigl( -\beta \mathcal{H}(s) \bigr)}{\exp\bigl( -\beta\mathcal{H}(s') \bigr)} = \exp(2\beta\, s_{k-1}s_{k}). $$ En particulier, $$ \frac{{\rm Prob}(\sigma_k=\sigma_{k-1})}{{\rm Prob}(\sigma_k = -\sigma_{k-1})} = \exp(2\beta). $$ Mais cela implique que $$ {\rm Prob}(\sigma_k=\sigma_{k-1}) = e^{2\beta}\, {\rm Prob}(\sigma_k = -\sigma_{k-1}) = e^{2\beta} \bigl( 1 - {\rm Prob}(\sigma_k = \sigma_{k-1}) \bigr), $$ et donc $$ (1+e^{2\beta})\, {\rm Prob}(\sigma_k=\sigma_{k-1}) = e^{2\beta}, $$ à partir de laquelle la formule pour $p$ suit immédiatement.