温度调到 2.269 磁铁突然失效：用 Monte Carlo 跑一遍 2D 伊辛模型相变

把温度滑杆从 1.4 拖到 2.3，屏幕上那片黑白格子会在几秒内从「几乎连成一整块」碎成「雪花」。这不是做出来的动画，是二维伊辛模型（2D Ising model）自己的相变，理论临界温度 T_c = 2 / ln(1 + √2) ≈ 2.2692。

实测的数字先摆在这里：T = 1.40 时 |M| = 0.9915，Onsager 精确解给的是 0.9914，两者差 0.02%；T = 2.00 时 0.9106 对 0.9113（差 0.08%）；温度过 2.2692 之后，T = 3.00 时实测只剩 0.0424，而理论已经精确归零。整个实验台是一个 index.html，零外部依赖，双击就能跑；页面里的 11 条自检全部通过，本地 node verify.js 可以复跑。

一、先看结果：2.2692 这条线，页面上自己画出来了

点一下「扫一遍温度」，页面会从 1.40 一路扫到 3.60，每 0.05 取一个温度点，先平衡 150 个 sweep 再采样 150 个 sweep，把实测的 |M| 画成红点，把 Onsager 精确解画成蓝线，再在 T_c = 2.2692 处拉一条黄色竖虚线。

关键就一句话：红点在低温段紧贴蓝线，跨过 2.2692 之后一起躺平到 0。把曲线上的点抠出来做成表，误差最大的地方正好在临界点附近：

低温段（1.40–2.20）误差都在 1.3% 以内，说明 Metropolis 采样和 Onsager 公式对得上；越靠近 T_c 误差越大（2.25 处 5.06%），这不是公式写错，而是有限尺寸效应的正常表现——64 × 64 的点阵在临界点附近关联长度已经接近整块板子，单次扫描根本平均不掉涨落。

二、为什么值得看：磁铁到底在临界点上发生了什么

多数人第一次听到「伊辛模型」是在统计物理课本里，配一张黑板上推出来的配分函数，看完还是不知道磁铁为什么会「突然」失磁。这个应用把课本里那条结论变成可以拖的滑杆：

• T 远低于 T_c：绝大多数自旋朝向一致，|M| 接近 1，画面上是一片连着的磁畴，翻转接受率只有十几个百分点。
• T 无限接近 T_c：磁畴被撕成大大小小的碎片，靠一次扫描测出来的 |M| 会随着采样窗口上下跳。
• T 高于 T_c：热涨落彻底压过交换作用，接受率冲到 40% 以上，|M| 掉到 0.05 以下，画面是一片均匀的噪点。

差值也很有意思：T = 2.10 时页面显示 E/N = -1.5439、|M| = 0.7578、翻转接受率 15.1%；T = 2.80 时 E/N 抬到 -0.9551、|M| = 0.0557、接受率 41.1%。能量密度往 0 靠、磁化往 0 掉、接受率往上走——三个读数一起变，比只背一句「临界温度之上没有自发磁化」直观得多。

三、准备环境：打开码道 Web，把需求说清楚

码道有三种使用方式：WebUI（浏览器对话）、TUI（终端命令行）和桌面 IDE（IDE 插件）。本文用的是 WebUI 版，浏览器打开码道 Web 版：https://devcloud.cn-north-4.huaweicloud.com/chat?source=dmzntgwsf&sourcead=dmzntgwwbwz，登录后就能在对话窗口里输入需求，不需要装软件。

四、任务描述：一句话把物理约束写全

在对话框里输入下面这段需求（和我在码道里实际发出去的原文一字不差）：

用单个 index.html（不要外部依赖、不要联网、不要构建工具，双击就能跑）做一个「二维伊辛模型相变实验室」网页应用，主题是 2D Ising model 的铁磁相变。

物理部分：
- 正方点阵，周期性边界条件；每个格点自旋 s = ±1；哈密顿量 H = -J Σ s_i s_j，取 J = 1、外场 h = 0。
- 用 Metropolis 单自旋翻转做 Monte Carlo 采样；一次 sweep 遍历全部格点；随机数用可复现的种子（默认 12345）。
- 点阵边长可选 32 / 64 / 96 / 128。

界面控件：
- 「温度 T」滑杆 1.00~4.00，步长 0.01，默认 2.00；旁边配一个能输入数字的输入框。
- 「点阵边长」下拉框。
- 「开始 / 暂停」按钮：以每秒若干 sweep 的速度实时演化并刷新画面。
- 「重置」按钮：回到随机初态。
- 「扫一遍温度」按钮：跑一次温度扫描测量，T 从 1.40 到 3.60、步长 0.05，每个温度先平衡再采样，记录磁化强度 |M| 并画成曲线。
- 「自检」按钮：在页面里跑一遍内置断言，把每条结果显示出来。

可视化：
- 主画布：把自旋画成格子，向上自旋用浅色、向下自旋用深色；画面上要能看出磁畴。
- 曲线面板：画 |M| vs T 的折线，一条是本次实测点，一条是 Onsager 精确解 |M| = (1 - sinh(2/T)^-4)^(1/8)（T < Tc 时有值，T ≥ Tc 时取 0），并在 T = 2.2692 处画一条竖直虚线。
- 数值面板：实时显示温度 T、能量密度 E/N、磁化强度 |M|、翻转接受率、已跑的 sweep 数。

「实测 vs 理论」并排面板：
- 同屏给出本次实测 |M|、Onsager 理论 |M|、相对误差百分比。
- 旁边注明 Tc = 2 / ln(1 + sqrt(2)) ≈ 2.2692。

自检面板（页面加载后自动跑一次，点按钮可重跑），至少 10 条断言，必须是真实计算而不是直接打印常量：
1. 全同向点阵的能量密度应等于 -2J，即 -2.0
2. 全同向点阵的 |M| 应等于 1
3. 棋盘格点阵的能量密度应等于 +2J，即 +2.0
4. Onsager 公式在 T = 1.0 时 |M| 应大于 0.99
5. Onsager 公式在 T = Tc 时应约等于 0
6. Onsager 公式在 T = 3.0 时应等于 0
7. 周期边界索引回绕正确（-1 变成 N-1，N 变成 0）
8. 同一个初始种子跑两次应得到完全相同的能量序列
9. 打到 200 次 sweep 后，T = 1.5 的 |M| 应明显大于 T = 3.0 的 |M|
10. Metropolis 接受率应落在 0 到 1 之间
11. 温度扫描结束时曲线点数应等于预期点数

失败路径与错误提示：
- 温度输入框里输入非法值（空、非数字、超出 1.00~4.00）时，弹出一条明确的中文提示，并保留上一次的有效温度，提示里写清「已保留上一次的有效值」。
- 点阵边长只接受给定的几档，其它值给提示。

性能：
- 128×128 点阵单次 sweep 要在 20 毫秒量级完成（用 Int8Array 存自旋，避免每帧大量对象分配）。
- 面板刷新要节流，不要每个 sweep 都重排 DOM。

完成后打开预览，让我直接看到运行效果。

写需求时我特意加了三类约束：物理量要写清定义（周期性边界、J = 1、h = 0），性能要给具体指标（128 × 128 单次 sweep 20 ms 量级），错误路径要写清文案格式（「已保留上一次的有效值」）。这三类写清楚了，验收时才有东西可对。

五、码道实际做了什么

这一轮任务跑了 11 分 57 秒，产出 5 个产物，交付物就是 /workspace/index.html，单文件 29.8 KB：

• 物理内核：IsingModel 用 Int8Array 存自旋，周期性边界用取模回绕，Metropolis 单自旋翻转；
• 随机数：自己实现 mulberry32，固定种子 12345，保证「同一种子跑两次结果完全一致」；
• 精确解：onsagerM(T) 按 |M| = (1 - sinh(2/T)^-4)^(1/8) 给出理论曲线，T_c 直接写成 2 / Math.log(1 + Math.sqrt(2))；
• 可视化：主画布画自旋（浅色向上、深色向下），第二块画布画 |M| vs T 的实测点和理论线；
• 自检：11 条断言全在页面里跑，结果直接渲染在「自检面板」上。

它没有停在「功能能跑」：一边写一边自己跑 node 脚本验证，还把两处真问题修掉了（见第七节和第八节）。

六、代码里最关键的几段

Metropolis 那一步是整个应用的性能瓶颈，页面里的写法是把 exp(-βΔE) 预先算好——因为 ΔE 只可能是 0、4、8 三种值：

IsingModel.prototype.sweep = function (T) {
  var beta = 1 / T;
  p4 = Math.exp(-4 * beta);
  p8 = Math.exp(-8 * beta);
  for (var i = 0; i < N; i++) {
    var sum = s[rowUp + x] + s[rowDown + x] + s[base + lx] + s[base + rx];
    var sval = s[i];
    var dE = 2 * sval * sum;          /* J = 1，外场 h = 0 */
    if (dE <= 0 || rng() < (dE === 8 ? p8 : p4)) {
      s[i] = -sval;
      e += dE;
      m -= 2 * sval;                  /* 翻转一个自旋，磁化增量的正确写法 */
      accepted++;
    }
  }
};

Onsager 精确解也没偷懒，T ≥ T_c 直接返回 0，避免了 sinh(2/T)^-4 在高温下变成负数再开 1/8 次方：

function onsagerM(T) {
  if (T < TC) {
    var sh = Math.sinh(2 / T);
    var v = 1 - Math.pow(sh, -4);
    if (v <= 0) return 0;
    return Math.pow(v, 1 / 8);
  }
  return 0;
}

七、验证：11 条自检全部通过，本地也能复跑

页面「自检」按钮背后是 11 条真算出来的断言，不是打印常量。下面这张表是从页面上抄下来的（括号里是计算值）：

同一批断言在本地也能跑，不用开浏览器：

node verify.js index.html
# 自检 11/11 通过
#   同参数对照（32×32、T=1.5、200 sweep、种子 12345）：
#     正确实现 |M| = 0.9902
#     符号写错 |M| = 6.0840  <- 超过 1，物理上不可能
#   128×128 单次 sweep 平均 4.389 ms

八、踩坑复盘：|M| 算出 6.08 的那次符号错误

这一轮最值钱的一段代码不是功能，而是码道自己抓到的 bug。它先按需求写完物理内核，跑性能基准时发现低温下的磁化强度不对劲：

/* 写错的样子：磁化增量的正负跟着 dE 的符号走 */
m -= 2 * (dE > 0 ? s[i] * -1 : s[i]);

翻转一个自旋，磁化强度的增量只跟翻转前的自旋值有关（ΔM = -2·s_old），跟这一步的能量变化 dE 无关。写成依赖 dE 的符号之后，误差会顺着 sweep 累积起来：同样是 32 × 32、T = 1.5、200 个 sweep、种子 12345，正确实现给出 |M| = 0.9902，写错的那版给出 |M| = 6.0840——磁化强度是「向上自旋占比的偏差」，物理上根本不可能超过 1。这条我自己在本地用 verify.js 复算了一遍，数字和码道对话里报的一致。

另一个坑跟物理无关，是「有限尺寸效应」被误读：64 × 64 的点阵在第一版扫描里，低温段（1.40–1.65）测出来只有 0.06–0.29，然后在 1.70 附近突然抬到 0.9 以上。原因不是模型错了，而是每个温度都从独立的随机初态开始、平衡步数不够，低温下 120 个 sweep 还没来得及长出全域磁化。改成「先高温热身、再连续升温扫描」之后，曲线才变成图 2 里那条光滑的样子。

九、亲手点三下（必须包含一次失败路径）

想自己确认一遍的话，按这三步点：

1. 看关键数值：把温度设成 2.10，点「开始」，等十几秒，右上角应该读到 |M| ≈ 0.75、接受率 15% 左右，画面上磁畴连成大片。
2. 换个状态：把温度拖到 2.80，再看 |M|，应该掉到 0.06 以下、接受率升到 40% 左右；再点「扫一遍温度」，红点会在 2.2692 附近掉头向下。
3. 走一遍失败路径：在温度输入框里输入 abc（或 5、留空后回车），页面不应该白屏，也不应该悄悄吞掉这次输入，而是弹出提示并保留上一次的有效温度：

温度输入无效：请输入 1.00 ~ 4.00 之间的数字（已保留上一次的有效值 2.80）

十、本地复现

git clone https://atomgit.com/deli007/demo_park.git
cd demo_park/codearts-ising-phase-transition
node verify.js index.html        # 复跑 11 条自检 + 复现那个符号 bug
python -m http.server 8000       # 或者直接双击 index.html，效果一样

verify.js 会把 <script id="sim-core"> 里的物理内核抽出来在 Node 里跑，不依赖浏览器；index.html 本身零外部依赖，没有 CDN、没有构建步骤。

十一、直接下载试玩

这个 Demo 已经完整开源（源代码、本文档、全部截图），托管在公开仓库，可直接下载来玩：

• Demo Park 仓库：https://atomgit.com/deli007/demo_park
• 本案例目录：https://atomgit.com/deli007/demo_park/tree/main/codearts-ising-phase-transition

打开上面的链接，进入 codearts-ising-phase-transition 文件夹，把 index.html 下载到本地，用浏览器打开即可运行，无需安装任何依赖。

十二、使用码道体会

• 物理类需求要把「可验收的量」写进需求里。这次我直接要求「实测值、理论值、相对误差同屏并排」，码道就把 Onsager 精确解接进来了；如果我只说「做一个伊辛模型演示」，回来的多半只有一张会动的格子图。
• 让它自己跑基准，比自己看代码有效。那个 |M| = 6.08 的符号错误不是靠读代码发现的，是它跑性能基准时看到「|M| 超过 1」才揪出来的。给需求加上「要有自检并显示结果」，等于给对方一根标尺。
• 性能上限要写具体数字。「128 × 128 单次 sweep 20 ms 量级」比「性能要好」有用，交付时我能在页面上直接读到 0.20 ms 的读数。
• 失败路径要写清文案格式。我写了「提示里要写清『已保留上一次的有效值』」，交付的提示原文就是可验收的。
• 有限尺寸效应这类「非 bug」要提前想好怎么解释。曲线在临界点附近抖，是模型特性不是实现错误，作者得自己分清楚，否则会把物理现象当成 bug 反复折腾。

十三、总结

一个 29.8 KB 的单文件网页，把二维伊辛模型的相变从公式变成了可以拖的滑杆：低温段实测与 Onsager 精确解差 0.02%–1.28%，跨过 T_c = 2.2692 后 |M| 掉到 0.03 附近，11 条自检全绿。想继续玩的话，可以试着把点阵开到 128 × 128 看临界点附近的涨落有多大，或者加一个零外场之外的外场 h，看相变线怎么移动。

说明：文中所有读数（自检值、扫描表、性能数字、错误提示原文）都是在本机实际运行页面和 verify.js 取到的，不是转述。