Algorithms Teaching Series · Advanced

模擬退火演算法完整教材

Simulated Annealing

Theory · Design · Applications


向冶金學借一個概念,解決最佳化最頑固的敵人——局部最優。 從 Boltzmann 分布與 Metropolis 準則出發,證明「為什麼偶爾往上爬是對的」與 「降溫多慢才保證收斂」(Hajek 定理),再把它變成一個通用引擎, 實戰 TSP、背包、N 皇后、圖著色、Max-Cut、排程與數獨——十三個可執行、可對照的 C++ 實例。

配套程式:annealing.cpp(C++17,單檔、泛型退火引擎、固定隨機種子、結果可重現)

2026 年 10 月

鐵匠要鍛一把好劍。他把鐵燒到通紅——高溫下原子劇烈振動、到處亂跑,什麼結構都可能出現;然後他慢慢讓它冷卻,原子一邊失去能量、一邊有足夠時間找到彼此最舒服的位置,最後排成整齊的晶格,這就是強度最高、能量最低的狀態。如果他把燒紅的鐵直接丟進水裡(淬火,quenching),原子來不及重排、被凍在混亂的位置,劍又脆又容易斷。

最佳化也一樣。「只往下走」的爬山法就是淬火——快,但凍在第一個遇到的谷底。模擬退火則給演算法一個「溫度」:高溫時容忍往上爬(走出谷底),然後慢慢降溫,讓它最終安定在最深的那個谷。

問題:局部最優

絕大多數有趣的最佳化問題——旅行推銷員、排程、晶片布局、蛋白質摺疊——都有一個共同的麻煩:解空間裡有天文數字個局部最優。以 TSP 為例,\(n\) 個城市有 \((n-1)!/2\) 條路線,隨便一個「鄰居都比我差」的路線就是局部最優,而這樣的路線多到任何從隨機起點出發、只往下走的方法幾乎必然卡住。

定義 1.1 (局部最優). 給定解空間 \(\mathcal{S}\)、能量(目標)函數 \(E:\mathcal{S}\to\mathbb{R}\) 與鄰域結構 \(\mathcal{N}(x)\subseteq\mathcal{S}\)。\(x\) 是局部最優若 \(E(x)\le E(y)\) 對所有 \(y\in\mathcal{N}(x)\);是全域最優若 \(E(x)\le E(y)\) 對所有 \(y\in\mathcal{S}\)。全域最優集合記為 \(\mathcal{S}^*\),最低能量記為 \(E^*\)。

逃離局部最優只有兩條路:換鄰域(讓它不再是局部最優,例如 2-opt 換成 3-opt)或容許暫時變差。模擬退火走第二條路,而且用物理學給了「該容許多少」一個漂亮的答案。

歷史

  • 1953,Metropolis 等人:為了在電腦上模擬液體的平衡態,提出以「往下一定走、往上以 \(e^{-\Delta E/T}\) 走」產生 Boltzmann 分布樣本的演算法——這是本教材的數學心臟(第 2 章)。

  • 1983,Kirkpatrick、Gelatt 與 Vecchi:在 IBM 做晶片布局時意識到「最佳化 = 找最低能量態」,把 Metropolis 演算法加上降溫,發表於 Science,命名「模擬退火」;論文同時展示了晶片布局、佈線與 TSP。

  • 1985,Černý:獨立提出同樣的想法(應用於 TSP)。

  • 1984,Geman & Geman:證明對數降溫 \(T_k\ge c/\ln k\) 保證收斂到全域最優(影像復原脈絡)。

  • 1988,Hajek:給出收斂的充要條件——\(c\) 必須至少等於最深的「非全域局部最優」的深度。這是模擬退火理論最深刻的結果。

  • 1985,Hinton & Sejnowski:Boltzmann 機——同一套 Metropolis + 溫度的思想進入神經網路。

  • 1990 年代至今:門檻接受法、自適應 SA、平行回火、量子退火(D-Wave,1999 年理論、2011 年商用);模擬退火至今仍是 VLSI 佈局、排程軟體與許多工業最佳化器裡的工作馬。

一句話定義與物理對照

在解空間上跑一條 Metropolis 馬可夫鏈(往下必走、往上以機率 \(e^{-\Delta E/T}\) 走),同時讓溫度 \(T\) 從高到低緩慢下降。固定溫度下鏈的平穩分布是 Boltzmann 分布 \(\pi_T(x)\propto e^{-E(x)/T}\);\(T\to0\) 時 \(\pi_T\) 集中在全域最優集合 \(\mathcal{S}^*\)。

物理退火 模擬退火 本教材的例子
物質的微觀狀態 一個可行解 \(x\in\mathcal{S}\) 一條 TSP 路線、一組皇后位置
能量 目標函數 \(E(x)\) 路線長度、衝突數、\(-\)價值
原子的微小位移 鄰域移動 \(x\to y\in\mathcal{N}(x)\) 2-opt 反轉、翻一個位元、交換兩格
溫度 \(T\) 控制參數 \(T\) 接受變差的意願
Boltzmann 分布 固定 \(T\) 時鏈的平穩分布 Ex 3 直接驗證
緩慢冷卻 降溫排程 \(T_0>T_1>\cdots\) 幾何、線性、對數、Lundy–Mees
基態(完美晶格) 全域最優 \(\mathcal{S}^*\) Held–Karp 最優、DP 最優
淬火(凍出缺陷) 降溫太快 \(\Rightarrow\) 卡在局部最優 Ex 4c 的 \(\alpha=0.8\)
比熱峰(相變) \(\operatorname{Var}(E)/T^2\) 的峰值 Ex 4 軌跡圖

在最佳化工具箱中的位置

模擬退火與遺傳演算法同屬元啟發式,但哲學不同:SA 是單點、靠溫度逃離;GA 是族群、靠重組探索。SA 的獨特之處在於它是元啟發式家族裡理論最完整的成員——有嚴格的收斂定理,而且參數(溫度)有清楚的物理與機率意義。

數學基礎:為什麼偶爾往上爬是對的

本章用四個定理回答四個問題:Metropolis 準則產生什麼分布(細緻平衡)?低溫時這個分布長什麼樣(集中定理)?降溫多慢才保證找到全域最優(Hajek 定理)?以及為什麼實務上不照定理做。

Boltzmann 分布與配分函數

定義 2.1 (Boltzmann(Gibbs)分布). 有限解空間 \(\mathcal{S}\)、能量 \(E\)、溫度 \(T>0\)。Boltzmann 分布為 \[\pi_T(x)=\frac{e^{-E(x)/T}}{Z_T},\qquad Z_T=\sum_{y\in\mathcal{S}}e^{-E(y)/T}\quad\text{(配分函數)}.\]

它是統計物理的基本定律:在溫度 \(T\) 達到平衡的系統,處於狀態 \(x\) 的機率正比於 \(e^{-E(x)/T}\)。對最佳化而言只需記住兩個性質:

性質 2.1.

  1. 能量越低機率越高,且兩狀態的機率比只與能量差有關、與 \(Z_T\) 無關: \[\frac{\pi_T(y)}{\pi_T(x)}=e^{-(E(y)-E(x))/T}.\]

  2. \(T\to\infty\) 時 \(\pi_T\to\) 均勻分布(溫度極高,能量無關緊要);\(T\to0\) 時 \(\pi_T\) 集中到 \(\mathcal{S}^*\)(下節定理)。

定理 2.1 (低溫極限). 令 \(E^*=\min E\),\(\mathcal{S}^*=\{x:E(x)=E^*\}\)。則對任意 \(x\), \[\lim_{T\to0^+}\pi_T(x)=\begin{cases}1/|\mathcal{S}^*| & x\in\mathcal{S}^*,\\ 0 & x\notin\mathcal{S}^*.\end{cases}\]

Proof. 分子分母同除 \(e^{-E^*/T}\): \[\pi_T(x)=\frac{e^{-(E(x)-E^*)/T}}{\sum_y e^{-(E(y)-E^*)/T}} =\frac{e^{-(E(x)-E^*)/T}}{|\mathcal{S}^*|+\sum_{y\notin\mathcal{S}^*}e^{-(E(y)-E^*)/T}}.\] 對 \(y\notin\mathcal{S}^*\) 有 \(E(y)-E^*>0\),故 \(e^{-(E(y)-E^*)/T}\to0\);分母 \(\to|\mathcal{S}^*|\)。分子在 \(x\in\mathcal{S}^*\) 時恆為 1,否則 \(\to0\)。 ◻

這就是模擬退火的全部策略:若能在每個溫度都從 \(\pi_T\) 抽樣,那麼把 \(T\) 降到接近 0,抽到的就是全域最優。剩下的問題是:怎麼從 \(\pi_T\) 抽樣?\(Z_T\) 是對整個解空間求和,TSP 的 \(|\mathcal{S}|=(n-1)!/2\),根本算不出來。

Metropolis 準則與細緻平衡

Metropolis 的洞見:不需要 \(Z_T\),只需要相鄰兩狀態的機率比。

目前解 \(x\),溫度 \(T\) 從鄰域 \(\mathcal{N}(x)\) 均勻隨機挑一個 \(y\) \(\Delta E\gets E(y)-E(x)\) \(x\gets y\) \(x\gets y\)

定義 2.2 (細緻平衡). 馬可夫鏈的轉移機率 \(K(x,y)\) 對分布 \(\pi\) 滿足細緻平衡若 \(\pi(x)K(x,y)=\pi(y)K(y,x)\) 對所有 \(x,y\)。

引理 2.1. 細緻平衡 \(\Rightarrow\) \(\pi\) 是平穩分布:\(\sum_x\pi(x)K(x,y)=\pi(y)\)。

Proof. \(\sum_x\pi(x)K(x,y)=\sum_x\pi(y)K(y,x)=\pi(y)\sum_xK(y,x)=\pi(y)\)。 ◻

定理 2.2 (Metropolis 鏈的平穩分布). 設鄰域對稱(\(y\in\mathcal{N}(x)\iff x\in\mathcal{N}(y)\))且大小一致(\(|\mathcal{N}(x)|=|\mathcal{N}(y)|\),或提議機率對稱 \(q(x,y)=q(y,x)\))。固定溫度 \(T\) 的 Metropolis 鏈對 Boltzmann 分布 \(\pi_T\) 滿足細緻平衡,故 \(\pi_T\) 為其平穩分布。

Proof. 對 \(x\ne y\),\(K(x,y)=q(x,y)\min\{1,e^{-(E(y)-E(x))/T}\}\)。不妨設 \(E(y)\ge E(x)\),則 \[\begin{align*} \pi_T(x)K(x,y)&=\frac{e^{-E(x)/T}}{Z_T}\,q(x,y)\,e^{-(E(y)-E(x))/T}=\frac{e^{-E(y)/T}}{Z_T}\,q(x,y),\\ \pi_T(y)K(y,x)&=\frac{e^{-E(y)/T}}{Z_T}\,q(y,x)\cdot1 . \end{align*}\] 由 \(q(x,y)=q(y,x)\) 兩者相等。\(x=y\) 時細緻平衡自動成立。 ◻

注意證明中 \(Z_T\) 在兩邊同時出現、從未被計算——這就是不需要配分函數的原因。再加上鏈不可約(任兩解可經鄰域互達)與非週期(有「留在原地」的機率),遍歷定理保證鏈的狀態分布收斂到 \(\pi_T\)。

範例 2.1 (Ex 3:鏈真的收斂到 Boltzmann 分布). 8 個狀態排成環,能量 \(E=(0,1,1,2,3,3,4,6)\),鄰域 = 左右兩鄰居。跑 \(2\times10^6\) 步 Metropolis,統計停留頻率:

T = 2.000   state:  0     1     2     3     4     5     6     7
          E:      0     1     1     2     3     3     4     6
   empirical P:   0.312 0.189 0.189 0.114 0.069 0.070 0.042 0.015
   Boltzmann  P:  0.311 0.189 0.189 0.115 0.069 0.069 0.042 0.015
T = 0.700   empirical P:   0.637 0.153 0.153 0.037 0.009 0.009 0.002 0.000
            Boltzmann  P:  0.638 0.153 0.153 0.037 0.009 0.009 0.002 0.000
T = 0.250   empirical P:   0.965 0.018 0.017 0.000 0.000 0.000 0.000 0.000
            Boltzmann  P:  0.964 0.018 0.018 0.000 0.000 0.000 0.000 0.000

三位小數全部吻合。同時看見定理 2.1:\(T\) 從 2 降到 0.25,全域最小的機率從 31% 升到 96%。

從抽樣到最佳化:為什麼要「慢慢」降溫

把定理 2.1 與 2.2 接起來:在每個溫度跑夠久的鏈 \(\Rightarrow\) 得到 \(\pi_T\) 的樣本 \(\Rightarrow\) 把 \(T\) 降到接近 0 \(\Rightarrow\) 樣本幾乎必為全域最優。問題在「跑夠久」:

  • 低溫時鏈的混合時間(到達平穩分布所需步數)隨 \(e^{d/T}\) 爆炸,\(d\) 是必須翻越的能量障礙。\(T\) 越低、要等越久。

  • 直接在低溫從隨機起點開始 \(=\) 淬火:鏈卡在第一個谷底,混合時間比宇宙年齡長。

  • 退火的巧思:高溫時混合快,先在高溫達到平衡;每次只降一點點溫度,前一個溫度的平衡分布已經很接近下一個的,所以只需少量步數就能「追上」。這叫準平衡(quasi-equilibrium)。

Hajek 定理:降溫多慢才夠

定義 2.3 (局部最優的深度). 狀態 \(x\) 是局部最優但非全域最優。從 \(x\) 出發、沿鄰域走到任何嚴格更低能量狀態的所有路徑中,「路徑上最高點減 \(E(x)\)」的最小值稱為 \(x\) 的深度 \(d(x)\)——即逃離 \(x\) 所必須翻越的最低障礙。令 \[d^*=\max\{d(x):x\text{ 為局部最優且 }x\notin\mathcal{S}^*\}.\]

定理 2.3 (Hajek 1988). 在弱可逆性假設下(任何能以高度 \(h\) 從 \(x\) 走到 \(y\) 的路徑,反方向也存在同高度的路徑),對降溫排程 \(T_k=\dfrac{c}{\ln(k+2)}\),模擬退火滿足 \[\lim_{k\to\infty}\mathbb{P}\big(x_k\in\mathcal{S}^*\big)=1\quad\Longleftrightarrow\quad c\ \ge\ d^*.\]

這是一個充要條件:\(c\) 太小(降溫太快),存在正機率永遠困在某個深度 \(>c\) 的局部最優;\(c\ge d^*\) 則保證收斂(機率意義)。直觀:在溫度 \(T\) 逃離深度 \(d\) 的谷需要大約 \(e^{d/T}\) 步;若 \(T_k=c/\ln k\),則 \(e^{d/T_k}=k^{d/c}\),而鏈到第 \(k\) 步已累積 \(k\) 步機會——當 \(d\le c\) 時 \(k^{d/c}\le k\),機會足夠;\(d>c\) 時 \(k^{d/c}\gg k\),永遠追不上。(Geman & Geman 1984 的早期版本是充分條件 \(c\ge\max E-\min E\),比 \(d^*\) 粗得多。)

對數降溫慢到荒謬:要讓 \(T_k\) 從 \(c\) 降到 \(c/100\),需要 \(k=e^{100}\approx10^{43}\) 步——超過宇宙原子數。而且 \(d^*\) 本身未知(知道它幾乎等於知道答案)。所以沒有人在實務上用對數降溫;定理的價值在於說明「慢」的方向是對的、以及為什麼過快會失敗。實務上用幾何降溫 \(T\gets\alpha T\),接受「以高機率得到好解」而非「以機率 1 得到最優解」——這是有限時間退火(finite-time annealing)的觀點。

比熱與「關鍵溫度」

統計物理告訴我們哪個溫度最重要。定義平均能量 \(\langle E\rangle_T=\sum_x E(x)\pi_T(x)\),則

定理 2.4 (能量漲落與比熱). \[C(T)\;\coloneqq\;\frac{d\langle E\rangle_T}{dT}=\frac{\operatorname{Var}_{\pi_T}(E)}{T^2}.\]

Proof. 寫 \(\langle E\rangle_T=\dfrac{\sum_x E(x) e^{-E(x)/T}}{Z}\),\(Z=\sum_x e^{-E(x)/T}\)。對 \(T\) 微分時用 \(\frac{d}{dT}e^{-E/T}=\frac{E}{T^2}e^{-E/T}\),得 \[\frac{d\langle E\rangle_T}{dT} =\frac{1}{T^2}\left(\frac{\sum_x E^2 e^{-E/T}}{Z}-\Big(\frac{\sum_x E e^{-E/T}}{Z}\Big)^2\right) =\frac{\langle E^2\rangle_T-\langle E\rangle_T^2}{T^2} =\frac{\operatorname{Var}_{\pi_T}(E)}{T^2}.\qedhere\] ◻

\(C(T)\) 大 \(\Rightarrow\) 溫度稍降、能量大幅下降 \(\Rightarrow\) 系統正在「結晶」(解的大尺度結構正在成形)。這是物理相變的訊號,也是最佳化中最該慢慢走的溫度區間。第 5.3 節的 TSP 軌跡會看到一個清楚的比熱峰。這個觀察催生了自適應降溫:\(C(T)\) 大時降溫慢、小時快。

演算法設計:把定理變成程式

完整演算法

初始解 \(x\)、初溫 \(T_0\)、終溫 \(T_{\min}\)、降溫率 \(\alpha\in(0,1)\)、每溫度步數 \(L\) \(T\gets T_0\),\(x^{\mathrm{best}}\gets x\) 隨機挑鄰居 \(y\in\mathcal{N}(x)\),計算 \(\Delta E=E(y)-E(x)\) \(x\gets y\) \(x^{\mathrm{best}}\gets x\) \(T\gets\alpha T\) \(x^{\mathrm{best}}\)

泛型引擎:問題與演算法分離

配套程式把演算法寫成一個模板函數 anneal(P, N, temp),任何「問題」只要提供五個方法就能退火:

方法 責任 TSP 的實作
energy() 回傳目前解的能量(只在開頭呼叫一次) 路線總長 \(O(n)\)
propose() 隨機產生鄰居、回傳 \(\Delta E\)(可暫定套用) 挑 \(i<j\),算兩條邊的差 \(O(1)\)
commit() 接受:把暫定移動正式化 反轉 tour[i..j]
undo() 拒絕:還原 無事可做
saveBest() 引擎通知「這是新的最佳解」 複製路線

引擎另外接受一個排程函數 temp(k, N)(第 \(k\) 步、共 \(N\) 步的溫度),所以降溫策略也是可抽換的。這個分離讓本教材 13 個實例共用同一份 40 行的引擎——也是工程上推薦的寫法:退火邏輯只寫一次、測一次。

鄰域設計:最重要的決策

  1. 連通:任兩解可經有限步鄰域移動互達(否則鏈不可約、定理失效、有些解永遠到不了)。

  2. 對稱:\(y\in\mathcal{N}(x)\iff x\in\mathcal{N}(y)\),Metropolis 的簡單比值才成立;非對稱提議要改用 Metropolis–Hastings 修正。

  3. 小步:一步只改一點點,\(\Delta E\) 的尺度小、接受率才有意義。「把整個解打亂重來」不是鄰域,是重啟。

  4. 增量 \(\Delta E\):永遠不要為了算 \(\Delta E\) 重算整個目標函數。TSP 的 2-opt 是 \(O(1)\) 而非 \(O(n)\),差 60 倍速度。

  5. 尊重約束:能用表示法保證可行(Sudoku 的區塊排列、N 皇后的排列)就不要用罰分;非用罰分不可時(背包)罰分係數要夠大。

初溫 \(T_0\)

初溫應讓鏈一開始「幾乎什麼都接受」(接受率 \(\chi_0\approx0.8\)),否則等於從低溫開始淬火。

隨機走 \(m\) 步,收集所有上坡移動的 \(\Delta E^+\) 平均值 \(\overline{\Delta E^+}\)。令初始接受率為 \(\chi_0\),則 \[\chi_0\approx e^{-\overline{\Delta E^+}/T_0}\quad\Longrightarrow\quad T_0=-\frac{\overline{\Delta E^+}}{\ln\chi_0}.\] \(\chi_0=0.8\) 時 \(T_0\approx4.5\,\overline{\Delta E^+}\)。

隨機初始解的鄰域 \(\Delta E\) 常常遠大於好解附近的 \(\Delta E\)(隨機 TSP 路線有很多長邊,2-opt 的 \(\Delta E\) 動輒上百;好路線附近只有個位數)。Ex 6 校準出 \(T_0=136\),實測 \(T_0\) 從 13.6 到 13614 效果幾乎相同(585–592),但 \(T_0=1.36\) 明顯變差(601)。教訓:\(T_0\) 寧高勿低——太高只是浪費前面幾個溫度層;太低則是淬火。

降溫排程

排程 公式 特性 Ex 5(TSP 60,同預算)
幾何 \(T_{k+1}=\alpha T_k\) 每個「數量級」花同樣時間;業界預設 633.8(最好之一)
線性 \(T_k=T_0-\beta k\) 高溫段太久、低溫段太短 640.1
Lundy–Mees \(T_{k+1}=\dfrac{T_k}{1+\beta T_k}\) 初期驟降(雙曲線),像淬火 668.1(差)
對數 \(T_k=c/\ln(k+2)\) 理論保證;實務幾乎是常溫 \(c{=}30\):627.1;\(c{=}3\):680.1
門檻接受 接受若 \(\Delta E<T_k\) 確定性,不需 \(\exp\) 與亂數 640.2
自適應 \(C(T)\) 大時慢降 把時間花在相變附近 見練習

鏈長、終止與其他細節

  • 每溫度步數 \(L\):至少要讓每個「自由度」有機會被動到——TSP 取 \(L\propto n\) 或 \(n^2\);或動態:接受次數達 \(L_{\text{acc}}\) 或嘗試次數達 \(L_{\max}\) 即降溫(Aarts & Korst)。

  • 終止:\(T<T_{\min}\);或連續數個溫度層最佳解不改善;或接受率低於 1%(鏈已凍結,繼續跑是浪費)。

  • 記錄最佳解:Metropolis 鏈會離開好解,一定要另存 \(x^{\mathrm{best}}\)。這是 Ex 5 中「幾乎常溫」的對數排程仍得好解的原因。

  • 重新加熱(reheating):凍結後若未達目標(Sudoku 衝突 \(>0\)),把 \(T\) 拉回中溫再退一次——比重頭來過便宜,因為保留了已成形的結構。

  • 重啟 vs. 長鏈:同預算下,幾條「夠長」的鏈取最好通常不輸一條長鏈且更穩(sd 小),但鏈短到還沒結晶就停會崩潰(Ex 7:1000 條 2000 步的鏈,解差 17%)。

  • 隨機種子:固定以利重現;正式實驗多種子取平均並報告標準差。

實驗 I:連續函數與理論驗證

一維多峰函數:爬山法 vs. 退火

\(f(x)=0.1x^2+4\sin2x+3\cos3x\),\(x\in[-10,10]\),全域最小 \(f^*=-6.5687\) 在 \(x^*=-0.9454\)(封面圖)。鄰域:\(y=x+0.5\,\varepsilon\),\(\varepsilon\sim\mathcal N(0,1)\)。1000 個隨機起點、各 20000 步:

hill climbing (T=0)     success  18.2%   mean final f = -2.452
SA  T0 = 1.00           success  27.8%   mean final f = -3.922
SA  T0 = 4.00           success  65.3%   mean final f = -5.424
SA  T0 = 16.0           success  63.6%   mean final f = -5.352

爬山法只在起點恰好落在全域最小的吸引盆(約 18% 的區間)時成功。\(T_0=1\) 太低——障礙高度約 5–10,\(e^{-8/1}\approx0\),與爬山法差不多。\(T_0=4\) 起成功率躍升到 65%,再高沒有更好(\(T_0=16\) 浪費了前段預算)。

Rastrigin:步長與溫度一樣重要

\(f(x,y)=20+x^2-10\cos2\pi x+y^2-10\cos2\pi y\),局部最小在所有整數格點、彼此相距 1、障礙約 10;全域最小 \(f=0\) 在原點。

sigma   method                steps     success     mean final f
0.15    hill climbing         20000     1%          17.352
0.15    SA (T0=5)             20000     4%          9.733
0.15    SA (T0=5), 10x longer 200000    11%         3.165
0.30    hill climbing         20000     84%         0.186
0.30    SA (T0=5)             20000     86%         0.144
0.30    SA (T0=5), 10x longer 200000    100%        0.002
0.50    hill climbing         20000     100%        0.005
0.50    SA (T0=5)             20000     100%        0.012
0.50    SA (T0=5), 10x longer 200000    100%        0.002

步長 \(\sigma=0.15\)(遠小於盆地間距 1)時,無論溫度多高、跑多久,成功率都很低——鄰域太小,每一步都在同一個盆地裡打轉,溫度只能幫它在盆地內亂跳。\(\sigma=0.5\) 時連爬山法都 100%——因為高斯步長有可觀機率直接跳進鄰近盆地,「局部」搜尋已經不局部了。真正展示退火價值的是 \(\sigma=0.3\):爬山法 84%、退火 100%。鄰域設計錯了,溫度救不了你;鄰域設計對了,溫度是臨門一腳。

理論驗證:Boltzmann 分布

見第 2 章範例——Metropolis 鏈的經驗頻率與 \(e^{-E/T}/Z\) 三位小數吻合,且 \(P(\text{全域最小})\) 隨 \(T\) 下降從 0.31 升到 0.96。

實驗 II:TSP 深入剖析

TSP 是模擬退火的經典戰場(Kirkpatrick 1983 的原始論文之一)。本章用一個 60 城隨機實例(座標 \([0,100]^2\))做五組實驗。鄰域為 2-opt 反轉、\(\Delta E\) 為 \(O(1)\)。

對照精確解與局部搜尋

=== Example 4a: TSP n=15, SA vs Held-Karp exact optimum (10 independent SA runs) ===
Held-Karp OPT = 368.64   SA best over runs = 368.64   mean = 368.64 +- 0.00   hit OPT in 10/10 runs

=== Example 4b: TSP n=60, SA vs nearest-neighbour + 2-opt local search ===
nearest neighbour        = 782.28
NN + 2-opt local search  = 671.42
random + 2-opt (best/20) = 653.77   mean 689.81
SA 2e6 steps (best/20)   = 645.34   mean 647.52 +- 2.45

15 城:10 次全部命中 Held–Karp 最優。60 城:SA 的平均(647.5)比 20 次隨機起點 2-opt 局部搜尋的最好(653.8)還好,且標準差只有 2.45。局部搜尋的 20 次結果散佈在 654–720 之間——這就是「卡在不同局部最優」的長相。

降溫速度掃描

每個溫度層 3000 步,\(T_0=30\to T_{\min}=0.01\),改變 \(\alpha\):

alpha   levels      steps       mean        sd        best
0.800   36          108000      659.80      9.26      647.15
0.950   157         471000      649.16      3.37      645.34
0.990   797         2391000     646.83      1.89      645.34
0.995   1598        4794000     646.97      1.24      645.34
0.999   8003        24009000    646.25      0.95      645.34

\(\alpha=0.8\) 是淬火:平均差 2%、sd 大。\(\alpha\ge0.99\) 之後平均只差 0.5 以內,再慢只是把 sd 壓小——報酬遞減。實務上「\(\alpha=0.95\)–\(0.99\)、\(L\propto n\)」是穩健的起點。

退火軌跡:能量、接受率與比熱

把 400 個溫度層(各 3000 步)的統計寫到 annealing_trace.dat,畫出來:

三件事值得看:

  1. 能量在 \(T\approx 20\to3\) 之間掉了大半(2100 \(\to\) 700),之後只剩細修。這段是「結晶期」。

  2. 比熱峰在 \(T\approx9\)——漲落最大、結構正在決定的溫度。自適應排程會在這附近放慢。

  3. \(T<0.5\) 後接受率 \(<0.1\%\)、最佳解不再改善:鏈已凍結。後面 200 個溫度層(一半的預算)是浪費的——這是「接受率低於門檻就終止」的理由。

排程比較、初溫與重啟

=== Example 5: cooling schedules under an identical budget (TSP n=60, 1.2e6 steps, 10 runs each) ===
schedule                    mean        sd        best      final T
hill climbing (T = 0)       691.65      36.63     627.05    1.0e-300
constant T = 1              654.47      22.04     627.05    1.0e+00
linear                      640.12      12.61     627.05    1.0e-02
geometric                   633.82      13.15     627.05    1.0e-02
Lundy-Mees                  668.14      22.74     638.06    1.0e-02
logarithmic c=30 (Hajek)    627.12      0.15      627.05    2.1e+00
logarithmic c=3             680.05      21.38     655.63    2.1e-01
threshold accepting         640.23      6.72      627.05    1.0e-02

對數排程 \(c=30\) 全程溫度都在 2–3 之間(1.2e6 步只讓 \(c/\ln k\) 從 30 降到 2.1),卻 10 次全部得到 627.05——這個實例的最佳解。原因:\(T\approx2.5\) 恰好是這個實例「結晶完成、但還能在好解之間遊走」的溫度,而引擎記錄最佳解,所以鏈在 1.2e6 步裡逛過的所有好解都被抓下來。這不是對數排程的勝利(\(c=3\) 就慘敗),而是「在對的溫度多待、並記錄最佳」的勝利。常溫 \(T=1\) 太冷(卡住,sd 22),爬山法最糟。

=== Example 6: calibrating T0 by target initial acceptance rate (TSP n=60, 1e6 steps, 10 runs) ===
mean uphill dE from random walk -> T0 for 80% acceptance = 136.14
T0                    mean        sd        best
1.361                 600.66      12.86     587.84
13.614                591.88      8.40      580.46
136.142               586.12      6.42      580.46
1361.418              589.03      7.63      580.46
13614.177             585.59      5.25      580.46

=== Example 7: one long chain vs restarts under the same budget (TSP n=60, 2e6 steps total, 5 trials) ===
strategy                  mean        sd        best
1 x 2000000 steps         579.27      3.54      577.63
10 x 200000 steps         577.89      0.52      577.63
100 x 20000 steps         578.19      1.11      577.63
1000 x 2000 steps         678.53      15.53     656.07

(Ex 5、6、7 各自產生新的隨機實例,數字不可跨實驗比較。)初溫:高 100 倍無害、低 100 倍有害。重啟:10 條 20 萬步的鏈比 1 條 200 萬步的鏈更穩(sd 0.5 vs 3.5);但 2000 步的鏈還沒走完結晶期就停了,解差 17%。

實驗 III:組合最佳化百寶箱

本章六個問題展示同一個引擎如何靠換表示法與鄰域處理完全不同的結構。每題給:表示法、能量、鄰域與 \(\Delta E\)、對照基準、結果。

0/1 背包:罰分能量

表示:位元字串 \(x\in\{0,1\}^n\)。能量:\(E=-\sum v_ix_i+\lambda\max\{0,\sum w_ix_i-W\}\),\(\lambda=2\max_i(v_i/w_i)\),保證「多放一單位重量的罰分 \(>\) 任何可能的價值增益」。鄰域:翻一個位元(50%)或一進一出交換(50%)——交換讓鏈能在「裝滿」的邊界上滑動而不必先經過一個很差的解。另外只記錄可行的最佳解。

W = 635, DP optimum = 1424
SA best feasible value: best = 1424.0, mean = 1422.9, hit optimum in 6/10 runs

平均離 DP 最優 0.08%。(第一版只用位元翻轉時平均差 1.3%、從未命中——鄰域加一個交換移動就解決了。)

N 皇后:表示法消滅約束

表示:排列 \(\sigma\),第 \(r\) 列的皇后在第 \(\sigma(r)\) 行——列、行約束由表示法自動滿足。能量:對角線衝突對數 \(\sum_{d}\binom{c_d}{2}\),\(c_d\) 為第 \(d\) 條對角線的皇后數(主、副對角線各 \(2n-1\) 條)。鄰域:交換兩列的皇后行位置。\(\Delta E\):只動四個對角線計數器,\(O(1)\)。

n = 8      initial conflicts = 5      final conflicts = 0  (solved within ~16 steps)
n = 50     initial conflicts = 43     final conflicts = 0  (solved within ~14200 steps)
n = 200    initial conflicts = 132    final conflicts = 0  (solved within ~60000 steps)
n = 1000   initial conflicts = 645    final conflicts = 0  (solved within ~332000 steps)

\(n=1000\) 在 33 萬步(不到一秒)解出——回溯法在這個規模已無望。這與 min-conflicts 啟發式(解過百萬皇后)同屬一家。

圖著色:固定 \(k\)、最小化衝突

表示:每點一個顏色 \(\in\{0,\dots,k-1\}\)。能量:同色邊數。鄰域:隨機點換隨機新色,\(\Delta E\) 只看該點的鄰居。策略:從 DSATUR 貪婪的 \(k\) 出發,每次 \(k\gets k-1\),直到 SA 找不到零衝突。

|E| = 505, DSATUR greedy uses 7 colours
k = 7   ->  proper colouring found
k = 6   ->  proper colouring found
k = 5   ->  failed, min conflicts = 2

SA 比貪婪少一色。\(k=5\) 失敗不證明 \(\chi(G)=6\)(SA 沒有保證),但「連 40 萬步 \(\times\) 5 次都剩 2 條衝突邊」是強烈暗示。

Max-Cut

表示:每點 \(\in\{0,1\}\)。能量:\(-\)割邊數。鄰域:翻一點,\(\Delta E=(\text{異側鄰居數})-(\text{同側鄰居數})\)。

n=20, |E|=76: brute-force optimum = 54, SA best = 54, SA mean = 54.0 (10 runs)
n=200, |E|=971: greedy local search mean = 668.6, best = 684.0   |   SA mean = 714.9, best = 715.0

\(n=20\) 對照 \(2^{19}\) 窮舉:10/10 最優。\(n=200\):貪婪局部搜尋(有改善就翻到不能翻)平均 669,SA 平均 715——高 7%,且 10 次結果幾乎一樣(714–715),顯示它們很可能逼近真正的最優。Max-Cut 是 NP-Hard 且 Goemans–Williamson 的 0.878 近似是理論最佳已知;SA 在隨機圖上的表現常常優於該保證。

平行機排程 \(P\|C_{\max}\):能量塑形

表示:每個工作指派一台機器。能量:\(C_{\max}=\max_i\text{load}_i\)。問題:\(C_{\max}\) 是平台狀的——移動非瓶頸機器上的工作 \(\Delta E=0\),鏈在平台上盲目漫遊。能量塑形:\(E=C_{\max}+\varepsilon\sum_i\text{load}_i^2\),\(\varepsilon=10^{-6}\) 小到不改變最優解的排序,卻給平台一個「往平衡方向」的坡度。鄰域:移動一個工作(50%)或交換兩台機器上的工作(50%)。

=== Example 12a: P||Cmax, n=12 jobs on m=3 machines - SA vs brute force (3^12 assignments) ===
brute-force optimum = 252, SA best = 252, hit optimum in 7/10 runs

=== Example 12b: P||Cmax, n=40 jobs on m=5 machines - LPT vs SA vs lower bound, energy shaping ===
sum p = 1978, lower bound = 396, LPT makespan = 398
SA with eps = 0e+00:  makespan best = 396, mean = 396.0, reached LB in 10/10 runs
SA with eps = 1e-06:  makespan best = 396, mean = 396.0, reached LB in 10/10 runs

40 個工作:SA 10/10 達到下界 \(\lceil\sum p/m\rceil=396\)(因此是最優),LPT 貪婪(有 \(4/3\) 近似保證)為 398。這個實例加了交換移動後平台問題已不嚴重,\(\varepsilon\) 沒有差別;在只有「移動」鄰域的早期版本中,\(\varepsilon\) 讓平均從 454.2 改善到 453.6——能量塑形是鄰域不夠好時的補救,而非替代。

數獨:表示法 + 重新加熱

表示:每個 \(3\times3\) 區塊內永遠是 \(1..9\) 的排列(提示格固定),所以區塊約束由表示法保證。能量:各列與各行「缺少的數字個數」之和,\(E=0\iff\) 解出。鄰域:同一區塊內交換兩個非提示格,\(\Delta E\) 只需重算兩列兩行。失敗處理:200000 步仍 \(E>0\) 則重新加熱再退一次。

solved 10/10 runs, reheats needed = 0, mean steps to solution = 16300
one solution found:
5 3 4 | 6 7 8 | 9 1 2
6 7 2 | 1 9 5 | 3 4 8
1 9 8 | 3 4 2 | 5 6 7
------+-------+------
8 5 9 | 7 6 1 | 4 2 3
4 2 6 | 8 5 3 | 7 9 1
7 1 3 | 9 2 4 | 8 5 6
------+-------+------
9 6 1 | 5 3 7 | 2 8 4
2 8 7 | 4 1 9 | 6 3 5
3 4 5 | 2 8 6 | 1 7 9

平均 1.6 萬步解出這個經典謎題(Wikipedia 範例)。SA 解數獨不是最有效率的方法(約束傳播 + 回溯快得多),但它展示了一個通用原則:不懂問題的任何解法,只要能定義「有多錯」與「怎麼小改」,退火就能跑。更難的謎題需要更多重新加熱——這也是 SA 在約束滿足問題上的典型行為。

六題總覽

問題 表示 能量 鄰域 約束處理 對照
背包 位元 \(-\)價值 \(+\) 罰分 翻位 / 交換 罰分 DP:差 0.08%
N 皇后 排列 對角衝突 交換兩列 表示法 \(n=1000\) 解出
圖著色 指派 同色邊 換一點色 固定 \(k\) 遞減 DSATUR:少 1 色
Max-Cut 位元 \(-\)割 翻一點 無 窮舉:10/10 最優
排程 指派 \(C_{\max}+\varepsilon\sum l^2\) 移 / 換 能量塑形 下界:10/10
數獨 區塊排列 列行缺數 區塊內交換 表示法 + 重加熱 10/10 解出

變體、親戚與延伸

接受準則的變體

  • 門檻接受(Threshold Accepting,Dueck & Scheuer 1990):接受若 \(\Delta E<\tau_k\),\(\tau_k\downarrow0\)。確定性、不需 \(\exp\) 與亂數,效果與 SA 相當(Ex 5:640 vs 634)。

  • 大洪水(Great Deluge,Dueck 1993):接受若 \(E(y)<\text{水位}_k\),水位隨時間下降。只有一個參數(降水速度)。

  • Record-to-Record Travel:接受若 \(E(y)<E^{\mathrm{best}}+\delta\)。

  • 熱浴 / Glauber 準則:接受機率 \(\dfrac{1}{1+e^{\Delta E/T}}\)(往下也可能拒絕)。同樣滿足細緻平衡,在 Ising 與 Boltzmann 機中常用。

溫度排程的變體

  • 自適應 SA(Ingber 1989, ASA / Very Fast SA):每個維度各自的溫度與步長、\(T_k=T_0e^{-ck^{1/d}}\) 的指數排程配合 Cauchy 式長尾提議,理論上允許比對數快的降溫。

  • 比熱導向:\(\Delta T\propto T/\sqrt{C(T)}\)(Huang et al. 1986)——在比熱峰附近慢、平坦區快。

  • 重新加熱與循環退火:\(T\) 週期性拉回,像金屬的「正火」;對多峰問題有效。

平行回火(Parallel Tempering / Replica Exchange)

同時跑 \(M\) 條不同溫度 \(T_1<\cdots<T_M\) 的鏈,定期以 Metropolis 機率 \[\min\Big\{1,\ \exp\big[(E_i-E_j)\big(\tfrac{1}{T_i}-\tfrac{1}{T_j}\big)\big]\Big\}\] 交換相鄰溫度鏈的狀態。高溫鏈負責跨越大障礙、低溫鏈負責精修,狀態在溫度階梯上下流動。它不降溫(每條鏈維持平穩分布),所以是正確的 MCMC 抽樣器,同時也是強力最佳化器——在自旋玻璃、蛋白質摺疊、貝氏多峰後驗中是標準工具。

量子退火

把熱漲落換成量子穿隧:系統哈密頓量 \(H(s)=(1-s)H_{\text{驅動}}+sH_{\text{問題}}\),\(s\) 從 0 到 1。絕熱定理說變化夠慢則系統停留在基態,終點即問題的最優解。D-Wave 以超導量子位元實作,解 Ising / QUBO 形式的問題(Max-Cut 直接就是 QUBO)。「夠慢」由能隙決定,與 Hajek 定理中 \(d^*\) 扮演類似的角色;量子穿隧對又高又窄的障礙有優勢,熱退火對寬而矮的障礙有優勢。是否有真正的量子加速仍是開放問題。

與其他元啟發式的比較

模擬退火 禁忌搜尋 遺傳演算法 蟻群
搜尋單位 單點 單點 族群 族群(構造式)
逃離局部最優 機率接受變差 記憶禁止回頭 交配 + 突變 費洛蒙蒸發
理論保證 收斂定理(Hajek) 弱 弱(schema 定理) 弱
參數 \(T_0,\alpha,L\) 禁忌長度 族群、交配率、突變率 多
易實作 最易 易 中 中
強項 平滑、鄰域自然的問題 有強結構、確定性 可重組的「基因」 圖上路徑構造

Memetic 組合很常見:GA 的每個子代出生後跑一小段 SA 或局部搜尋(本系列第 13 篇的練習)。No Free Lunch 定理提醒:沒有一種元啟發式在所有問題上最好——SA 的優勢是最少的設計決策、最清楚的理論、最快的原型。

實務應用一覽

  • VLSI 實體設計:Kirkpatrick 的起點。標準單元佈局(TimberWolf)、FPGA 佈局(VPR 至今仍用 SA)、佈線。

  • 排程與時刻表:航班機組、護理師排班、大學課表、鐵路時刻——約束雜、目標多,SA 是 OR 軟體的常客。

  • 影像處理:Geman & Geman 的 MAP 影像復原——以 Markov 隨機場為能量,SA 求最大後驗。

  • 蛋白質結構預測與分子對接:能量 = 力場,鄰域 = 扭轉角擾動;AutoDock 的核心搜尋之一。

  • 機器學習:Boltzmann 機訓練、超參數搜尋、神經架構搜尋的早期方法;模擬退火的「溫度」概念也出現在 softmax 溫度、知識蒸餾與語言模型取樣中。

  • 統計物理:自旋玻璃基態、Ising 模型——SA 的「回娘家」。

  • OR 經典問題:TSP、VRP、二次指派(QAP)、圖劃分、設施選址——幾乎每個 NP-Hard 問題都有 SA 的文獻基準。

完整 C++ 程式

以下為配套程式 annealing.cpp 全文:泛型退火引擎、五種降溫排程、初溫校準、13 個實例(每個都有精確解或強基準對照)。單一檔案、僅用標準函式庫、固定隨機種子。編譯與執行:

g++ -std=c++17 -O2 -Wall -Wextra -o annealing annealing.cpp
./annealing > annealing_output.txt     # 約 35 秒;另產生 annealing_trace.dat
實例 函數 主題 對照基準
1 ex1_hill_vs_sa_1d 一維多峰、\(T_0\) 的影響 格點搜尋全域最小
2 ex2_rastrigin 步長 \(\times\) 排程 已知最小 \(f=0\)
3 ex3_boltzmann 固定 \(T\) 鏈 \(\to\) Boltzmann \(e^{-E/T}/Z\) 解析值
4 ex4_tsp 2-opt、\(\alpha\) 掃描、軌跡檔 Held–Karp、NN+2-opt
5 ex5_schedules 五種排程 + 門檻接受 同預算互比
6 ex6_initial_temperature \(T_0\) 校準 同預算互比
7 ex7_restarts 長鏈 vs. 重啟 同預算互比
8 ex8_knapsack 罰分能量、混合鄰域 DP
9 ex9_nqueens 排列表示、\(O(1)\) \(\Delta E\) 零衝突
10 ex10_coloring 固定 \(k\) 遞減 DSATUR
11 ex11_maxcut 翻點鄰域 \(2^{19}\) 窮舉、貪婪
12 ex12_scheduling 能量塑形、交換鄰域 \(3^{12}\) 窮舉、LPT、下界
13 ex13_sudoku 區塊排列表示、重新加熱 已知唯一解

完整執行結果

以下為 ./annealing 一次執行的完整輸出(Apple clang 21,-O2,約 35 秒)。

所有隨機性來自單一 std::mt19937_64 與固定種子,同一平台重跑結果完全相同。不同標準函式庫的 normal_distribution 實作可能不同,使 Ex 1、2 的數字略異——但每個實驗都報告了多次執行的平均與標準差,結論建立在這些統計量上,不在單一數字上。

練習題

  1. (理論)證明熱浴準則 \(\alpha(x\to y)=\dfrac{1}{1+e^{(E(y)-E(x))/T}}\) 對 Boltzmann 分布滿足細緻平衡。比較它與 Metropolis 準則在 \(\Delta E<0\) 時的差異,哪一個混合較快?

  2. (理論)對 Ex 3 的 8 狀態環,手算每個非全域局部最優的深度與 \(d^*\)。以 \(T_k=c/\ln(k+2)\) 分別取 \(c=d^*/2\) 與 \(c=2d^*\) 模擬 \(10^6\) 步,觀察末端停留在狀態 0 的比例是否符合 Hajek 定理的預測。

  3. (理論)證明定理 2.1 的收斂是指數快的:對 \(x\notin\mathcal{S}^*\),\(\pi_T(x)\le\dfrac{|\mathcal{S}|}{|\mathcal{S}^*|}e^{-(E(x)-E^*)/T}\)。

  4. (設計)把 TSP 的鄰域從 2-opt 改成「Or-opt」(把一段 1–3 個城市的子路徑搬到別處),推導 \(O(1)\) 的 \(\Delta E\),與 2-opt 比較同預算下的解品質。

  5. (設計)實作接受率終止:連續 5 個溫度層接受率 \(<1\%\) 即停。對 Ex 4 的軌跡估算可省下多少預算,把省下的預算用在更慢的 \(\alpha\) 上是否更好?

  6. (排程)實作比熱導向的自適應排程 \(T_{k+1}=T_k\exp\!\big(-\lambda T_k/\sqrt{C(T_k)}\big)\)(Huang et al.),\(C\) 由該層的 \(\operatorname{Var}(E)/T^2\) 估計。與幾何排程同預算比較。

  7. (背包)把罰分係數 \(\lambda\) 改為 \(0.5\max(v_i/w_i)\)(不足以保證可行解占優),觀察 SA 如何被「超重但高價值」的解吸引,以及只記錄可行最佳解的策略還救不救得回來。

  8. (圖著色)把能量改為 Johnson et al. (1991) 的 \(E=-\sum_c|C_c|^2+2\sum_c|B_c|\)(\(C_c\) 為第 \(c\) 色的點集、\(B_c\) 為第 \(c\) 色內的壞邊),它允許 \(k\) 在退火中自動變化。與固定 \(k\) 遞減法比較。

  9. (排程)在 Ex 12 中只保留「移動一個工作」的鄰域,重跑 \(\varepsilon=0\) 與 \(\varepsilon=10^{-6}\),確認能量塑形的效果;再試 \(\varepsilon=10^{-2}\),說明 \(\varepsilon\) 太大會如何改變最優解。

  10. (數獨)找一個「極難」等級的謎題(提示格 \(\le22\)),統計需要幾次重新加熱;嘗試改為 \(T\) 固定在 0.3 的常溫 Metropolis + 記錄最佳,與退火比較。

  11. (平行回火)實作 8 條鏈的平行回火解 Rastrigin 10 維,溫度幾何分布於 \([0.05,20]\),每 100 步嘗試交換相鄰鏈。與單鏈 SA 同預算比較成功率。

  12. (挑戰)二次指派問題(QAP):\(n\) 個設施放到 \(n\) 個位置,成本 \(\sum_{i,j}f_{ij}d_{\sigma(i)\sigma(j)}\)。排列表示、交換鄰域、推導 \(O(n)\) 的 \(\Delta E\)。在 QAPLIB 的 nug12(已知最優 578)上測試。

小結與延伸閱讀

一張表看完全教材

問題 答案 出處
為什麼往上爬是對的 固定 \(T\) 的鏈收斂到 \(\pi_T\propto e^{-E/T}\) 定理 2.2、Ex 3
為什麼降溫會找到最優 \(T\to0\) 時 \(\pi_T\) 集中在 \(\mathcal{S}^*\) 定理 2.1
降溫多慢才保證 \(T_k=c/\ln k\),\(c\ge d^*\)(最深局部最優) Hajek 定理 2.3
為什麼實務不照做 \(e^{100}\) 步、\(d^*\) 未知;用幾何降溫換「高機率好解」 第 2 章
最重要的設計決策 鄰域(連通、對稱、小步、\(O(1)\) \(\Delta E\)) Ex 2、8、12
初溫 接受率 80% 校準;寧高勿低 Ex 6
降溫率 \(\alpha=0.95\)–\(0.99\),報酬遞減 Ex 4c
該在哪個溫度多待 比熱峰附近 Ex 4 軌跡
何時該停 接受率 \(<1\%\) 即凍結 Ex 4 軌跡
約束怎麼處理 表示法 \(>\) 罰分 \(>\) 能量塑形 Ex 8、9、12、13
一條長鏈還是多次重啟 夠長的鏈多跑幾條更穩;太短全毀 Ex 7

三個帶走的觀念

  1. 溫度是「容忍變差的匯率」。\(e^{-\Delta E/T}\) 把能量差換成接受機率;高溫寬容、低溫嚴格。這個單一參數把「探索」與「利用」連成一條連續的路,而不是二選一。

  2. 理論告訴你方向,實驗告訴你數字。Hajek 定理證明「慢」是對的、太快必失敗;但 \(\alpha\)、\(T_0\)、\(L\) 的具體值必須靠實驗——而且每次都要報告平均與標準差,因為 SA 的輸出是隨機變數。

  3. 退火是一種思考方式:遇到卡在局部最優的搜尋,先問「我能不能偶爾接受變差?接受多少?多久之後該收斂?」——從 1953 年的液體模擬到今天的量子退火機,這個問句的答案始終是同一個指數函數。

延伸閱讀

  • Metropolis, Rosenbluth, Rosenbluth, Teller, & Teller (1953). Equation of State Calculations by Fast Computing Machines. J. Chem. Phys.

  • Kirkpatrick, S., Gelatt, C. D., & Vecchi, M. P. (1983). Optimization by Simulated Annealing. Science 220(4598)——開山之作,仍值得全文閱讀。

  • Černý, V. (1985). Thermodynamical Approach to the Traveling Salesman Problem. J. Optim. Theory Appl.

  • Geman, S., & Geman, D. (1984). Stochastic Relaxation, Gibbs Distributions, and the Bayesian Restoration of Images. IEEE PAMI——對數降溫收斂的首次證明。

  • Hajek, B. (1988). Cooling Schedules for Optimal Annealing. Math. Oper. Res. 13(2)——充要條件。

  • Aarts, E., & Korst, J. Simulated Annealing and Boltzmann Machines(Wiley 1989)——理論最完整的專書。

  • van Laarhoven, P. J. M., & Aarts, E. Simulated Annealing: Theory and Applications(1987)。

  • Johnson, D. S., Aragon, C. R., McGeoch, L. A., & Schevon, C. (1989, 1991). Optimization by Simulated Annealing: An Experimental Evaluation, Parts I & II. Oper. Res.——圖劃分、圖著色、數字劃分的大規模實驗,實務調參的聖經。

  • Dueck, G., & Scheuer, T. (1990). Threshold Accepting. J. Comput. Phys.

  • Ingber, L. (1989). Very Fast Simulated Re-annealing. Math. Comput. Modelling.

  • Earl, D. J., & Deem, M. W. (2005). Parallel Tempering: Theory, Applications, and New Perspectives. PCCP.

  • Kadowaki, T., & Nishimori, H. (1998). Quantum Annealing in the Transverse Ising Model. Phys. Rev. E.

  • Henderson, D., Jacobson, S. H., & Johnson, A. W. (2003). The Theory and Practice of Simulated Annealing. Handbook of Metaheuristics——簡潔的綜述。



Algorithms Teaching Series · Advanced · Simulated Annealing