J = Table[If[Random[] ≤ 0.5, 0, 1], {i, 1, 8}, {j, 1, 8}]; 
    MatrixForm[J]

inimigos = 
      Function[{i, j}, 
        Module[{diagonal1, diagonal2, linha, coluna}, 
          expelir = 
            Function[l, 
              1 ≤ l[[1]] ≤ 8 && 
                1 ≤ l[[2]] ≤ 
                  8 && (l ≠ {i, j})]; 
          diagonal1 = Table[{i - x, j - x}, {x, -7, 7}]; 
          diagonal2 = Table[{i - x, j + x}, {x, -7, 7}]; 
          linha = Table[{i, k}, {k, 1, 8}]; 
          coluna = Table[{k, j}, {k, 1, 8}]; 
          Select[Join[diagonal1, diagonal2, linha, coluna], 
            expelir]]]

inimigos[1, 1]

mostraSistema = 
      Function[M, 
        FrameBox[
            GridBox[Table[
                If[M[[i, j]] == 
                    0, " ", "Ψ"], {i, 
                  1, Dimensions[M][[1]]}, {j, 
                  1, Dimensions[M][[2]]}], RowLines -> True, 
              ColumnLines -> True]] // DisplayForm]

mostraSistema[J]

estadoInimigos = 
      Function[{M, i, j}, 
        Module[{estado}, 
          estado = Function[ponto, M[[ponto[[1]], ponto[[2]]]]]; 
          Map[estado, inimigos[i, j]]]]

energiaDeUm = 
      Function[{T, b, p, i, j}, 
        Module[{numerario}, 
          numerario = 
            Function[p, T[[p[[1]], p[[2]]]]]; -1/
              2*(b*T[[i, j]] - 
                Apply[Plus, 
                  p*T[[i, j]] 
                    Map[numerario, inimigos[i, j]]])]]

energia = 
      Function[{M, b, p}, 
        Apply[Plus, 
          Flatten[Table[
              energiaDeUm[M, b, p, i, j], {i, 
                1, Dimensions[M][[1]]}, {j, 
                1, Dimensions[M][[2]]}]]]]

probabilidade = 
      Function[{delta, T}, 1/(1 + Exp[1. *delta/T])]

reforco = 
      Function[{T, b, p, i, j}, 
        Module[{numerario}, 
          numerario = Function[p, T[[p[[1]], p[[2]]]]]; 
          If[T[[i, j]] == 0, 
            0, -1/
                2*(b*T[[i, j]] - 
                  2*Apply[Plus, 
                      p*T[[i, j]] 
                        Map[numerario, inimigos[i, j]]])]]]

sistemaBoltzmann = 
      Function[{M, b, p, T}, 
        Module[{i, j, NS, delta}, 
          i = Random[Integer, {1, Dimensions[M][[1]]}]; 
          j = Random[Integer, {1, Dimensions[M][[2]]}]; NS = M; 
          NS[[i, j]] = If[M[[i, j]] == 0, 1, 0]; 
          delta = reforco[NS, b, p, i, j] - reforco[M, b, p, i, j]; 
          If[Random[] < probabilidade[delta, T], NS, M]]]

boltzmannAux = 
      Function[{M, T}, sistemaBoltzmann[M, 1, 5, T]]

Anealing = 
      Function[{T0, a, b, p, n}, 
        Module[{M, T}, 
          boltzmannAux = Function[{M, T}, sistemaBoltzmann[M, b, p, T]]; 
          M = Table[If[Random[] ≤ 0.5, 0, 1], {i, 8}, {j, 8}]; 
          T = T0; While[T ≥ 1, 
            M = Last[FoldList[boltzmannAux, M, Table[T, {i, 1, n}]]]; 
            Print[mostraSistema[M]]; 
            Print["Temperatura=", T, " Energias=", energia[M, b, p]]; 
            T = a*T]]]

Anealing[20, 0.9, 175, 145, 64]