function E = energy(s,J,B);
    E = -0.5*J*sum(sum(s.*(s(:,[2:end,1])+s(:,[end,1:end-1])+s([2:end,1],:)+s([end,1:end-1],:)))) - B*sum(sum(s));
end