模拟退火简记

概述#

模拟退火是一种随机算法，一般用于最优化状态的问题中，主要思路是随机调整状态，如果调整后的状态更优则接受，否则以一定概率接受。

设状态 $s$ 得到的答案为 $C(s)$，我们需要找到最小的答案。

设置初始状态 $s_0$，初始温度 $T_0$（一般取在 $[2000,3000]$ 中，可以尽量取大），温度变化量 $\Delta T$（一般取 $0.9996$）。

设当前状态为 $s$，温度为 $T$。随机调整 $s$，随机值固定时调整幅度需要与温度成正比，得出新的状态 $s'$。然后与已知的最优解 $w$ 比较。如果 $C(s')<w$，直接接受（将 $s$ 赋值为 $s'$），并将 $w$ 赋值为 $C(s')$；否则以 $\exp(\frac{w-C(s')}{T})$ 的概率接受。

重复以上步骤，每做一次都要降温（$T\gets T\cdot\Delta T$），直到温度接近 $0$（在代码中表现为小于 eps）。

下面通过两个例题来理解。

例题#

UVA10228 A Star not a Tree?#

题目大意#

平面上有 $n$ 个点，求出到所有点距离和最小的点到所有点的距离和。

题解#

模拟退火，状态为平面上的一个点，答案函数为状态到所有点的距离和，初始状态可以设置为所有点两个坐标的平均值。随机调整时，可以首先随机一个 $[-7,7]$ 之间的数（或者别的什么，只要保证每个点都可能跑到），再乘以温度，作为坐标的变化量。

参考代码#

#include<bits/stdc++.h>
using namespace std;
typedef double lf;
constexpr int N=105;
constexpr lf eps=1e-8,T0=3000,Td=0.9996;
mt19937 gen(chrono::system_clock::now().time_since_epoch().count());
lf rd(lf l,lf r){return uniform_real_distribution<>(l,r)(gen);}
lf dis(lf X,lf Y,lf x,lf y){return sqrt((X-x)*(X-x)+(Y-y)*(Y-y));}
int n,x[N],y[N];
lf f(lf X,lf Y){
	lf ret=0;for(int i=1;i<=n;i++)ret+=dis(X,Y,x[i],y[i]);
	return ret;
}
lf SA(lf X,lf Y){
	lf x=X,y=Y,s=f(X,Y);for(lf T=T0;T>=eps;T*=Td){
		lf nx=x+rd(-7,7)*T,ny=y+rd(-7,7)*T,u=f(nx,ny);
		if(u<s)x=nx,y=ny,s=u;else if(rd(0,1)<=exp((s-u)/T))x=nx,y=ny;
	}
	return s;
};
int main(){
	ios::sync_with_stdio(false);
	cin.tie(nullptr);
	int T;cin>>T;while(T--){
		cin>>n;lf X=0,Y=0;
		for(int i=1;i<=n;i++){
			cin>>x[i]>>y[i];
			X+=x[i];Y+=y[i];
		}
		X/=n;Y/=n;
		printf("%.0lf\n",SA(X,Y));
		if(T)putchar('\n');
	}
	return 0;
}

P5544 [JSOI2016]炸弹攻击1#

题目大意#

平面上有 $N$ 个圆，$M$ 个点。找到一个圆心坐标任意，半径不超过 $R$，且与那 $N$ 个圆不相交（可以相切）的圆，求这个圆最多能够覆盖 $M$ 个点中的多少个。

题解#

在保证不会超时的情况下，可以适当把初始温度和温度变化量调高一点，并且多跑几次模拟退火，增大得到最优解的概率。

这道题要求的是最大值，这个时候当新解不优时应该以 $1-\exp(\frac{w-C(s')}{T})$ 的概率接受。

参考代码#

#include<bits/stdc++.h>
using namespace std;
typedef double lf;
constexpr int N=1e3;
constexpr lf eps=1e-8,T0=2500,Td=0.9997;
mt19937 gen(chrono::system_clock::now().time_since_epoch().count());
lf rd(lf l,lf r){return uniform_real_distribution<>(l,r)(gen);}
lf dis(lf X,lf Y,lf x,lf y){
	return sqrt((X-x)*(X-x)+(Y-y)*(Y-y));
}
int n,m,R,x[N],y[N],r[N],p[N],q[N];
int f(lf X,lf Y){
	lf rr=R;int ret=0;
	for(int i=0;i<n;i++){
		rr=min(rr,dis(X,Y,x[i],y[i])-r[i]);
	}
	if(rr<0)return 0;
	for(int i=0;i<m;i++){
		if(dis(X,Y,p[i],q[i])<rr+eps)ret++;
	}
	return ret;
}
int SA(lf x,lf y){
	int s=f(x,y);for(lf T=T0;T>eps;T*=Td){
		lf nx=x+rd(-10,10)*T,ny=y+rd(-10,10)*T;
		int u=f(nx,ny);if(u>s)x=nx,y=ny,s=u;
		else if(rd(0,1)>exp((s-u)/T))x=nx,y=ny;
	}
	return s;
}
int main(){
	ios::sync_with_stdio(false);
	cin.tie(nullptr);
	cin>>n>>m>>R;
	for(int i=0;i<n;i++){
		cin>>x[i]>>y[i]>>r[i];
	}
	lf X,Y;
	for(int i=0;i<m;i++){
		cin>>p[i]>>q[i];
		X+=p[i];Y+=q[i];
	}
	X/=m;Y/=m;int ans=0;auto t=clock();
	for(int i=0;i<5;i++)ans=max(ans,SA(X,Y));
	printf("%d\n",ans);
	return 0;
}