Production 設定Integrator, thermostat, constraints
timestep、控溫、約束：怎麼換到速度而不換掉物理。
01一句話只有一句
每個 production 參數都是一次換匯：物理正確性換牆鐘時間。
02Why should I care?我為什麼要讀這頁
因為這一頁的參數決定了你能跑幾條 replica，而 replica 數才是你整套可靠性論證的地基。
你的 5F9R 完整三元系統是 535,265 個原子，4 條 replica、每條約 87 ns。timestep 從 2 fs 換到 4 fs，同一批 GPU 時數大約買到兩倍的軌跡 —— 你可以把它花在長度上，也可以花在 replica 數上。這件事不是效能問題，是統計力的預算問題：在 2 條 replica 的時候，逐殘基 top-8 接觸殘基看起來有 6/8 共享，非常漂亮；跑到 4 條掉到 3/8（詳見取樣與收斂）。負不負擔得起第三條與第四條，決定了你會不會在 n=2 就停手。 production 設定買到的不是速度，是「有沒有本錢做對照」。
用你熟的世界對照一下，這裡有一個 MD 特有的陷阱值得先講：
在 VASP 裡，ENCUT 和 k-mesh 是純粹的收斂旋鈕 —— 調大只是變慢，答案往同一個極限靠。production MD 的參數沒有一個是這種旋鈕。它們分成三類，而你必須知道手上轉的是哪一類：
- 真的是收斂旋鈕：timestep。 往小調，積分誤差往零走。這是唯一一個可以照 DFT 直覺處理的參數。
- 看起來像收斂旋鈕、其實是參數化的一部分：real-space cutoff。 DFT 直覺會說「cutoff 調大一定更準」，在這裡是錯的。力場的 Lennard-Jones 參數是在某個 cutoff 底下配出來的（amber14 這一系的慣例是 0.8 到 1.0 nm 搭 PME），把它推到 1.2 nm 不會更接近真值，只會離開參數化的條件。它是敏感度測試的旋鈕。
- 完全不是旋鈕，是定義：thermostat 與 barostat。 選了錯的 thermostat，你跑一微秒也不會收斂到正確答案，因為你從頭到尾在取樣另一個系綜。
DFT 世界裡沒有一個 flag 是「跑越久、錯得越穩定」的，MD 裡有。這就是為什麼這一頁不能跳過。
03What這是什麼
Production 就是平衡完之後、拿掉所有 restraint、開始把數據寫進硬碟的那一段。 它由八個決定組成，每一個都可以單獨毀掉結果：
- Integrator —— 怎麼把 F = ma 往前推一步。
- Timestep Δt —— 一步多長。
- Constraints —— 哪些鍵長／鍵角被鎖成剛體。演算法的名字依套件而異：Amber 走 SHAKE／RATTLE、GROMACS 走 LINCS、水幾乎一律走解析解 SETTLE，而 OpenMM 內部用的是 CCMA。你在 OpenMM 文件裡找不到 SHAKE 這個字，不是你漏設了什麼。
- Hydrogen mass repartitioning（HMR） —— 把質量從重原子搬到氫上。
- Thermostat —— 怎麼控溫，以及因此在取樣哪個系綜。
- Barostat —— 要不要控壓、怎麼控。
- 長程靜電 —— PME 加上 real-space cutoff。
- 輸出節奏 —— 多久寫一幀、多久寫一次 checkpoint。
timestep 的物理上限由系統裡最快的運動決定。 這條規則很好記：Δt 大約要小於最快週期的十分之一，否則數值積分會抓不到那個振盪。把波數換成週期就看得出各層級差在哪：
- X–H 伸縮 ≈ 3000 cm⁻¹ → 週期 ≈ 11 fs → 允許 Δt ≈ 1 fs
- 含氫的角彎曲 ≈ 1400 cm⁻¹ → 週期 ≈ 24 fs → 允許 Δt ≈ 2 fs
- 重原子伸縮 ≈ 1000 cm⁻¹ → 週期 ≈ 33 fs → 允許 Δt ≈ 3 fs
（換算就是 T = 1/(c·ν̃)。這幾個波數是官能基的典型量級，不是你系統的實測值。）
所以標準做法是一層一層把最快的運動拿掉：先用 constraints 把 X–H 鍵長鎖死，最快的東西就變成含氫的角彎曲，Δt 進到 2 fs；再用 HMR 把氫變重，角彎曲也慢下來，Δt 進到 4 fs。
但這個十分之一法則是給決定性積分器的，而你用的不是。 OpenMM 使用者手冊把話講得很明白：同樣是 constraints=HBonds，Verlet 大約 2 fs，Langevin 大約 4 fs。差別在於隨機 thermostat 對構形分布比對軌跡本身寬容得多 —— 而這件事有代價，代價寫在 When NOT 第 3 條：它同時也讓 timestep 太大時不會崩潰給你看。
HMR 到底在做什麼。 振動頻率 ω = √(k/μ)。含氫運動的約化質量幾乎就是氫本身，所以把氫加重就直接把頻率壓下來。把氫從 1.008 amu 加到 3.024 amu（Amber 生態的慣例值 —— ParmEd 的 HMassRepartition 預設就是 3.024 Da），如果只看氫的質量，最快的含氫模式頻率粗估降到 1/√3 ≈ 0.58 倍，週期拉長約 1.7 倍。這個 1.7 倍是樂觀上界：多出來的質量是從相連的重原子扣掉的（系統總質量不變），重原子變輕會吃掉一部分約化質量的增幅，對 C–H 實際算出來大約是 1.6 倍。取哪個數都好，量級足以支撐 2 fs 進到 4 fs。關鍵是：位能面完全沒動，動到的只有質量矩陣。
改 timestep 改的是精度，改 thermostat 改的是你在算什麼。 前者是收斂問題，後者是定義問題。這兩件事必須分開處理，混在一起就是本頁所有錯誤的來源。
04Why為什麼重要
因為你要報的量，全都是「漲落」，而系綜決定漲落。
看一下你的可靠性表在量什麼：BSA、接觸數、鹽橋數、RMSF。這些沒有一個是平均位置，全部是分布的寬度或分布的一階矩。統計力學裡，漲落的大小由系綜直接決定 —— 熱容來自能量漲落、壓縮率來自體積漲落、RMSF 來自座標漲落。所以一個「溫度看起來對、但漲落錯」的 thermostat，會讓你的每一個數字都偏移，而且溫度計上完全看不出來。
Berendsen thermostat 是最經典的例子。它用速度縮放把溫度拉回設定值，效率很高，但它產生的不是 canonical 系綜：動能漲落被系統性壓窄。更糟的是所謂 flying ice cube —— 動能會從高頻模式（鍵伸縮）漏到低頻模式（整體質心平移），最後你得到一塊在水盒裡飛的冰。文獻上的建議已經很直白：production 不要再用 Berendsen，改用真的產生 canonical 分布的方法。Berendsen barostat 也一樣，它給你正確的平均密度，但不給你正確的體積漲落。
Constraints 的代價比較好講。 把 X–H 鍵長鎖死，等於把那些自由度從系綜裡拿掉。你會失去什麼？X–H 的振動光譜與零點能相關的東西。你在乎嗎？不在乎 —— 古典力場本來就算不對 X–H 振動（那是 DFT 的工作），而你要的是構形取樣。這是一筆划算的交易，而且划算的理由必須寫得出來，因為 reviewer 會問。
HMR 的代價就沒那麼單純了，而且對你的專案特別敏感。 HMR 不動位能面，所以平衡態的構形分布原則上不變 —— 熱力學量守住了。但質量變了，動力學就變了：擴散係數、自相關時間、任何速率常數都會被重新縮放。你這個專案的核心主張是「Cas9 的專一性是 kinetics 不是平衡 ΔG」。如果你哪天要從軌跡裡讀出時間尺度，HMR 就是必須先處理掉的第一個 confounder。目前的用法（讀構形與介面接觸的平衡性質）是安全的，但這條界線要自己記住。
05How怎麼做
這個專案實際跑的設定長這樣（拼自建系統與 production 兩支腳本，不是逐字的單一片段）：
# 建 system
system = ff.createSystem(topology,
nonbondedMethod=PME, nonbondedCutoff=1.0*nanometer,
constraints=HBonds, rigidWater=True,
hydrogenMass=1.5*amu)
# 掛 barostat、開 integrator（每條 replica 換 SEED）
baro = MonteCarloBarostat(1*bar, 300*kelvin, 25); baro.setRandomNumberSeed(SEED)
system.addForce(baro)
integ = LangevinMiddleIntegrator(300*kelvin, 1/picosecond, 2*femtosecond)
integ.setRandomNumberSeed(SEED)
...
integ.setStepSize(4*femtosecond) # 平衡走 2 fs，只有密度平衡與 production 用 4 fs
逐項拆開：
LangevinMiddleIntegrator—— OpenMM 的 LFMiddle 離散化。官方文件說它與 BAOAB 產生相同的軌跡，差別只在 LFMiddle 回傳半步（leapfrog）速度、BAOAB 回傳整步速度，而前者對熱系綜的取樣更準；同一份文件也明說它在構形取樣上比舊的LangevinIntegrator更準。它同時是 integrator 也是 thermostat，而且是產生正確 canonical 分布的那一類。用它，不要用 Berendsen。- 摩擦係數 1/ps —— 弱耦合。太大會把動力學拖慢成擴散、太小會控溫無力。1/ps 是蛋白質模擬的慣用值。注意這個數字會影響動力學觀測量，要談 kinetics 就得記得它。
constraints=HBonds+rigidWater=True—— 鎖住所有含氫鍵長，水用剛體（水走 SETTLE 那類解析解，比迭代解快且穩）。hydrogenMass=1.5*amu+ 4 fs —— 這一項的定位比直覺微妙，見 When NOT 第 4 條。先講結論：它不是超速，但它落在兩個生態的慣例中間，所以要準備好被問。MonteCarloBarostat(1 bar, 300 K, 25)—— NPT，每 25 步嘗試一次體積移動（25 就是 OpenMM 建構子的預設值，不是這個專案特別挑的數字，方法段照實寫成預設即可）。Monte Carlo 型 barostat 用 Metropolis 判準接受或拒絕體積變動，給出正確的 isothermal-isobaric 體積分布；Berendsen barostat 則是把體積漲落壓窄，平均密度對、分布不對。- PME + 1.0 nm cutoff —— 標準組合。cutoff 只切 real-space 那一段，倒空間由 PME 補回來，所以 1.0 nm 不代表你在 1.0 nm 之外就不算靜電。
平衡到 production 的階梯（Skeens 式分段釋放，這段是專案實際用的排程）：
重點在restraint 是逐步放掉的，而且 4 fs 只在最後才打開。加熱與早期 NPT 全程 2 fs，等到密度平衡那 5 ns 才切到 4 fs，production 才維持 4 fs。這個順序不是儀式，是因為系統最不穩定的時候（剛加熱、還有 restraint、密度還沒收）恰好最禁不起大 timestep。
GPU 上的實務：
- 精度用
mixed。 單精度算力、雙精度累加。全雙精度在消費級卡上是災難級的慢，全單精度會有可見的能量漂移。 - 每條 replica 換種子，而且要換全部三個。 你的腳本裡
SEED同時餵給 integrator、barostat、以及setVelocitiesToTemperature。嚴格講，光換 integrator 種子就已經是不同的 Langevin 噪音實現，軌跡本來就會分開；三個一起換的理由不是物理，是方法段可以一句話講完、reviewer 沒有縫可以鑽。真正限制獨立性的從來不是種子，是四條共用同一個起始座標，那件事換再多種子也解決不了。 - 一定要 checkpoint。 你的排程是每 100 ps 寫一次，腳本偵測到
prod_<seed>.chk就續跑。共享叢集上同時可跑的 job 數與單次 walltime 通常都有上限，沒有續跑能力就等於每次被砍都從頭來。 - 輸出節奏要先算硬碟。 DCD 每 20 ps、state log 每 10 ps。這個系統的 DCD 就已經 22 GB 等級，而家目錄的配額通常遠小於軌跡體積 —— 軌跡放大容量檔案系統、分析在叢集跑完只傳小 CSV 回來。存太密不只吃硬碟，還會拖慢 GPU（reporter 是同步 I/O）。
- 成本量級： 你自己在 H100 上量過一個 272k 原子的系統跑到 207 到 210 ns/day（4 fs、HMR）。5F9R 是 535k 原子，會明顯更慢；確切數字看你的
prod_<seed>.log裡的speed欄，那是最誠實的來源。
06When什麼時候用
- 4 fs + HMR： 當你要的是構形與介面的平衡性質（接觸、BSA、鹽橋、RMSF），而且需要多條 replica 才有統計力時。這就是你目前的情境。
- 2 fs +
constraints=HBonds、不用 HMR： 當你要比較時間尺度、算擴散、或要跟不含 HMR 的文獻結果對齊時。也適合當 4 fs 的對照組。 - 1 fs、不加任何 constraint： 當你要看含氫的振動行為，或在做力場開發／驗證。這種需求在你的 workflow 裡幾乎不會出現。
- NPT： 平衡密度階段，以及任何你在乎體積、密度、或系統含有會改變體積的過程時。溶劑化的生物分子 production 用 NPT 是安全的預設。
- NVT： 密度已經收好、而且你想要拿掉 barostat 這個額外的隨機來源時（例如做嚴格的 replica 比較）。
- 動 cutoff（例如試 1.2 nm）： 只當敏感度測試，不當升級。cutoff 是力場參數化的一部分，amber14 這一系配的是 0.8 到 1.0 nm 搭 PME；1.2 nm 是 CHARMM 那邊的慣例，換力場才跟著換。要測就固定其他所有參數只動這一個，然後比較你真正要報的量。
07When NOT什麼時候別用
-
不要用 Berendsen thermostat 或 Berendsen barostat 跑 production。 為什麼會壞：它不產生 canonical 分布，動能漲落被壓窄，而且動能會從高頻模式漏進質心平移（flying ice cube）。你所有靠漲落定義的量 —— RMSF、接觸數的變異、任何誤差棒 —— 都會系統性偏掉。怎麼看出它壞了：把 state log 裡的瞬時溫度畫直方圖，跟理論的相對寬度 √(2/Ndf) 比；Berendsen 下這個分布會明顯太窄。注意分母是自由度數不是 3N —— 你用了
constraints=HBonds加剛體水，實際自由度明顯少於 3N，直接跟 OpenMM 要（integrator算溫度時用的就是它），不要自己用 3N 湊。另外檢查質心速度有沒有隨時間單調長大。 -
不要在沒有 constraints 也沒有 HMR 的情況下用 4 fs。 為什麼會壞：X–H 伸縮週期約 11 fs，4 fs 一步等於一個週期只取樣不到三個點，積分器直接發散。怎麼看出它壞了：這個最好認 —— log 裡溫度飆到幾千 K，或者直接吐
Particle coordinate is NaN。它會在幾百皮秒內就爆掉，不會等到你分析才發現。 -
不要以為「沒爆掉」就等於 timestep 是安全的 —— 這是 Langevin 動力學特有的陷阱。 為什麼會壞：隨機 thermostat 會持續把能量拉回設定溫度，所以 timestep 太大造成的積分誤差不會表現成能量漂移或崩潰，而是表現成構形分布的靜默偏移。溫度計永遠是 300 K，但你的鍵角分布、二面角族群比例已經偏了。怎麼看出它壞了：唯一可靠的方法是做對照 —— 拿同一個平衡好的結構，跑一小段 2 fs（不含 HMR）與 4 fs，比較你真正要報的那幾個量（BSA、接觸數）的分布。差異在誤差棒內才算通過。另一個便宜的檢查是關掉 thermostat 跑一小段 NVE，看能量漂移速率。
-
HMR 只加到 1.5 amu 卻用 4 fs，不要當成無代價 —— 但也不要當成超速。 為什麼會壞：這一條的風險比直覺小、比零大，而把它講清楚本身就是方法段的一段。先把兩個生態的慣例攤開：Amber 那邊 4 fs 幾乎一律配 3.024 Da（ParmEd
HMassRepartition的預設值），OpenMM 使用者手冊的 heavy-hydrogen 範例用的是 4 amu 搭constraints=AllBonds。而同一份 OpenMM 手冊也寫著：光是constraints=HBonds，Langevin 動力學就已經可以用「約 4 fs」（Verlet 才是約 2 fs）。也就是說你的 4 fs 落在 OpenMM 的建議範圍內，1.5 amu 是額外的安全邊際，不是不足的 HMR。真正的麻煩是你卡在兩個慣例中間：√(1.008/1.5) ≈ 0.82 倍的頻率壓縮遠不到 3.024 Da 給的 1/√3 ≈ 0.58 倍，所以遇到 Amber 背景的 reviewer，他很可能直接說你超速，而要反駁就需要一段自己的對照數據。怎麼看出它壞了：不會崩潰，只會安靜地增加積分誤差 —— 也就是上一條講的那種陷阱，只有 2 fs 對照跑得出來。不要順手把它當成別的問題的解釋：你專案裡 DNA RMSF 的高 CV，自己的分析已經指向建模補出來的 5′ 單股末端在甩，groove 接觸計數的 145.5% 則主要是 CV 分母趨近零的假象（見指標可靠性）—— 把這些算到 HMR 頭上是沒有證據的歸因。這一點應該誠實寫進方法段與 limitation，正確的處理是補一段 2 fs 對照，不是在 paper 裡把它藏起來。 -
用了 HMR 就不要直接報動力學量。 為什麼會壞：HMR 不動位能面所以平衡性質守得住，但質量矩陣變了，擴散係數、自相關時間、速率常數全部被重新縮放。你這個專案的整個論點是「專一性由 kinetics 決定」，在含 HMR 的軌跡上談時間尺度，reviewer 一定會打。怎麼看出它壞了：跑一條不含 HMR 的 2 fs 對照，比較同一個量的自相關時間；兩者若差一個顯著的倍率，你的時間軸就不能直接引用。
-
不要在 production 還留著 position restraint 或忘記歸零。 為什麼會壞：restraint 會把你要量的柔軟度直接壓掉，RMSF 全面偏低、接觸數偏穩定，看起來反而「更收斂」—— 這是最陰險的一種假象。怎麼看出它壞了：RMSF 曲線異常平坦、連末端 loop 都很硬；或者直接檢查
context.getParameter("k")是不是 0。 -
不要對帶大量淨電荷的系統忽略中和與離子平衡。 為什麼會壞：Cas9 加 sgRNA 加完整 R-loop DNA（sgRNA 把 target 股搶走、non-target 股被擠成單股的那個三股結構）是一個高度帶電的組合。PME 對非中性系統會隱含一個均勻背景電荷，跟你的淨電荷交互作用產生系統性偏差；離子分布本身也需要時間才凝聚到 DNA 周圍。怎麼看出它壞了：印出系統淨電荷確認是 0；把 DNA 周圍的 Na⁺ 徑向分布對時間畫出來，看它在 production 開始前有沒有 plateau。沒 plateau 就是你的「平衡」還沒到。
這是最容易把方法段寫錯的一個細節，因為它會改變水的轉動與擴散，進而影響所有溶劑主導的過程。OpenMM 的 createSystem 文件寫得很直接：若用 rigidWater 把水鎖成剛體，水的氫不會被更動。你的設定正是 rigidWater=True，所以水的氫維持 1.008 amu，方法段就照這樣寫。（Amber 那邊一樣，ParmEd 的 HMassRepartition 預設不動水，要動得加 dowater。）但驗證比引用文件便宜：把 system 載進來，印幾個水分子裡氫原子的 system.getParticleMass(i)，確認是 1.008 而不是 1.5。一行程式的事。
08Project Lens跟我的 project 多相關
為什麼是四星而不是五星： 因為 production 設定不是這個專案的瓶頸所在。你自己的可靠性分析已經指認出真正的元兇：逐殘基 top-8 共享率從 Chai-1 起始的 0/8、ρ≈−0.03（隨機），換成 5F9R 晶體起始就跳到 3/8、Spearman +0.77 —— 元兇是起始結構，不是 integrator 設定。把 timestep 從 2 fs 調到 4 fs 讓你多跑一倍 replica，這是很實在的收穫，但它改變的是統計力，不是系統誤差。真正決定答案對不對的那一層在系統建立與結構來源。
為什麼不是更低：
- 「便宜」是這個專案的硬約束，而這一頁就是便宜的來源。 在共享叢集上，同時可跑的 job 數與可用儲存空間都有硬性上限。4 fs 加 HMR 加菱形十二面體盒（535k 原子，同樣的最小鏡像距離換成立方盒是 762k）是讓「4 replica × 87 ns」這件事在這個配額下成立的直接原因。少了這幾個決定，同樣的預算大概只夠一半的軌跡 —— 而在 n=2 的時候，那個 6/8 看起來就是可以用的。
- 一個錯的 thermostat 可以無聲地作廢整張可靠性表。 你報的每一個 CV 都是漲落的統計量。系綜錯了，表上每一格都錯，而且錯得看起來很正常。這種「一失全失」的風險值得四星。
- 它是 methodology paper 方法段必須寫清楚的部分。 這一段寫得夠精確，reviewer 才會往下讀你的可靠性論證。
具體怎麼用： 把這一頁的設定當成 workflow 的基礎設施層，一次定案、全系統統一、寫進方法段就不要再動。要比較 WT 與變體時，production 參數必須逐字相同，否則差異來源說不清。唯一還掛著的是 HMR 的質量 —— 它在 OpenMM 的建議範圍內、卻在 Amber 慣例之外，所以該補的是一段短的 2 fs 對照（成本以小時計，不是以 GPU-day 計），而不是等 reviewer 問。要注意這個對照回答的是「積分誤差有多大」，不是「逐位點量為什麼不穩」—— 後者的元兇你已經另外指認過了（起始結構、建模補出的區段、稀疏計數的分母）。
09Reviewer Thinkingreviewer 會問什麼
你用 4 fs 但氫只重分配到 1.5 amu，有沒有對照 2 fs 驗證過？點開看參考答案
分兩層答。第一層先把慣例講清楚：OpenMM 使用者手冊寫的是「constraints=HBonds 之下，Langevin 動力學約可用 4 fs」，所以這個 timestep 本身就在建議範圍內，1.5 amu 是額外邊際；Amber 生態則習慣 3.024 Da 配 4 fs，粗估 1.5 amu 只把可用餘裕放大 √(1.5/1.008) ≈ 1.22 倍而不是 √3 ≈ 1.73 倍。第二層是把可防守的對照講清楚：拿一條 replica 的平衡結構，跑短段 2 fs 無 HMR 對照，比較你實際報的聚合量（BSA、總接觸數、鹽橋數）的分布是否落在跨 replica 誤差棒之內。聚合量對積分誤差的敏感度遠低於逐位點量，所以這個對照很可能會通過 —— 而通過的證據本身就是論文裡的一段。
你聲稱在取樣哪個系綜？thermostat 和 barostat 分別是什麼？點開看參考答案
NPT。LangevinMiddleIntegrator（LFMiddle 離散化，摩擦 1/ps，300 K）加 MonteCarloBarostat（1 bar、300 K、每 25 步）。兩者都是產生正確系綜分布的方法，不是 Berendsen 那類會壓抑漲落（恆溫器壓窄動能漲落、恆壓器壓窄體積漲落）、因此根本沒在取樣正確系綜的弱耦合方案。方法段要把摩擦係數與 barostat 頻率寫出來（頻率 25 是 OpenMM 的建構子預設值，照實標明），因為前者會影響動力學觀測量、後者會影響體積漲落。
你的主張是專一性由 kinetics 決定，但你用了 HMR。怎麼談時間？點開看參考答案
這是對這個專案最致命的一問，因為它把方法選擇和核心論點對撞。誠實的答案分兩層：第一，目前的 MD 只用來讀平衡性質（介面接觸、BSA、構形分布），HMR 不動位能面，這些量原則上不受影響；第二，任何要引用時間尺度的分析都必須另外跑不含 HMR 的軌跡，或明確標示時間軸已被質量重分配縮放。先自己畫出這條界線，比被 reviewer 畫出來好。
你的四條 replica 有沒有真的獨立？點開看參考答案
種子那一層沒有爭議：integrator、barostat、初始速度指派各自吃同一個 SEED，而四條的 SEED 不同，所以是四組不同的 Langevin 噪音實現與初速。真正的限制在上一層，要主動講出來 —— 它們共享同一個 minimize 過的起始座標、同一次溶劑化與離子放置、同一套力場、同一組質子化狀態。這件事要在方法段講明白，因為它決定了你的跨 replica CV 在量什麼 —— 它量的是「同一個起始盆地內的取樣散布」，不是「起始結構不確定性」。後者你剛好有獨立證據：Chai 起始與 5F9R 起始給出完全不同的逐殘基結果。
系統帶大量淨電荷，你怎麼中和？離子平衡到了嗎？點開看參考答案
0.15 M NaCl，neutralize=True，所以系統淨電荷為零，PME 的均勻背景電荷問題不成立。離子平衡是另一回事：DNA 周圍的反離子凝聚需要時間，要證明它到了，就把 Na⁺ 對 DNA 磷酸的徑向分布或配位數對時間畫出來，看有沒有 plateau。這個檢查很便宜，而且很少人報 —— 主動報反而是加分。
10Common Mistakes最常犯的錯
- 用 Berendsen 控溫跑 production。 後果：系綜不對，所有漲落量與誤差棒失真，而且溫度看起來完全正常。正解：用 Langevin（LFMiddle／BAOAB）或 Nosé-Hoover 這類產生正確 canonical 分布的方法；barostat 用 Monte Carlo 型。
- 把 timestep 當純效能參數，看沒崩潰就放心。 後果：隨機 thermostat 掩蓋了積分誤差，構形分布靜默偏移。正解：做 2 fs 對照，比較你真正要報的量；必要時跑一小段 NVE 看能量漂移。
- 把 timestep、constraints、HMR 三個數字拆開，各自引用不同來源的慣例。 後果：組合出來的設定沒有任何一篇文獻背書，被問就答不出來。正解：整組引用同一個生態的慣例（Amber：3.024 Da 配 4 fs；OpenMM：
HBonds+ Langevin 本身已經是約 4 fs），任何偏離都自己補一段對照數據。 - 在含 HMR 的軌跡上直接讀速率或擴散。 後果：時間軸被縮放過，數字不能跟實驗或其他模擬比。正解：動力學分析另跑無 HMR 對照，或在文中明確標註縮放。
- WT 與變體用了不同的 production 參數。 後果：差異來源無法歸因，整組比較作廢。正解：參數一次定案、全系列統一，並在補充材料放完整的輸入檔。
- 軌跡存太密。 後果：硬碟爆掉、reporter 的同步 I/O 把 GPU 拖慢，而且相鄰幀高度相關、多存的幀不增加有效樣本數。正解：先估你要分析的量的自相關時間，用它決定 cadence；20 ps 對介面接觸類的量是合理起點。
- 忘記把 restraint 歸零就開始 production。 後果：柔軟度被壓掉、結果看起來異常穩定，最像「收斂」的假象。正解：production 開始前顯式
setParameter("k", 0.0)，並用 RMSF 曲線的末端是否柔軟當作 sanity check。
11Further Reading讀哪幾篇
12Summary三句話
Production 設定是一組換匯決定：constraints 換掉振動自由度、HMR 換掉質量分布、大 timestep 換掉積分精度，換回來的都是牆鐘時間與 replica 數。這些交易絕大多數划算，但 thermostat 是例外 —— 它改變的不是精度而是你在取樣的系綜，選錯了跑再久也不會對。你目前的 1.5 amu 配 4 fs 落在 OpenMM 手冊的建議範圍內、卻在 Amber 慣例之外，補一段短的 2 fs 對照就能把它從「會被質疑的點」變成「方法段的一句話」。