
clear all

% value of transition cost
beta = 1.0;

% Size of image = (NxN) 
N = 16;

% number of iterations = M
M = 100;

% number of displayed frames = M0
M0 = 10;
samp_rate = ceil(M/M0);

color = [0,0,0;1,1,1];
colormap(color);
%M_data = moviein(M0);

%Set initial condition of random pixel values
x = rand(N)>0.5;
image(x+1)
pause(1)
%M_data(:,1) = getframe;

for n=2:M
  for i=1:N
  for j=1:N
    ip= rem(N+i,N)+1;
    im= rem(N+i-2,N)+1;
    jp= rem(N+j,N)+1;
    jm= rem(N+j-2,N)+1;

    v = (1~=x(ip,j)) + (1~=x(i,jp)) + (1~=x(im,j)) +(1~=x(i,jm));

    prob = 1/( 1+exp(2*beta*(v-2)) );
    x(i,j) = rand(1)<prob;
  end
  end
  if rem(n,samp_rate)==0
    image(x+1);
    S = sprintf('Ising model: Beta = %f, Iteration = %d',beta,n);
    title(S)
    pause(1)
%    M_data(:,n) = getframe;
  end
end

%movie(M_data,10,6)
