# Speeding up Module

I am working on the Ising 2D model in Mathematica and this is currently my Metropolis sampling function:

flip[B_, J_, lat_, num_, T_] :=
Module[{i, j, down, latbuf = lat , left, right, up, s, x},
Do[i = Random[Integer, {1, num}];
j =  Random[Integer, {1, num}];
If[i == num, down = 1, down = i + 1];
If[j == num, right = 1, right = j + 1];
If[i == 1, up = num, up = i - 1];
If[j == 1, left = num, left = j - 1];

s = latbuf [[down, j]] + latbuf [[up, j]] + latbuf [[i, right]] +
latbuf [[i, left]];
x = latbuf [[i, j]];

If[ (2 x (B + s J)) < 0 || Random[] < Exp[-(2 x (B + s J))/T],
latbuf [[i, j]] = -latbuf [[i, j]],
latbuf[[i, j]] = latbuf[[i, j]],
latbuf[[i, j]] = latbuf[[i, j]]], 1000];
latbuf ]


but it is ridiculously slow, especially when I want to calculate magnetization:

 num = 16;
lattice = 2 Table[Random[Integer], {num}, {num}] - 1;

magnetization =
Evaluate[Table[{x,
N @(Sum[Total[Total[flip[0, 1, lattice, num, x ]]], {i, 1, 100}]/(
50 num^2)) }, {x, 0.001, 6, 0.1} ]]


Is there any way to speed it up?

## 1 Answer

num = 16;
a = 2 RandomInteger[{0, 1}, {num, num}] - 1;
ker = CrossMatrix;
ker[[2, 2]] = 0;
B = 0;
J = 1;
T = 0.1;


In principle, this should do the same as your function flip but more efficient.

Do[
a *= (2 UnitStep[Exp[-2 a (B + ListCorrelate[ker, a, {2, 2}] J)/T] - RandomReal[{0, 1}, {num, num}]] - 1);
, {1000}]; // AbsoluteTiming // First


0.041588

This should give you an idea. I haven't tested it thoroughly, yet, because I did not understand what exactly you try to achieve.