Markov ising.py

From Werner KRAUTH

(Difference between revisions)
Jump to: navigation, search
Revision as of 22:24, 4 March 2024
Werner (Talk | contribs)
(Program)
← Previous diff
Current revision
Werner (Talk | contribs)
(Version)
Line 1: Line 1:
-This page presents the program markov_ising.py, a Markov-chain algorithm for the Ising model on an LXL square lattice in two dimensions.+This page presents the Python3 program markov_ising.py, a Markov-chain algorithm for the Ising model on an LXL square lattice in two dimensions, with periodic boundary conditions. The program uses the Metropolis algorithm.
 + 
__FORCETOC__ __FORCETOC__
=Description= =Description=
- +The program is described in my 2024 Oxford Lecture No 8, and also in my book. This version of the program estimates the energy per particle, and the specific heat.
=Program= =Program=
- +
import random, math import random, math
- +
L = 6 L = 6
N = L * L N = L * L
Line 32: Line 33:
E2_av = E2_tot / float(nsteps) E2_av = E2_tot / float(nsteps)
c_V = beta ** 2 * (E2_av - E_av ** 2) / float(N) c_V = beta ** 2 * (E2_av - E_av ** 2) / float(N)
- print(E_av / N,c_V)+ print(E_av / N, c_V)
=Version= =Version=
See history for version information. See history for version information.
-[[Category:Python]] [[Category:Honnef_2015]] [[Category:MOOC_SMAC]]+[[Category:Python]] [[Category:Oxford_2024]] [[Category:MOOC_SMAC]]

Current revision

This page presents the Python3 program markov_ising.py, a Markov-chain algorithm for the Ising model on an LXL square lattice in two dimensions, with periodic boundary conditions. The program uses the Metropolis algorithm.


Contents

Description

The program is described in my 2024 Oxford Lecture No 8, and also in my book. This version of the program estimates the energy per particle, and the specific heat.

Program

import random, math

L = 6
N = L * L
nbr = {i : ((i // L) * L + (i + 1) % L, (i + L) % N,
            (i // L) * L + (i - 1) % L, (i - L) % N) \
                                    for i in range(N)}
nsteps = 10000000
T = 1.0
beta = 1.0 / T
S = [random.choice([1, -1]) for k in range(N)]
E = -0.5 * sum(S[k] * sum(S[nn] for nn in nbr[k]) \
                                for k in range(N))
E_tot, E2_tot = 0.0, 0.0
for step in range(nsteps):
    k = random.randint(0, N - 1)
    h = sum(S[nn] for nn in nbr[k])
    Sk_old = S[k]
    delta_E = 2.0 * S[k] * h
    if random.uniform(0.0, 1.0) < math.exp(-beta * delta_E):
        S[k] *= -1
        E -= 2.0 * h * S[k]
    E_tot += E
    E2_tot += E ** 2
E_av  = E_tot / float(nsteps)
E2_av = E2_tot / float(nsteps)
c_V = beta ** 2 * (E2_av - E_av ** 2) / float(N)
print(E_av / N, c_V)

Version

See history for version information.

Personal tools