关于模拟退火

感谢 caq 的友情讲解 orz

模拟退火是一种靠随机化解决最优解问题的一种方法。

它通常是设定一个初温、降温系数以及极值，在不断降温的过程中随机找最优解。

对于每次随机到的一个新的方案，如果它的答案比当前更优，那么我们更新当前答案；

如果不比当前更优，则我们有一定的几率更新答案。这个几率的计算也有特定公式。

以下设初温为 (T)，降温系数为 (d)，当前最优解为 ( ext{ans})，随机到的新答案为 ( ext{new})。

那么 (e^{Large{frac{ ext{ans}- ext{new}}{T imes k}}}) 就是我们在遇到非最优解时的更新概率。其中 (k) 大多情况下为 (1)。

具体流程如图所示：

假设我们求的是最小值，那么概率接受答案的代码如下：

if(exp(-del/T)>(double)rand()/RAND_MAX) ... ; //如果满足条件进行更新。

因为是随机化算法，要运行多次才能使答案尽量优。

为了使答案能更优，我们要尽可能地多去进行随机比较，所以我们可以这么写：

while((double)clock()/CLOCKS_PER_SEC<0.9) SA();//表示在 0.9s 内不断进行随机比较

随机比较的效果也取决于我们选择的参数大小。

通常情况下，降温系数是个 (< 1) 但尽量接近于 (1) 的实数，最低极值是个 (>0) 但尽量接近于 (0) 的实数，初温是个比较大的整数。

随着温度的不断降低，答案的更新频率也越来越低。

具体降温的模板如下：

void SA(){
  double T=2021;//初温
  while(T>lim){//lim 为极值，这里表示若温度够则不断降温
    ...//此处搭配具体题目写的计算函数进行随机答案
    del=now-ans;//
    if(del<0) ans=now;//这里求的是最小值，若更优则直接更换
    else if(exp(-del/T)>(double)rand()/RAND_MAX) ... ;
    //若有几率更换则进行相应操作（不能直解赋值更新，因为它不优）
    T*=d;//d 为降温系数，这里表示不断降温
  }
}

当然了，毕竟是随机化算法，所以有的时候要多交几遍才能 AC 也是很正常的事。

所以调参数，随 rp 也是模拟退火十分重要的一点（

这个时候稳住心态才是最重要的（

当然了，模拟退火本身也不是很难，所以我们不扯别的概念啥的了，直接上几个典型题就懂了。

[]

UVA10228 A Star not a Tree?

首先明确，我们要随机的是点，即一个横坐标和一个纵坐标。

我们就通过不断随机这个点，来暴力计算到每一个顶点的距离，并按上面说的方法更新答案即可。

当然我们也可以对此进行优化。我们运用人类智慧，可以猜测这个点在多边形中间位置。

当然我们没法直接找出来，所以可以先把每个点的横纵坐标分别加起来，然后除以 (n)，之后在这个点周围不断随机即可。

然后就做完了。要注意点的横纵坐标是实数，以及 UVA 的输出格式。

#include<ctime>
#include<cmath>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
#define maxn 100010
#define INF 0x3f3f3f3f
//#define int long long

using namespace std;

int t,n;
double nowx,nowy,Now,del;
double limx,limy,ansx,ansy,ans;
//ansx,ansy 是最优横纵坐标
//limx,limy 存可能最优的横纵坐标
const double lim=1e-10,d=0.996;
struct node{double x,y;}a[maxn];

int read(){
  int s=0,w=1;char ch=getchar();
  while(ch<'0'||ch>'9'){if(ch=='-')w=-1;ch=getchar();}
  while(ch>='0'&&ch<='9')s=(s<<1)+(s<<3)+ch-'0',ch=getchar();
  return s*w;
}

double Dis(double sx,double sy,double tx,double ty){
  return sqrt(pow(sx-tx,2)+pow(sy-ty,2));
}

double calc(double sx,double sy){
  double res=0.0;
  for(int i=1;i<=n;i++)
    res+=Dis(sx,sy,a[i].x,a[i].y);
  return res;
}

void SA(){
  double T=2021;
  while(T>lim){
    nowx=limx+((rand()<<1)-RAND_MAX)*T;
    nowy=limy+((rand()<<1)-RAND_MAX)*T;
    //在靠中点的点周围进行随机。
    //nowx,nowy 是当前随机到的点。
    Now=calc(nowx,nowy);del=Now-ans;
    if(del<0)ans=Now,ansx=limx=nowx,ansy=limy=nowy;
    else if(exp(-del/T)>(double)rand()/RAND_MAX)limx=nowx,limy=nowy;
    T*=d;
  }
}

int main(){
  t=read();
  for(int q=1;q<t;q++){
    n=read();
    for(int i=1;i<=n;i++)
      scanf("%lf%lf",&a[i].x,&a[i].y),
      limx+=a[i].x,limy+=a[i].y;
    limx/=(double)n;limy/=(double)n;
    ansx=limx;ansy=limy;ans=calc(limx,limy);
    SA();SA();SA();SA();SA();SA();SA();SA();
    printf("%.0lf

",ans);limx=limy=ansx=ansy=ans=0.0;
  }
  n=read();
  for(int i=1;i<=n;i++)
    scanf("%lf%lf",&a[i].x,&a[i].y),
    limx+=a[i].x,limy+=a[i].y;
  limx/=(double)n;limy/=(double)n;
  ansx=limx;ansy=limy;ans=calc(limx,limy);
  SA();SA();SA();SA();SA();SA();SA();SA();
  printf("%.0lf
",ans);
  return 0;
}

[]

[JSOI2016]炸弹攻击1

和上个题很像，不过要注意几个地方。

当我们随机到一个点时，我们要钦定它的半径。

由于边缘的敌人可以被消灭，而边缘的建筑不会受损，所以我们可以一直扩展到边缘与建筑边缘重合为止。

所以我们当前点的半径就是 到建筑中心的距离减去建筑半径 的最小值。

由于数据范围很小，所以我们随机完每个点后，暴力枚举每个建筑找这个值即可。、

注意这题求的是最大值，所以符号与上题不同。剩下的与上面基本相似。

#include<ctime>
#include<cmath>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
#define maxn 100010
#define INF 0x3f3f3f3f
//#define int long long

using namespace std;

int n,m;
double nowx,nowy,now,del;
const double lim=1e-10,d=0.997;
double R,limx,limy,ansx,ansy,ans;
struct enemy{double x,y;}e[maxn];
struct build{double x,y,r;}a[maxn];

int read() {
  int s=0,w=1;char ch=getchar();
  while(ch<'0'||ch>'9'){if(ch=='-')w=-1;ch=getchar();}
  while(ch>='0'&&ch<='9')s=(s<<1)+(s<<3)+ch-'0',ch=getchar();
  return s*w;
}

double Dis1(double sx,double sy,double dx,double dy) {
  return sqrt(pow(sx-dx,2)+pow(sy-dy,2));
}

double Dis2(double sx,double sy,build a) {
  return sqrt(pow(sx-a.x,2)+pow(sy-a.y,2))-a.r;
}

int calc(double sx,double sy) {
  int res=0.0;
  double dis=R;
  for(int i=1; i<=n; i++)
    dis=min(dis,Dis2(sx,sy,a[i]));
  for(int i=1; i<=m; i++)
    res+=(Dis1(sx,sy,e[i].x,e[i].y)<=dis);
  return res;
}

void SA() {
  double T=2020;
  while(T>lim) {
    nowx=limx+((rand()<<1)-RAND_MAX)*T;
    nowy=limy+((rand()<<1)-RAND_MAX)*T;
    now=calc(nowx,nowy);del=now-ans;
    if(del>0) ans=now,ansx=limx=nowx,ansy=limy=nowy;
    else if(exp(-del/T)<(double)rand()/RAND_MAX)limx=nowx,limy=nowy;
    T*=d;
  }
}

int main() {
  n=read();m=read();
  scanf("%lf",&R);
  for(int i=1; i<=n; i++)
    scanf("%lf%lf%lf",&a[i].x,&a[i].y,&a[i].r);
  for(int i=1; i<=m; i++)
    scanf("%lf%lf",&e[i].x,&e[i].y),
  limx+=e[i].x,limy+=e[i].y;
  limx/=(double)m;ansx=limx;
  limy/=(double)m;ansy=limy;ans=calc(limx,limy);
  while((double)clock()/CLOCKS_PER_SEC<0.40)SA();
  printf("%.0lf
",ans);
  return 0;
}

[]

[HAOI2006]均分数据

我们肯定是要使每一组的值都尽量相似，才能保证方差更小，标准差才能更小。

手动枚举每一组分确实麻烦，所以我们考虑换个简单方法。同时要结合它的随机性这一特点来处理。

“你们想想小学生怎么分啊！！！（” ——caq。

我们可以贪心地来分，我们每找到一个数，就将它加到当前值最少的组合里，然后计算贡献。

这样看起来很无厘头，但是我们是随机啊！

我们随机这个序列的顺序，然后用这个方法来做。通过多遍的随机化，我们有很大地概率可以分到正确的方案。

至于随机打乱顺序，我们可以用 random_shuffle 来实现。

注意这题的一个坑点，我们分完后只有 (m) 组，但题目中给出的式子被除数是 (n)。

#include<ctime>
#include<cmath>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
#define maxn 25
#define INF 3e5
//#define int long long

using namespace std;

int n,m;
double a[maxn],x[maxn];
double sum,ans=1000.0,now,del;

int read(){
  int s=0,w=1;char ch=getchar();
  while(ch<'0'||ch>'9'){if(ch=='-')w=-1;ch=getchar();}
  while(ch>='0'&&ch<='9')s=(s<<1)+(s<<3)+ch-'0',ch=getchar();
  return s*w;
}

double calc(){
  double res=0.0;
  for(int i=1;i<=m;i++)
    res+=pow(double(x[i])-sum,2);
  res/=(double)m;res=sqrt(res);
  return res;
}

void SA(){
  for(int i=1;i<=INF;i++){
    random_shuffle(a+1,a+n+1);
    memset(x,0,sizeof x);
    int Minx=INF;
    for(int i=1;i<=n;i++){
      int now=1;
      for(int j=1;j<=m;j++)
        if(x[j]<x[now]) now=j;
        x[now]+=a[i];
    }	
    now=calc();ans=min(ans,now);
  }
}

int main(){
  n=read();m=read();
  for(int i=1;i<=n;i++)
    scanf("%lf",&a[i]),sum+=(double)a[i];
  sum/=(double)m;SA();
  printf("%.2lf
",ans);
  return 0;
}

至于这题为什么不用模拟降温来实现……因为这题答案如果不优就直接会被遗弃，不需要计算概率，所以直接循环也可以啦。

[]

[SCOI2008]城堡

是一道有了思路几分钟出代码的黑题呢。

首先不管怎样，我们计算贡献需要用到两个点之间的最短距离，所以跑最短路。(nle 50) 可以 Floyd。

另外我们利用人类智慧发现，城堡越多越优，所以我们最多放 (k) 个那就放 (k) 个。

计算贡献的时候，我们暴力枚举每个点与其它的城堡点，算出他们的最远距离即可。

然后我们随机每个城堡放在哪里。

我们可以钦定前 (k) 个城市有城堡，算出贡献然后比较并撤回操作。之后再随机交换两个无城堡的点，同样钦定前 (k) 个城市有城堡。如此多次。

通过多次随机，我们大概率可以找出正确的城市。

注意处理重边。然后就完成了，没啥坑点，也不难写，不知道为啥评黑。

#include<ctime>
#include<cmath>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
#define maxn 1010
#define INF 0x3f3f3f3f
//#define int long long

using namespace std;

bool vis[maxn],flag[maxn];
int n,m,k,r[maxn],limx,limy,pos[maxn];
int ans,now,del;
int Dis[maxn][maxn];
const double d=0.998,lim=1e-10;

int read(){
  int s=0,w=1;char ch=getchar();
  while(ch<'0'||ch>'9'){if(ch=='-')w=-1;ch=getchar();}
  while(ch>='0'&&ch<='9')s=(s<<1)+(s<<3)+ch-'0',ch=getchar();
  return s*w;
}

void Floyd(){
  for(int k=1;k<=n;k++)
    for(int i=1;i<=n;i++)
      for(int j=1;j<=n;j++)
        Dis[i][j]=min(Dis[i][j],Dis[i][k]+Dis[k][j]);
}

int calc(){
  int res=-1;
  for(int i=1;i<=n;i++){
    int noww=INF;
    for(int j=1;j<=n;j++){
      if(!vis[j]) continue;
      noww=min(noww,Dis[i][j]);
    }
    res=max(res,noww);
  }
  return res;
}

void SA(){
  double T=2020;
  while(T>lim){
    T*=d;int sum=0;
    limx=rand()%n+1,limy=rand()%n+1;
    for(int i=1;i<=n;i++)
      if(!flag[i]&&vis[i])vis[i]=0;
    swap(pos[limx],pos[limy]);
    for(int i=1;i<=n;i++){
      if(vis[pos[i]]) continue;
      vis[pos[i]]=1;sum++;
      if(sum==k) break;
    }
    now=calc();del=now-ans;if(del<0) ans=now;
    else if(exp(-del/T)>(double)rand()/RAND_MAX){T*=d;continue;}
    else swap(pos[limx],pos[limy]);
    T*=d;
  }
}

int main(){
  memset(Dis,0x3f,sizeof Dis);
  n=read();m=read();k=read();
  for(int i=1;i<=n;i++)r[i]=read()+1,pos[i]=i,Dis[i][i]=0;
  for(int i=1;i<=n;i++)
    Dis[i][r[i]]=Dis[r[i]][i]=min(Dis[i][r[i]],read());
  for(int i=1,x;i<=m;i++)x=read()+1,vis[x]=flag[x]=1;
  Floyd();int sum=0;
  for(int i=1;i<=n;i++){
    if(vis[i]) continue;
    vis[i]=1;sum++;
    if(sum==k) break;
  }
  ans=calc();
  while((double)clock()/CLOCKS_PER_SEC<0.4) SA();
  printf("%d
",ans);
  return 0;
}

[]

Coloring

真的难题他来了。

这题我调参 + 调细节调了接近 2h，交了 16 遍

首先我们看看怎么计算一个矩阵的贡献。

可以发现，如果对每个格子都算四个方向的贡献，那么除了矩阵最外面的一圈之外，其余的都被算了两次。

所以我们可以对于每个格子，只算右方和下方的贡献。这样可以保证不重不漏。

然后我们第一次先随便填这个矩阵，之后对其进行随机微调。每次随机两个位置，将他们的颜色交换，同时统计贡献。

我们算出交换前两个位置四个方向的贡献，再算出交换后的，就可以比较出答案是否变优。

基本思路就是这样了。具体实现的时候，注意要开两个矩阵，一个存最优答案，一个存与当前随机比较的可能最优的答案。

另外这题调参很恶心，但是 caq 人很好，给了我他的参数，然后我在此基础上调了调就过了。

#include<ctime>
#include<cmath>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
#define maxn 1010
#define INF 0x3f3f3f3f
//#define int long long

using namespace std;

int ans,del;
int limx,limy,Limx,Limy;
int n,m,c,p[maxn],val[maxn][maxn];
int Ans[maxn][maxn];
int Minx=INF;
const double d=0.999996,lim=1e-15;

int read(){
  int s=0,w=1;char ch=getchar();
  while(ch<'0'||ch>'9'){if(ch=='-')w=-1;ch=getchar();}
  while(ch>='0'&&ch<='9')s=(s<<1)+(s<<3)+ch-'0',ch=getchar();
  return s*w;
}

int Changex(int x,int y){
  return (val[x-1][y]!=val[x][y]&&x!=1&&x!=n+1);
} 

int Changey(int x,int y){
  return (val[x][y-1]!=val[x][y]&&y!=1&&y!=m+1);
}

int Change(int x,int y){
  return ((val[x-1][y]!=val[x][y]&&x!=1&&x!=n+1)+(val[x][y-1]!=val[x][y]&&y!=1&&y!=m+1));
}

int change(int sx,int sy,int tx,int ty){
  int res=0;
  res-=Changex(sx+1,sy)+Changex(tx+1,ty)+Changey(sx,sy+1)+Changey(tx,ty+1)+Change(sx,sy)+Change(tx,ty);
  swap(val[limx][limy],val[Limx][Limy]);
  res+=Changex(sx+1,sy)+Changex(tx+1,ty)+Changey(sx,sy+1)+Changey(tx,ty+1)+Change(sx,sy)+Change(tx,ty);
  return res;
}

void SA(){
  double T=500;
  while(T>lim){
    limx=rand()%n+1;limy=rand()%m+1;
    Limx=rand()%n+1;Limy=rand()%m+1;
    del=change(limx,limy,Limx,Limy);
    if(ans+del<Minx){
      for(int i=1;i<=n;i++)
        for(int j=1;j<=m;j++)
          Ans[i][j]=val[i][j];
      ans=ans+del;Minx=ans;
    }
    else if(exp(-del/T)*RAND_MAX>(double)rand())ans+=del;
    else swap(val[limx][limy],val[Limx][Limy]);
    T*=d;
  }
}

int main(){
  srand(200030103);
  n=read();m=read();c=read();
  for(int i=1;i<=c;i++)p[i]=read();int now=1;
  for(int i=1;i<=n;i++)for(int j=1;j<=m;j++)
    {if(p[now])p[now]--;else now++,p[now]--;val[i][j]=Ans[i][j]=now;}
  for(int i=1;i<=n;i++)for(int j=1;j<=m;j++)ans+=Change(i,j);
  Minx=ans;while(clock()/CLOCKS_PER_SEC<1.5) SA();
  for(int i=1;i<=n;i++){for(int j=1;j<=m;j++)printf("%d ",Ans[i][j]);puts("");}
  return 0;
}

[]

写在最后

感觉这东西适合处理一些看上去不好解决，但数据范围很小的问题。

其实本质也是暴力，只不过是随机多次而已。最重要的不是退火的过程，而是计算函数的写法。

不过也不能漫无目的地随机，要善于运用人类智慧（雾），比如上面的钦定选择方案或者点的大致位置等等。

考场上骗分感觉还是可以的。

例题

P2210 Haywire
SP4587 FENCE3 - Electric Fences
P3878 [TJOI2010]分金币