模拟退火(simulated annealing)算法详解
模拟退火算法来源于固体退火原理,得益于材料统计力学的研究成果,并且该算法也是一种基于概率的算法。该算法主要用于求解最优解问题,如巡航问题、函数极值等。
在材料统计力学中,材料中的粒子的不同结构对应粒子的不同能量水平,高温条件下粒子能量高,可以自由运动和重新排列;低温条件下粒子能量低。如果从高温开始缓慢降温(即退火),粒子就可以在每个温度下达到热平衡。当系统完全被冷却时,最终形成低能状态晶体。
算法原理
设材料在状态i的能量为E(i),那么材料在温度T时从状态i进入状态j的规律如下:
(1)如果E(j)<=E(i):接受该状态j
(2)如果E(j)>E(i),则以以下概率被接受。其中K是玻尔兹曼常数,算法中我们通常让K=1。T是材料温度,算法中我们通常设初始温度T=1。
当温度降至很低时,我们认为材料以很大概率进入能量最小状态。
算法流程
(1)选任一状态s0作为初始解s(0)=s0,设初始温度为T0,T0可以是1;设终止温度e,可以让e=0.1^30;设降温速度alpha,可以是0.999,迭代次数M,迭代次数,次数足够就可以,初始时i=0。
(2)i=0时当前解s(0)=s0,T0=1。每次迭代后T(i)=alpha*T(i-1)。并且按照某一规定方式(例如随机产生数)产生当前解s(i)。
(3)取评价函数C
如果该式小于零,就接受s(i)作为新的解;如果大于零,就以概率:
接受s(i)作为新的解。
(4)令i=i+1,进行下一次迭代。如果i达到迭代次数M或温度T(i)已经小于最终温度e,则停止迭代。最终得到结果s(i)。
算法实现
本文以求解函数最小值为例。设函数为f(x)=x^2+2*x+1。我们求函数在区间(0,10)的极小值点和极小值。
class SAA:
'''模拟退火算法'''
def __init__(self,T=1,e=0.1**30,alpha=0.999,M=5000):
self.T=T #初始温度,这里设为1
self.e=e #终止温度,这里设为0.1*830
self.alpha=alpha #alpha为降温速度
self.M=M #M为迭代次数,这里是5000次
def f1(self,x):
'''所需要求的函数'''
return x**2+2*x+1
def main(self):
resule_list=[] #存储每次迭代后的解
x1=np.random.rand()*10 #在区间0到10随机选取初始值
for k in range(self.M): #迭代
x=np.random.rand()*10 #随机选取新的值
dx=self.f1(x)-self.f1(x1) #公式里面的C(s(i))-C(s(i-1))
if dx<0:
x1=x #接受新解
elif dx>0:
if np.exp(-dx/self.T)>=np.random.rand():
x1=x #以概率接受新的解
self.T*=self.alpha #降温
if self.T<self.e: #判断是否降到所设定的温度
break
resule_list.append(x1)
return resule_list #返回每次迭代后的结果
if __name__=="__main__":
SAA=SAA()
result_list=SAA.main()
M=range(5000)
plt.plot(M,result_list,color='red')
plt.plot(M,SAA.f1(np.array(result_list)),color='green')
plt.show()
print(result_list)
结果如下图所示:
图中红色的是x值随迭代次数的变化,绿色是函数值随迭代次数的变化。我们发现迭代3000次之后结果趋于稳定。函数达到最小值,而红色的极小值点也趋近于1。如果我们把区间扩大,比如把-1包含进去,那么最小值就将趋近于0,极小值点趋近于-1。。
总结
模拟退火算法相对其它算法而言比较容易实现,理论上该算法具有概率的全局优化性能,可以趋于全局最优,并且在诸多领域得到广泛应用,如机器学习、神经网络、信号处理等。使用的时候需要学会灵活使用其中的参数。例如:温度下降过快,就不一定得到最优解,下降过慢又会增加运行时间,终止温度也不宜过高或过低,迭代次数适当等等。