模拟退火学习记录

模拟退火学习记录

没错，这篇博客连笔记都算不上，因为我还没有确切地理解这个算法的奥妙，但是略有所感，就想写下来。争取用最简单的语言描述下来，让小白能容易地入门。下面省略无数个“个人认为”……

一、简介

模拟退火是一种随机化算法。当一个问题的方案数量极大（甚至是无穷的）而且不是一个单峰函数时，我们常使用模拟退火求解
——摘自OI WIKI

这个算法是将一个物理过程的模拟转化成一个数学模型，然后应用到一些不易解决的问题。由于是随机化算法，所以需要多次执行退火过程来确保最优解的质量更高。

---

二、前置知识

（1）模拟退火需要3个量：起始温度T0，降温系数d，终止温度Tk。每个时刻的温度决定着枚举下一个状态的随机化程度，降温系数决定了降温的快慢。
（2）这里还需要引入能量的概念来刻画与目标解的差距 。能量越大，与最优答案相差越远；能量越小，越接近最优答案。
相应地，能量差 便表示 当前随机化枚举出来的答案，相对于前一个答案而言，与最优解的接近程度。能量差为负，表示更接近最优解，反之则更远离；能量差的绝对值大小表示接近程度的差距大小。
注意，下面用D表示能量差，E表示能量。
（3）我们对于新枚举的一个状态，该怎么决定是否跳到新的状态？这个判断依据非常奇特：如果D<0，则跳到该状态；反之则以P=e^-D/T 的概率跳到新状态。
引用OI WIKI的一个比喻：有一个兔子喝醉了，它会尽可能地朝山顶跳，它有可能跳到山顶，但是也有可能一下子越过山顶。随后兔子慢慢清醒，跳跃也更加冷静。
不难发现，我们在不断地朝着最优解的方向行进，虽然时有偏差，但总趋势是不断逼近最优解的。

---

三、流程

退火过程的框架如下：

• 设定T0,Tk,d，通常d的选取在[0.985~0.999]范围内选取。
• 设计随机化函数得到目标状态，并计算当前状态与目标状态的能量差，通过上文介绍的判断方法决定是否跳到下一状态
• 令T=T*d,重复上述过程，直至T<Tk,此时算法结束

就像这样：

for(double T=你给定的值,d=你给定的值;T>你给定的值;T*=d){
	随机化得到新的状态
	计算出能量差
	比较，如果更优则接受，更劣则一定概率接受 
} 得到答案 

四、例题

干说还是有点空洞难懂，这里举一道例题

P1337 [JSOI2004]平衡点 / 吊打XXX
这道题需要一个前置知识：绳结最终稳定在∑d[i]w[i]最小的位置，其中d[i]表示绳结距离i号点的欧氏距离。我们用模拟退火解决。绳结的状态可用(x,y)表示，那么我们随机化地找到新的状态(ex,ey)并计算出此时的能量ew=∑d[i]w[i]，计算能量差，比较即可。值得注意的是，为了降低枚举量，我们取各个点x值、y值的平均值作为初始状态中的x,y。

#include<bits/stdc++.h>
using namespace std;
#define il inline
#define re register
il int read(){
	int s=0,w=1;char c=getchar();
	while(c<'0'||c>'9'){ if(c=='-') w=-1;c=getchar();}
	while(c>='0'&&c<='9'){ s=(s<<1)+(s<<3)+c-'0';c=getchar();}
	return s*w;
} 
const int N=1e4+10;
const double Down=0.996;
int n,x[N],y[N],w[N];
double ansx,ansy,answ;
il double Rand(){ return rand()*2-RAND_MAX; }
//生成[-RAND_MAX,RAND_MAX]范围内的随机数 
il double calc(double xx,double yy){//计算新状态能量值 
	double res=0,dx,dy;
	for(re int i=1;i<=n;++i){
		dx=x[i]-xx,dy=y[i]-yy;
		res+=sqrt(dx*dx+dy*dy)*w[i];
	}
	return res;
}
il void sa(){
	for(re double T=8000,ex,ey,ew,de;T>=1e-15;T*=Down){
		ex=ansx+Rand()*T,ey=ansy+Rand()*T,ew=calc(ex,ey),de=ew-answ;
		//(ex,ey):随机化枚举出的下一状态
		//ew:由(ex,ey)确定的下一状态的能量值
		//de：能量差，初态-末态 
		if(de<0){//如果此答案让能量更低，就接受 
			ansx=ex,ansy=ey,answ=ew;
		} 
		else if(exp(-de/T)*RAND_MAX>rand()){//防止精度误差，除法变乘法 
			ansx=ex,ansy=ey; 
		}
	}
}
il void solve(){
	sa();sa();sa();sa();
}
int main()
{
	srand(0);//实际上不应该选择确定的种子，但是我看到题解用这个种子直接就过了
	//实际上应该是这个：srand(time(0)) 
	n=read();
	for(re int i=1;i<=n;++i)
		x[i]=read(),y[i]=read(),w[i]=read(),ansx+=x[i],ansy+=y[i];
	ansx/=n,ansy/=n,answ=calc(ansx,ansy);
	solve();
	printf("%.3lf %.3lf",ansx,ansy);
	return 0;
	
}