G
N
I
D
A
O
L

模拟退火算法详解

别着急,干货在最后面!!!

(本文用c++实现,可以在评论区讨论,后面还有情况的话还会更新,有问题欢迎指正哦~)

可以在右上角看目录,左下角点歌哦(不行的话刷新一下就好了~)

本文章也介绍了模拟退火的使用情景,以免误入歧途(本蒟蒻就是)。

很多人都学过贪心,但是贪心在一些情况并不适用,比如:

已知我们从黄色出发,找最小值。

贪心策略当然是一直往函数大小减小的地方偏移——但是,万一不是单谷呢?我们会陷入如图的蓝色中无法自拔。

肯能你会想到:随机找一个点出发,然后贪心找最小值?多随机几遍,然后求全局最小值?

你会发现复杂度暴增!!!!!!!!——所以如何处理这种问题呢?

模拟退火:啊,对对对~

是的!模拟退火就是一种类似于随机化贪心的一个算法,在OI界也小有名气(冥器)!(如题[NOIP2021] 方差

原理图:

如图:在物理应用中分子排布可能是紊乱的,如果我们将它升温然后缓慢降温,就可以生成完美的晶形!

而对于我们求解的:

 

怎么形象描述它呢?

一个有自己一定卡路里的人爬山,想翻山找远方的草药给自己心爱的妻子,但是他并不知道山的那头是什么。所以他会在卡路里多的时候尽量去远方探险,但是每次会花费卡路里以至于他后面不能翻过太高的山丘,而且他的背包蛮大的,装填着无数爱的草药而芳香四溢。

所以我们立刻(啊,对对对~)能设定模拟退火的参数:

1.初始温度 T (1000-7000)

2.末尾温度 P(1e-6~1e-15)

3.降温系数 K (0.91~0.9975)

4.状态空间(被降温物体) S

5.当前能量 E(new)

6.全局能量 E(old)

一:Metropolis准则

以概率接受新状态:

 

 

 

这就是物理(化学)方面类似的推论——一定概率的更新。

什么意思呢?

我们已知:当前能量 E ( new ) , 全局能量 E ( old ),那么我们的目标是什么,不就是减少目前的能量吗?

所以:当当前能量少于全局能量(即更新前的能量),那么我们有概率为 1 的更新概率;

      当当前能量大于全局能量(即更新前的能量),那么我们有概率为 exp( - (E( new )-E( old ))/T) 的更新概率 ( T为当前温度)

[exp(x)函数:e的x次方的函数  如 exp(1)表示e的1次方=e=2.718281828…  exp(0)表示e的0次方=1  exp(2)表示e的平方=7.3890561…  e是一个常数,等于2.718281828…];

注意:有时候也不一定以以上方式更新,这只是比较妥的做法,概率方面是可以自己定的,但是一定以当前能量与全局能量的关系来设定的。(除非直接暴力的随机算法)

二:生成新温度

那么怎么生成新的当前温度呢?,以生成小数为例:

当前将更新温度=全局温度+(rand()*2-RAND_MAX)*t;
if(不在状态空间内){
    当前将更新温度=fmod(当前将更新温度,状态空间大小)
}

  即:在当前状态的邻域结构内以一定概率方式(均匀分布、正态分布、指数分布等)产生。

三:温度更新函数

若固定每一温度,算法均计算至平稳分布,然后下降温度,则称为时齐算法;

若无需各温度下算法均达到平稳分布,但温度需按一定速率下降,则称为非时齐算法。

本人用的:

T*=K;

四:外循环终止准则 

本人使用的:

(t>1e-15)//可以改大一点

其他常用方法:

    (1)设置终止温度的阈值。

    (2)设置外循环迭代次数。

    (3)算法搜索到的最优值连续若干步保持不变。

    (4)概率分析方法。

五:实现流程图:

六:关于其他类似算法的优缺点比较

遗传算法:其优点是能很好地处理约束,跳出局部最优,最终得到全局最优解。缺点是收敛速度慢,局部搜索能力弱,运行时间长,容易受到参数的影响。

模拟退火:具有局部搜索能力强、运行时间短的优点。缺点是全局搜索能力差,容易受到参数的影响。

爬山算法:显然爬山算法简单、效率高,但在处理多约束大规模问题时,往往不能得到较好的解决方案。

七:退火口诀:

初始温度小心设(1000-3000),又粗又大wa一脸

多次sa更保险,忘了卡时直接T[if((double)clock()/CLOCKS_PER_SEC>=0.993)](七遍模拟退火也行)

退火系数大胆设,不过0.9975会很厄

全局、状态不一样,全局必须菊部优

百年骗分一场空,不开srand见祖宗

退火需谨慎,退火不规范,灵封两行泪

八:使用条件:

我们可以看下这道题:https://www.luogu.com.cn/problem/P6879

第一次做看到最优解我直接退火了,要不是捆绑数据还以为自己正确了。

代码:

#include <bits/stdc++.h>
using namespace std;
const int N=500;
const double dw=0.9975;
int resx,n,ans,l,T[N];
long long a[N];
bool st[N];
int em(int x){
	memset(st,false,sizeof st);
	int res=0,pd=0;
	pd=min(abs(x-0),abs(x-l));
	for(int i=1;i<=n*2;i++){
		if(abs(a[i]-x)+pd<=T[i]&&!st[i]){
			res++;st[i%(2*n)]=st[(i+n)%(2*n)]=true;
		}
	}
	return res+(rand()-RAND_MAX+1)%2;
}
void sa(){
	double t=1500;
	while(t>1e-15){
		int x=abs(resx+rand()*2-RAND_MAX+1)%l;
		int e=em(x),dt=ans-e;
		if(dt<0){
			ans=e;resx=x;
		}else if(exp((double)(-dt/t))*RAND_MAX<rand()){
			resx=x;
		}
		t*=dw;
	}
}
int main(){
	srand((unsigned)time(0));
	cin>>n>>l;
	for(int i=1;i<=n;i++){
		scanf("%lld",&a[i]);
		a[i+n]=a[i]+l;
	}
	for(int i=1;i<=n;i++){
		scanf("%d",&T[i]);
		T[i+n]=T[i];
	}
	l*=2;
	ans=em(0);
	sa();sa();sa();sa();sa();sa();sa();
	cout<<ans; 
}

  

发现可以A掉一部分?其实是不对滴——计算当前物品能量的函数写的并不是正确的。

因为当你到了圆上的一点后,发现你仍然需要下一步决策走向最优解,而并不是直接可以计算出来可以得到的贡献值。

这不禁让我反思——退火可以解DP题吗?

我们可以观察一下这道题:大的最优解是由更小的最优解转移过来的,就像一个树形结构,由子节点向父节点转移,像这样的题是不可以用退火的。

但是比如01背包, [NOIP2021]方差 和等类贪心题目是可以做的,因为你可以通过概率水掉局部最优解对全局最优解的不可转移。

所以,每次做题要看退火的当前能量计算复不复杂,有没有子问题的限制!

然后是[NOIP2021]方差的实现(玄学万岁):

 1 #include <bits/stdc++.h>
 2 using namespace std;
 3 #define LL long long
 4 const int N=1e5+10;
 5 const double dw=0.9975;
 6 int a[N],n,c[N];
 7 long long ans;
 8 bool cmp(int a,int b){
 9     return a>b;
10 }
11 LL en(){
12     LL em=0;
13     LL ranss=0;
14     for(int i=2;i<=n;i++){
15         a[i]=a[i-1]+c[i];
16     }
17     for(int i=1;i<=n;i++){
18         em+=(long long)a[i]*a[i];
19     }em=(long long)em*n;
20     for(int i=1;i<=n;i++){
21         ranss=(long long)ranss+a[i];
22     }ranss=(long long)ranss*ranss;
23     return (long long)(em-ranss);
24 }
25 void sa(){
26     double t=1000;
27     while(t>1e-15){
28         if((double)clock()/CLOCKS_PER_SEC>=0.993){
29             cout<<ans;
30             exit(0);
31         }
32         int x=rand()%(n-1)+2,y=rand()%(n-1)+2;
33         while(x==y)x=rand()%(n-1)+2;
34         swap(c[x],c[y]);
35         LL m=en(),dt=ans-m;
36         if(dt>0){
37             ans=m;
38         }else if((double)rand()>=(double)RAND_MAX*(double)exp((double)dt/t)){
39             swap(c[x],c[y]);
40         }
41         t*=dw;
42     }
43 }
44 int main(){
45     srand((unsigned)time(0));
46     cin>>n;
47     for(int i=1;i<=n;i++){
48         scanf("%d",&a[i]);
49         c[i]=a[i]-a[i-1];
50     }
51     sort(c+2,c+n/2+1,cmp);
52     sort(c+n/2+1,c+n+1);
53     ans=en();
54     while(1)sa();
55 }

 附赠导论:骗分

posted @ 2022-01-25 22:50  中指半仙  阅读(3842)  评论(4编辑  收藏  举报