玄学算法——模拟退火

引入

有时我们需要解决一些看似无法解决的问题，比如这题：P1337 [JSOI2004] 平衡点 / 吊打XXX - 洛谷。

总不能把每个坐标都枚举过去吧。（当然这道题也有许多其他优秀的算法。）

这时就需要玄学登场了。

模拟退火

什么是退火？

退火是一种金属热处理工艺，指的是将金属缓慢加热到一定温度，保持足够时间，然后以适宜速度冷却。广义上说，退火是一种对材料的热处理工艺，包括金属材料、非金属材料。而且新材料的退火目的也与传统金属退火存在异同。

——百度百科

模拟退火要做的就是在当前解的基础上 按当前温度的幅度来随机情况 并判断此解的优劣。形式化地说：

有参数 \(T_0\) 表示初始温度，\(T_k\) 表示终止温度，\(d\) 表示降温系数。\(S\) 表示当前最优解，\(T\) 表示当前温度。注意 \(T_0>T_k\)（即整个退火是降温操作）。

每次通过随机化生成一个解 \(S'\)，那么接受这个解的概率为：

\[\begin{cases} 1&\text{S' is better than S}\\ e^{-\frac{\Delta}{T}}&\text{otherwise} \end{cases} \]

其中 \(\Delta\) 表示两个解的差。由于 \(e^{-x}<1~(x>0)\) 且递减，所以 \(\Delta\) 越大（解越差），接受概率越小。同时 \(T\) 越大，接受概率越大（更狂、更不理智、更赌）。（模拟退火含金量最高的地方。好像与玻尔兹曼分布有关。其实不用了解太多。）

那么怎么退火？

• \(T\gets T_0\)。
• 随机化：
  • 随机解并按概率接受解。随机解时应按 \(T\) 的大小控制随机幅度（具体看代码）。
  • \(T\gets T\times d\)。
• 重复上一步（随机化），如果 \(T<T_k\)，算法结束。

#include <iostream>
#include <cmath>
#include <ctime>
using namespace std;
 
typedef double d;
 
// 分别为：降温系数、开始温度、结束温度、用于卡时的。
const d down=0.997,Beg=1e5,End=1e-15,MAX_TIME=0.8;
// 降温系数一般接近 1（徐徐降温），开始温度一般足够大（看数据范围），结束温度足够小（看答案精度），卡时常数略小于时限（单位 秒）。
 
const int N=2145;
int n;
d ansx,ansy,answ;
 
struct node
{
	int x,y,w;
};
 
node a[N];
 
d calc(d x,d y)// 计算解的优劣
{
	d r=0,dx,dy;
	for(int i=1;i<=n;i++)
	{
		dx=x-a[i].x;
		dy=y-a[i].y;
		r+=sqrt(dx*dx+dy*dy)*a[i].w;
	}
	return r;
}
 
void SA()
{
	d t=Beg;
	while(t>End)
	{
		d ex=ansx+(rand()*2-RAND_MAX)*t,ey=ansy+(rand()*2-RAND_MAX)*t;// 按当前温度随机化解
		// printf("%.5lf %.5lf\n",ex,ey);
		d ew=calc(ex,ey);// 计算优劣
		d de=ew-answ;// 计算优劣差异
		if(de<0)// 随机的解更好
		{
			ansx=ex,ansy=ey,answ=ew;// 接受
		}
		else if(exp(-de/t)*RAND_MAX>rand())// 否则以概率接受
		{
			ansx=ex,ansy=ey;
		}
		t*=down;// 降温
	}
}
 
inline void init()
{
	srand(time(NULL));
	srand(rand()+114514);
	srand(rand()+1919810);
	srand(rand());
}
 
int main()
{
	init();
	scanf("%d",&n);
	for(int i=1;i<=n;i++)
	{
		scanf("%d %d %d",&a[i].x,&a[i].y,&a[i].w);
		ansx+=a[i].x;
		ansy+=a[i].y;
	}
	ansx/=n,ansy/=n;// 以平均数为初始解
	answ=calc(ansx,ansy);
	while((double)clock()/CLOCKS_PER_SEC<MAX_TIME)SA();// 卡时，有时间就一直模拟退火
	printf("%.3lf %.3lf",ansx,ansy);// 输出解
	return 0;
}

卡时寄巧

定一个略小于时限的卡时常数（必须满足：卡时常数 \(+\) 单次模拟退火用时 \(<\) 时限）。一般定 \(0.8\sim0.9\)（时限 \(1~\text s\)）。

这样只要还有时间，就一直跑模拟退火（增加找到解的概率）。

后记

初学模拟退火，以上就是我的心得体会。以后有新东西会补上。

EOF