Nested for loops in python for Ising Model

Question

Im working on statistical mechanics currently, and trying to apply some programming to it since they fit so well together! Im working on finding the partition function for a finite number of particles. However..the partition function is defined as a sum of a sum! I guess we could write this as a list of a list, so we would use nested for-loops, but i just cant quite figure out the correct way of writing it.

Z=\sum_{s_1}^{s_N}e^(s_1s_2+...+s_(N-1)s_N) is the partition function.

the possible values of s_i are -1,+1.

Effectively the ising model(1D) is a chain with N points on it and each point can have s_i=-1 or +1. The energy of the system depends on the values of s_i, and each possible combination is called a state. the total sum of these states is called Z, the partition fucntion.

So for a chain of length N=5(hence 2^5=32 possible states) how would i calculate this Z? I dont really have any code to show, but i know from the formula the result should be something like e^(+1+1+1+1+1)+e^(-1+1+1+1+1)+...+e^(-1-1-1-1-1). The question is..how on earth do I go about doing that? Ive generate the set of possible states:

import itertools
counting=0
for state in itertools.product([1,-1],repeat=5):
    print(state)
    counting+=1
print('the total possible number of states is',counting).

but how can i use this to get to a value for Z?

user_7867968 · Answer

I'd use a function to calculate the sum for each state, then do the overall sum afterwards:

import itertools
from math import exp

def each_state(products):
     for state in products:
         yield sum(state)


Z = sum(exp(x) for x in each_state(itertools.product([1,-1],repeat=5)))

The benefit of this approach is that it is in keeping with the spirit of itertools: to not aggregate everything into memory at once. So while a numpy solution might be faster, say you wanted to calculate Z for many states, a numpy implementation would start to hit memory issues whereas the generator expression will not:

from itertools import product
import numpy as np
from math import exp

# this will yield a single number, and product will yield
# each state one at a time, never aggregating the
# full set of objects into memory (even though it might seem slow)
x = sum(exp(sum(x)) for x in product([1,-1], repeat=500))


# On my 16GB MacBook, this process will be killed because
# we collect all of the states into memory
x = np.array(list(product([1, -1], repeat=500))
[1]    7743 killed     python

The general rule of thumb is that list(giant_iterable) runs out of space whereas for item in giant_iterable will run out of time

user_7469434 · Answer

Based on your description of the problem, you can calculate it using numpy as follows:

import itertools
import numpy as np

states = np.array([state for state in itertools.product([1,-1], repeat=5)])
print("There are %d states" % states.shape[0])  # 32 states

# calculate the sum for each state
sum_over_each_state = np.sum(states, axis=1)  
print(sum_over_each_state)

# calculate e^(sum(state)) for each state
exp_of_all_states = np.exp(sum_over_each_state)
print(exp_of_all_states)

# sum up all exponentials
Z = np.sum(exp_of_all_states)
print("Z:", Z)

This gives Z = 279.96.

Nested for loops in python for Ising Model

2 Answers