模拟退火学习笔记

%你退火，听上去一个十分高级的算法，实际上他只是贪心的一种随机化，他的思想很简单，一般是随机几个数，然后用这几个数来算出答案，并更新答案。

但是，如果只是单纯的随机几个数，那么他的也不会叫做%你退火，所以这个算法也有他自己的奥妙。

%你退火实际上是爬山算法的一种升级版，爬山算法就是要不断的选择更优解，如果到了一个高峰（低谷），他就会满足于这个局部最优解，但是外面的世界他并不会去探索。

而模拟退火，就很好地解决了这一个问题，他好像一个喝醉酒的兔子，在群山之间乱窜，走着走着，就慢慢醒来了，窜动的频率也就越来越小，这就是模拟退火。

模拟退火实际上是仿照金属退火，但是为了方便理解，我们还是用以上兔子醉酒的例子来理解：

一只兔子喝了很多酒，他的体内有很多酒精也不知道喝了酒为什么会醉，所以就用酒精的量来表示酒醉的程度好了）

然后这只兔子就开始乱跳，找到了一个下坡就向下滚，陷入了低古，也会借着酒劲往上爬。

但随着酒精量的减少，兔子会渐渐地回复理智，不想爬山，于是他就只是借着下坡的力量或者渐渐减少的酒劲前进，知道体内没有酒精。

好了，扯淡结束，但是这和我们\(OI\)有什么关系吗？

其实模拟退火算法就是以上兔子醉酒的实现，我们先设定一个初始酒精量（我一般设为\(1926\)），然后再设定一个每一次爬山兔子损失的酒精量（由于酒精消耗越来越慢，所以我们用乘法来表示消耗，一般的乘数略小于\(1\)就行），然后兔子的酒精量越来越少，知道一个值我们就认为他体内没有酒精了（我一般设为\(1e-18\)）

然后我们就可以用来操作我们的程序

例题：
一道最经典的模拟退火的题目就是P1337 [JSOI2004]平衡点 / 吊打TBR，代码如下：

#include<bits/stdc++.h>
using namespace std;
#define il inline
#define re register
il int read()
{
    re int x=0,k=1;re char c=getchar();
    while(c<'0'||c>'9'){if(c=='-') k=-1;c=getchar();}
    while(c>='0'&&c<='9')x=(x<<3)+(x<<1)+(c^48),c=getchar();
    return x*k;
}
#define D double
#define T 1926//初始酒精量
#define EPS 1e-18//可以理解为酒精级限量
#define delta 0.9986//每一次的酒精减少量
#define maxn 1005
#define inf 123456789101112.0
struct node
{
    int x,y,w;
}e[maxn];
int n;
il D getans(D x,D y)
{
    D ans=0;
    for(re int i=1;i<=n;++i)
    {
        ans+=(sqrt((e[i].x-x)*(e[i].x-x)+(e[i].y-y)*(e[i].y-y)))*e[i].w;
    }
    return ans;
}//求出当一个点的坐标是(x,y)的答案，一般按照题意模拟即可
D xans,yans,ans=inf;
il void get()
{
    re D xx=xans,yy=yans,t=T;
    while(t>EPS)//只要现在的酒精量还大于极限酒精量
    {
        re D xt=xx+(rand()*2-RAND_MAX)*t;
        re D yt=yy+(rand()*2-RAND_MAX)*t;
        //随机两个数，但是为什么不是直接rand呢？
        //我们想想兔子是怎么跑的，因为他喝醉了酒，所以它可能向前也可能向后，体现出来就是(rand()*2-RAND_MAX)
        //RAND_MAX是一个常量，表示rand的最大值，rand()*2-RAND_MAX的意思是：rand*2可能比RAND_MAX小也可能大，而且是等概率的，所以我们用这个式子可以很好的表示向前和向后
        //为什么要*t呢？因为兔子体内酒精越少，越不想动，所以它动的频率也就越低，由于t是单调递减，所以我们*t可以模拟出酒精越少越不动的样子
        re D now=getans(xt,yt);
        re D deta=now-ans;
        if(deta<0)//如果新的两个点比答案大，就更新答案，也就是相当于往下滚
        {
            xx=xans=xt;
            yy=yans=yt;
            ans=now;
        }
        else if(RAND_MAX*exp(-deta/t)>rand())//这个柿子据说是科学家推出来的，意思是如果当前的答案不优，兔子也会凭着酒劲乱跑
        {
            xx=xt,yy=yt;
        }
        t*=delta;
    }
}
il void SA()
{
	//由于是随机，所以我们要尽量的多做几遍来保证正确率
    re int pax=11;
    while(pax--) get();
}
signed main()
{
#ifndef ONLINE_JUDGE
    freopen("a.in","r",stdin);
#endif
    srand(19260817);
    n=read();
    for(re int i=1;i<=n;++i)
    {
        e[i].x=read(),e[i].y=read(),e[i].w=read();
    }
    SA();
    printf("%.3lf %.3lf",xans,yans);
    return 0;
}

如果你有疑问：如果这份代码总是WA怎么办？

法1：调大T（用处不大）

法2：调大delta（一般就是调这个）

法3：调小eps（这个也很有用）

法4：洗把脸

例题2：P1429 平面最近点对（加强版）

这道题目和上一题没有什么本质区别，唯一注意的一点就是我们在随机两个点的时候，我们可以先把所有点按照x排好序，因为最近点对的x肯定不会太远，所以我们再算点的时候可以加一个限制，也就是不让两个点的x差太远。

代码如下：

#include<bits/stdc++.h>
using namespace std;
#define il inline
#define re register
il int read()
{
    re int x=0,k=1;re char c=getchar();
    while(c<'0'||c>'9'){if(c=='-') k=-1;c=getchar();}
    while(c>='0'&&c<='9')x=(x<<3)+(x<<1)+(c^48),c=getchar();
    return x*k;
}
#define inf 1234567890.0
#define debug printf("Now is Line %d\n",__LINE__)
#define D double
#define T 1926
#define delta 0.99992
#define eps 1e-18
#define bd 11//这是一个常量，表示允许两个点的x排名的差值
struct node
{
    D x,y;
}e[200005];
il D far(int a,int b)
{
    re D x=e[a].x-e[b].x,y=e[a].y-e[b].y;
    return sqrt(x*x+y*y);
}//计算两点的距离
int n,now;
D ans=inf;
il int p(int x)
{
    return (rand()&1)?x:-x;
}
il void get()
{
    re D t=T;
    now=rand()%n+1;
    while(t>eps)
    {
        re int tnow=now+p(rand()%bd);
        if(tnow==now||tnow<1||tnow>n) continue;//越界了就continue
        re D nans=far(tnow,now);
        re D data=nans-ans;
        if(data<0)
        {
            now=tnow;
            ans=nans;
        }
        else if(RAND_MAX*exp(data/t)>rand())
        {
            now=tnow;
        }
        t*=delta;
    }
}
il void SA()
{
    re int pax=12;
    while(pax--) get();
}
il bool cmp(node a,node b)
{
    return a.x<b.x;
}
int main()
{
#ifndef ONLINE_JUDGE
    freopen("a.in","r",stdin);
#endif
    srand(19260817);
    n=read();
    for(re int i=1;i<=n;++i)
    {
        e[i].x=read(),e[i].y=read();
    }
    sort(e+1,e+n+1,cmp);
    SA();
    printf("%.4lf",ans);
    return 0;
}

例题3：P3959 宝藏

这道题的模拟退火实际上不是很标准，我们只是随机生成一个顺序，用Prim的贪心思想跑一边答案就行。

代码如下：

#include<bits/stdc++.h>
using namespace std;
#define il inline
#define re register
#define debug printf("Now is Line %d\n",__LINE__)
il int read()
{
    re int x=0,k=1;re char c=getchar();
    while(c<'0'||c>'9'){if(c=='-') k=-1;c=getchar();}
    while(c>='0'&&c<='9')x=(x<<3)+(x<<1)+(c^48),c=getchar();
    return x*k;
}
#define T 1926
#define delta 0.993
#define eps 1e-18
#define inf 1234567890
#define maxn 20
int n,m,a[maxn][maxn],dis[maxn],ans=inf,id[maxn];
il int doit()
{
    memset(dis,-1,sizeof(dis));
    dis[id[1]]=0;
    re int spend=0;
    for(re int i=2;i<=n;++i)
    {
        re int minn=inf;
        for(re int j=1;j<i;++j)
        {
            if(a[id[j]][id[i]]<inf&&(dis[id[j]]+1)*a[id[j]][id[i]]<minn)
            {
                dis[id[i]]=dis[id[j]]+1;
                minn=(dis[id[j]]+1)*a[id[j]][id[i]];
            }//模拟Prim的思想贪心
        }
        if(minn==inf) return minn;//如果这个顺序不合法（也就是前面的点无法到达后面的点）
        spend+=minn;
    }
    return spend;
}
il void SA()
{
    re int t=T;
    while(t>eps)
    {
        re int x=rand()%n+1,y=rand()%n+1;
        swap(id[x],id[y]);
        re int now=doit();
        if(now<ans) ans=now;//这个其实没有完全按照兔子喝醉的思想，为什么呢？我们后面会讲。
        t*=delta;
    }
}
signed main()
{
#ifndef ONLINE_JUDGE
    freopen("a.in","r",stdin);
#endif
    n=read(),m=read();
    memset(a,127,sizeof(a));
    for(re int i=1;i<=n;++i) id[i]=i,a[i][i]=0;
    for(re int i=1;i<=m;++i)
    {
        re int u=read(),v=read(),w=read();
        a[u][v]=a[v][u]=min(a[u][v],w);
    }
    re int pax=233;
    while(pax--) SA();
    printf("%d",ans);
    return 0;
}

为什么说他是不完整的模拟退火呢？因为以下代码也可以过

#include<bits/stdc++.h>
using namespace std;
#define il inline
#define re register
#define debug printf("Now is Line %d\n",__LINE__)
il int read()
{
    re int x=0,k=1;re char c=getchar();
    while(c<'0'||c>'9'){if(c=='-') k=-1;c=getchar();}
    while(c>='0'&&c<='9')x=(x<<3)+(x<<1)+(c^48),c=getchar();
    return x*k;
}
#define T 1926
#define delta 0.993
#define eps 1e-18
#define inf 1234567890
#define maxn 20
int n,m,a[maxn][maxn],dis[maxn],ans=inf,id[maxn];
il int doit()
{
    memset(dis,-1,sizeof(dis));
    dis[id[1]]=0;
    re int spend=0;
    for(re int i=2;i<=n;++i)
    {
        re int minn=inf;
        for(re int j=1;j<i;++j)
        {
            if(a[id[j]][id[i]]<inf&&(dis[id[j]]+1)*a[id[j]][id[i]]<minn)
            {
                dis[id[i]]=dis[id[j]]+1;
                minn=(dis[id[j]]+1)*a[id[j]][id[i]];
            }
        }
        if(minn==inf) return minn;
        spend+=minn;
    }
    return spend;
}
il void SA()
{
    random_shuffle(id+1,id+n+1);
    ans=min(ans,doit());//每一次随机生成一个序列，模拟答案即可
}
signed main()
{
#ifndef ONLINE_JUDGE
    freopen("a.in","r",stdin);
#endif
    n=read(),m=read();
    memset(a,127,sizeof(a));
    for(re int i=1;i<=n;++i) id[i]=i,a[i][i]=0;
    for(re int i=1;i<=m;++i)
    {
        re int u=read(),v=read(),w=read();
        a[u][v]=a[v][u]=min(a[u][v],w);
    }
    re int pax=23333;
    while(pax--) SA();
    printf("%d",ans);
    return 0;
}

其实以上两个代码的本质是一样的，所以其实这个算法只能叫做随机化乱搞

例题4：P3382 【模板】三分法

其实这个模板也是可以用模拟退火来做的，只是lr比较小，所以我们把初始酒精量调低一点就行了

代码如下：

#include<bits/stdc++.h>
using namespace std;
#define file(a) freopen(#a".in","r",stdin),freopen(#a".out","w",stdout)
#define re register
#define il inline
#define debug printf("Now is line %d\n",__LINE__)
#define D double
#define maxn 20
#define inf 123456789.0
#define T 1
#define delta 0.997
#define eps 1e-18
#define rand_max 1
int n,da;
D l,r,ans=-inf,xx,a[maxn],zb;
il D f(D u)
{
    re D ans=0;
    for(re int i=0;i<=n;++i)
    {
        ans=ans*u+a[i];
    }
    return ans;
}//统计答案的时候用了秦九韶公式：
//ax^n+bx^n-1+……yx+z=z+x(y+x(……+(b+ax)))
//其实不理解用快速幂就行了
il int random(int x)
{
    return rand()%x;
}
il void SA()
{
    re D t=T;
    xx=(l+r)/2.0;
    while(t>eps)
    {
        re D now=xx+(rand()*2-RAND_MAX)*t;
        if(now<l) now=l;
        if(now>r) now=r;
        re D nans=f(now);
        re D data=nans-ans;
        if(data>0)
        {
            xx=zb=now;
            ans=nans;
        }
        else if(RAND_MAX*exp(data/t)>rand())
        {
            xx=now;
        }
        t*=delta;
    }
}
signed main()
{
#ifndef ONLINE_JUDGE
    freopen("a.in","r",stdin);
#endif
    srand(19260817);
    scanf("%d%lf%lf",&n,&l,&r);
    for(re int i=0;i<=n;++i)
    {
        scanf("%lf",&a[i]);
    }
    re int pax=20;
    while(pax--)
    {
        SA();
    }
    printf("%.5lf\n",zb);
    return 0;
}