- 抽象來源:模仿冶金過程中的退火原理。
- 核心思想:在冶金退火過程中,隨着溫度的下降,系統內部分子的平均動能逐漸降低,分子在自身位置附近的擾動能力也隨之下降,即分子自身的搜索范圍隨着溫度的下降而下降。利用該特性,我們可以對給定狀態空間(待求解空間)內的某個狀態產生函數(待求解函數)的最值進行求解。在高溫狀態下,由於分子的擾動能力較強,對較差狀態(遠離最值所對應的狀態)的容忍性高,因此可以在給定狀態空間內進行全局的隨機搜索,從而有較高概率跳出局部極值。隨着溫度的逐漸下降,分子的擾動能力減弱,對較差狀態的容忍性隨之降低,導致此時的全局隨機搜索能力下降,相應地對局部極值的搜索能力上升。綜合整個退火過程,在理想情況下,最終解應該對應於給定狀態空間內的最值。
- 迭代公式:
系統溫度為$T$時,出現能量差為$dE$的降溫概率為:
\begin{equation}\label{eq_1}
p(dE) = e^{dE/kT}
\end{equation}
其中,$dE<0$,因此$0<p(dE)<1$。
以該降溫概率為基礎,采用狀態轉移概率$p(\Delta f)$來表示對較差狀態的容忍性:
\begin{equation}\label{eq_2}
p(\Delta f) = e^{\Delta f/T}
\end{equation}
其中,$f$為狀態產生函數,因此$\Delta f = f(x_{new}) - f(x_{old})$。常數$k$通過改變溫度$T$的取值范圍而被忽略,即式\ref{eq_1}中的$kT$等效於式\ref{eq_2}中的$T$。
注意,在求取最值的過程中,存在兩類可接受解:1)更優解 --- $BS$(較當前狀態更加接近最值的狀態);2)容忍解 --- $TS$(較當前狀態更加遠離最值的解)。前者的接受概率1;后者的接受概率為狀態轉移概率$p(\Delta f)$,且必須保證$p(\Delta f)\in (0, 1)$。- 求取最小值:
\begin{equation}
\begin{cases}
\Delta f < 0 \rightarrow BS \rightarrow 1\\
\Delta f \geq 0 \rightarrow TS \rightarrow p(-\Delta f)
\end{cases}
\end{equation} - 求取最大值:
\begin{equation}
\begin{cases}
\Delta f > 0 \rightarrow BS \rightarrow 1\\
\Delta f \leq 0 \rightarrow TS \rightarrow p(\Delta f)
\end{cases}
\end{equation}
- 求取最小值:
- Python代碼實現:
1 import numpy as np 2 import matplotlib.pyplot as plt 3 import random 4 5 class SA(object): 6 7 def __init__(self, interval, tab='min', T_max=10000, T_min=1, iterMax=1000, rate=0.95): 8 self.interval = interval # 給定狀態空間 - 即待求解空間 9 self.T_max = T_max # 初始退火溫度 - 溫度上限 10 self.T_min = T_min # 截止退火溫度 - 溫度下限 11 self.iterMax = iterMax # 定溫內部迭代次數 12 self.rate = rate # 退火降溫速度 13 ############################################################# 14 self.x_seed = random.uniform(interval[0], interval[1]) # 解空間內的種子 15 self.tab = tab.strip() # 求解最大值還是最小值的標簽: 'min' - 最小值;'max' - 最大值 16 ############################################################# 17 self.solve() # 完成主體的求解過程 18 self.display() # 數據可視化展示 19 20 def solve(self): 21 temp = 'deal_' + self.tab # 采用反射方法提取對應的函數 22 if hasattr(self, temp): 23 deal = getattr(self, temp) 24 else: 25 exit('>>>tab標簽傳參有誤:"min"|"max"<<<') 26 x1 = self.x_seed 27 T = self.T_max 28 while T >= self.T_min: 29 for i in range(self.iterMax): 30 f1 = self.func(x1) 31 delta_x = random.random() * 2 - 1 32 if x1 + delta_x >= self.interval[0] and x1 + delta_x <= self.interval[1]: # 將隨機解束縛在給定狀態空間內 33 x2 = x1 + delta_x 34 else: 35 x2 = x1 - delta_x 36 f2 = self.func(x2) 37 delta_f = f2 - f1 38 x1 = deal(x1, x2, delta_f, T) 39 T *= self.rate 40 self.x_solu = x1 # 提取最終退火解 41 42 def func(self, x): # 狀態產生函數 - 即待求解函數 43 value = np.sin(x**2) * (x**2 - 5*x) 44 return value 45 46 def p_min(self, delta, T): # 計算最小值時,容忍解的狀態遷移概率 47 probability = np.exp(-delta/T) 48 return probability 49 50 def p_max(self, delta, T): 51 probability = np.exp(delta/T) # 計算最大值時,容忍解的狀態遷移概率 52 return probability 53 54 def deal_min(self, x1, x2, delta, T): 55 if delta < 0: # 更優解 56 return x2 57 else: # 容忍解 58 P = self.p_min(delta, T) 59 if P > random.random(): return x2 60 else: return x1 61 62 def deal_max(self, x1, x2, delta, T): 63 if delta > 0: # 更優解 64 return x2 65 else: # 容忍解 66 P = self.p_max(delta, T) 67 if P > random.random(): return x2 68 else: return x1 69 70 def display(self): 71 print('seed: {}\nsolution: {}'.format(self.x_seed, self.x_solu)) 72 plt.figure(figsize=(6, 4)) 73 x = np.linspace(self.interval[0], self.interval[1], 300) 74 y = self.func(x) 75 plt.plot(x, y, 'g-', label='function') 76 plt.plot(self.x_seed, self.func(self.x_seed), 'bo', label='seed') 77 plt.plot(self.x_solu, self.func(self.x_solu), 'r*', label='solution') 78 plt.title('solution = {}'.format(self.x_solu)) 79 plt.xlabel('x') 80 plt.ylabel('y') 81 plt.legend() 82 plt.savefig('SA.png', dpi=500) 83 plt.show() 84 plt.close() 85 86 87 if __name__ == '__main__': 88 SA([-5, 5], 'max')
筆者所用示例函數為 :
\begin{equation}
f(x) = (x^2 - 5x)sin(x^2)
\end{equation}
- 結果展示:
- 參考:https://blog.csdn.net/AI_BigData_wh/article/details/77943787?locationNum=2&fps=1