Clusters_Ising.py

import numpy as np


def find_clusters(sample):

    N = sample.shape[0]
    mask = np.ones(N**2, dtype=bool)

    done = set()
    todo = set()

    cluster_sizes = []
    while mask.sum() > 0:

        s0 = np.random.choice(np.arange(N**2)[mask])
        i0 = s0//N; j0 = s0%N

        sign = sample[i0,j0]
        todo.add(s0)
        count = 1

        while len(todo) > 0:
            s = todo.pop()
            i, j = s//N , s%N
            
            possible = {(1, 0), (-1, 0), (0, 1), (0, -1)}

            for (delta_i, delta_j) in possible:
                ii = (i+delta_i)%N
                jj = (j+delta_j)%N
                ss = N*ii + jj
                if sample[ii,jj] != sign:
                    continue 
                
                if (ss not in done) and (ss not in todo) and mask[ss]:
                    todo.add(ss)
                    count += 1 
            
            done.add(s)
            mask[s] = False

        cluster_sizes.append(count)

    return cluster_sizes