我写过一个用mathematica的程序,比较简单.要的话可以给你.
再问： 好的，麻烦上传一下
再答： Clear[n, JT, s, M, Z, A, de]; n = 6; JT = .5; s = Table[1, {n}, {n}];(*选取初始configuration为6*6的自旋全部向上*) M = n^2.;(*M标记update过程中的总的磁化强度的绝对值，初始为36*) Z = 0;(*标记随机求和所进行的次数，初始为0*) data = Table[(A = 0; k = 10000; Do[ Module[{k1, k2, \[Xi]}, k1 = RandomInteger[{1, 6}]; k2 = RandomInteger[{1, 6}]; \[Xi] = .1 Random[]; de = 0; If[k1 > 1, de = de + s[[k1 - 1, k2]], de]; If[k1 < 6, de = de + s[[k1 + 1, k2]], de]; If[k2 > 1, de = de + s[[k1, k2 - 1]], de]; If[k2 < 6, de = de + s[[k1, k2 + 1]], de]; R = Exp[-2 JT de s[[k1, k2]]]; If[R > \[Xi], s[[k1, k2]] = -s[[k1, k2]], s[[k1, k2]]]; If[R > \[Xi], M = M + 2 s[[k1, k2]], M]; A = A + Abs[M]; ], {k}]; N[A/k]), {30}] 做图： Show[{Plot[Mean[data], {x, 0, 30}, PlotStyle -> Red], ListLinePlot[data]}, PlotRange -> {4, 7}, Frame -> True]
