模拟退火算法(simulated annealing,SA)来源于固体退火原理,是一种基于概率的算法。
模拟退火算法(SA)来源于固体退火原理,是一种基于概率的算法。将固体加温至充分高的温度,再让其徐徐冷却,加温时,固体内部粒子随温升变为无序状,内能增大,分子和原子越不稳定。而徐徐冷却时粒子渐趋有序,能量减少,原子越稳定。在冷却(降温)过程中,固体在每个温度都达到平衡态,最后在常温时达到基态,内能减为最小。
模拟退火算法从某一较高初温出发,伴随温度参数的不断下降,结合概率突跳特性在解空间中随机寻找目标函数的全局最优解,即在局部最优解能概率性地跳出并最终趋于全局最优。模拟退火算法是通过赋予搜索过程一种时变且最终趋于零的概率突跳性,从而可有效避免陷入局部极小并最终趋于全局最优的串行结构的优化算法。
(相关资料图)
模拟退火算法从某一较高初温出发,伴随温度参数的不断下降,结合一定的概率突跳特性在解空间中随机寻找目标函数的全局最优解,即在局部最优解能概率性地跳出并最终趋于全局最优。
这里的 “一定的概率” 的计算参考了金属冶炼的退火过程,这也是模拟退火算法名称的由来。将温度 T 当作控制参数,目标函数值 f 视为内能 E ,而固体在某温度 T 时的一个状态对应一个解 \(x_i\),然后算法试图随着控制参数 T 的降低,使目标函数 f (内能 E )也逐渐降低,直至趋于全局最小值(退火中低温时的最低能量状态),就像金属退火过程一样。
关于爬山算法与模拟退火,有一个有趣的比喻:
爬山算法:兔子朝着比现在高的地方跳去。它找到了不远处的最高山峰。但是这座山不一定是珠穆朗玛峰。这就是爬山算法,它不能保证局部最优值就是全局最优值。
模拟退火:兔子喝醉了。它随机地跳了很长时间。这期间,它可能走向高处,也可能踏入平地。但是,它渐渐清醒了并朝最高方向跳去。这就是模拟退火。
从上面我们知道,会结合概率突跳特性在解空间中随机寻找目标函数的全局最优解,那么具体的更新解的机制是什么呢?如果新解比当前解更优,则接受新解,否则基于Metropolis
准则判断是否接受新解。接受概率为:
假设开始状态在A,随着迭代次数更新到B局部最优解,这时发现更新到B时,能量比A要低,则说明接近最优解了,因此百分百转移,状态到达B后,发现下一步能量上升了,如果是梯度下降则是不允许继续向前的,而这里会以一定的概率跳出这个坑,这个概率和当前的状态、能量等都有关系,如果B最终跳出来了到达C,又会继续以一定的概率跳出来,直到到达D后,就会稳定下来。
(1) 初始化:初始温度 T (充分大),初始解状态S(是算法迭代的起点),每个 T 值的迭代次数 L
(2) 对 k=1, …, L 做第(3)至第(6)步:
(3) 产生新解 S′
(4) 计算增量 ΔT = C(S′)-C(S),其中 C(S) 为目标函数, C(S) 相当于能量
(5) 若 ΔT<0 则接受 S′ 作为新的当前解,否则以概率 exp(-ΔT/T) 接受S′作为新的当前解.
(6) 如果满足终止条件则输出当前解作为最优解,结束程序。
(7) T 逐渐减少,且 T->0 ,然后转第2步。
其中有几个需要注意的点:
可以大概分为这四个步骤:
# -*- coding: utf-8 -*-"""Created on Mon Apr 3 19:17:28 2023@author: steve"""from random import randomimport mathimport matplotlib.pyplot as pltmax_iter = 100 # 每一次温度降低的迭代次数alpha = 0.99 # 降温系数T_f = 0.01 # 温度的终值T_n = 100 # 当前的温度,也是初始温度x, y = [random() * 10 - 5 for i in range(max_iter)], [random() * 10 - 5 for i in range(max_iter)] # 进行数据的初始化f = lambda x, y : (4 * x ** 2 - 2.1 * x ** 4 + x ** 6 / 3 + x * y - 4 * y ** 2 + 4 * y ** 4) # 我们需要求的函数result = { "f": [], "t": [] } # 用来存放每一次下降的最优解
def metrospolis(T_n, old, new): """ 进行 Metrospolis 准则的判断 Parameters ---------- T_n : int 当前的温度. old : double 函数扰动之前的值. new : double 函数扰动之后的值. Returns ------- int 是否需要重新寻找值. """ if old >= new: return 1 else: p = math.exp((old - new) / T_n) if random() < p: return 1 else: return 0
def generate_new(T_n, x, y): """ 其为扰动的过程 Parameters ---------- T_n : int 当前的温度. x : double DESCRIPTION. y : double DESCRIPTION. Returns ------- list 返回生成的新值. """ while True: x_new = x + T_n * (random() - random()) y_new = y + T_n * (random() - random()) if (-5 <= x_new <= 5) & (-5 <= y_new <= 5): break #重复得到新解,直到产生的新解满足约束条件 return x_new, y_new
def best(max_iter, f, x, y): """ 计算这个温度下的最优值 Parameters ---------- max_iter : int 最大的迭代次数. f : function 我们需要求的函数. x : double DESCRIPTION. y : double DESCRIPTION. Returns ------- list 返回最优值,以及它的索引. """ min_val, min_inx = f(x[0], y[0]), 0 for i in range(max_iter): val = f(x[i], y[i]) if min_val > val: min_val, min_inx = val, i return [min_val, min_inx]
def plot(result): plt.plot(result["t"], result["f"]) plt.title("SA") plt.xlabel("t") plt.ylabel("f") plt.gca().invert_xaxis() plt.show()def main(max_iter, alpha, T_f, T_n, x, y, f, result): count = 0 # 统计迭代了多少次 while T_n > T_f: # 外层循环,当前温度小于最低温度时,终止循环 for i in range(max_iter): # 内层循环 x_new, y_new = generate_new(T_n, x[i], y[i]) # 产生新值 if metrospolis(T_n, f(x[i], y[i]), f(x_new, y_new)): # 将原值与新产生的值进行比较 x[i] = x_new # 如果接收新值,则存入数组中 y[i] = y_new # 迭代 max_iter 后,记录该温度下的最优解 [ft, _] = best(max_iter, f, x, y) result["f"].append(ft) result["t"].append(T_n) T_n = T_n * alpha # 进行降温操作 count += 1 # 得到最优解 [f_best, inx] = best(max_iter, f, x, y) print(f"F={f_best}, x_1={x[inx]}, x_2={y[inx]}") # 进行图像表示 plot(result)
# -*- coding: utf-8 -*-"""Created on Mon Apr 3 19:17:28 2023@author: steve"""from random import randomimport mathimport matplotlib.pyplot as pltdef metrospolis(T_n, old, new): """ 进行 Metrospolis 准则的判断 Parameters ---------- T_n : int 当前的温度. old : double 函数扰动之前的值. new : double 函数扰动之后的值. Returns ------- int 是否需要重新寻找值. """ if old >= new: return 1 else: p = math.exp((old - new) / T_n) if random() < p: return 1 else: return 0 def generate_new(T_n, x, y): """ 其为扰动的过程 Parameters ---------- T_n : int 当前的温度. x : double DESCRIPTION. y : double DESCRIPTION. Returns ------- list 返回生成的新值. """ while True: x_new = x + T_n * (random() - random()) y_new = y + T_n * (random() - random()) if (-5 <= x_new <= 5) & (-5 <= y_new <= 5): break #重复得到新解,直到产生的新解满足约束条件 return x_new, y_new def best(max_iter, f, x, y): """ 计算这个温度下的最优值 Parameters ---------- max_iter : int 最大的迭代次数. f : function 我们需要求的函数. x : double DESCRIPTION. y : double DESCRIPTION. Returns ------- list 返回最优值,以及它的索引. """ min_val, min_inx = f(x[0], y[0]), 0 for i in range(max_iter): val = f(x[i], y[i]) if min_val > val: min_val, min_inx = val, i return [min_val, min_inx]def plot(result): plt.plot(result["t"], result["f"]) plt.title("SA") plt.xlabel("t") plt.ylabel("f") plt.gca().invert_xaxis() plt.show()def main(max_iter, alpha, T_f, T_n, x, y, f, result): count = 0 # 统计迭代了多少次 while T_n > T_f: # 外层循环,当前温度小于最低温度时,终止循环 for i in range(max_iter): # 内层循环 x_new, y_new = generate_new(T_n, x[i], y[i]) # 产生新值 if metrospolis(T_n, f(x[i], y[i]), f(x_new, y_new)): # 将原值与新产生的值进行比较 x[i] = x_new # 如果接收新值,则存入数组中 y[i] = y_new # 迭代 max_iter 后,记录该温度下的最优解 [ft, _] = best(max_iter, f, x, y) result["f"].append(ft) result["t"].append(T_n) T_n = T_n * alpha # 进行降温操作 count += 1 # 得到最优解 [f_best, inx] = best(max_iter, f, x, y) print(f"F={f_best}, x_1={x[inx]}, x_2={y[inx]}") # 进行图像表示 plot(result) max_iter = 100 # 每一次温度降低的迭代次数alpha = 0.99 # 降温系数T_f = 0.01 # 温度的终值T_n = 100 # 当前的温度,也是初始温度x, y = [random() * 10 - 5 for i in range(max_iter)], [random() * 10 - 5 for i in range(max_iter)] # 进行数据的初始化f = lambda x, y : (4 * x ** 2 - 2.1 * x ** 4 + x ** 6 / 3 + x * y - 4 * y ** 2 + 4 * y ** 4) # 我们需要求的函数result = { "f": [], "t": [] } # 用来存放每一次下降的最优解main(max_iter, alpha, T_f, T_n, x, y, f, result)
最后的运行结果为:
关键词: