膜模魔 模拟退火
一. 爬山算法 ( Hill Climbing )
介绍模拟退火前,先介绍爬山算法。爬山算法是一种简单的贪心搜索算法,该算法每次从当前解的临近解空间中选择一个最优解作为当前解,直到达到一个局部最优解。
爬山算法实现很简单,其主要缺点是会陷入局部最优解,而不一定能搜索到全局最优解。如图1所示:假设C点为当前解,爬山算法搜索到A点这个局部最优解就会停止搜索,因为在A点无论向那个方向小幅度移动都不能得到更优的解。
二. 模拟退火(SA,Simulated Annealing)思想
爬山法是完完全全的贪心法,每次都鼠目寸光的选择一个当前最优解,因此只能搜索到局部的最优值。模拟退火其实也是一种贪心算法,但是它的搜索过程引入了随机因素。模拟退火算法以一定的概率来接受一个比当前解要差的解,因此有可能会跳出这个局部的最优解,达到全局的最优解。以图1为例,模拟退火算法在搜索到局部最优解A后,会以一定的概率接受到E的移动。也许经过几次这样的不是局部最优的移动后会到达D点,于是就跳出了局部最大值A。
模拟退火算法描述:
若J( Y(i+1) )>= J( Y(i) ) (即移动后得到更优解),则总是接受该移动
若J( Y(i+1) )< J( Y(i) ) (即移动后的解比当前解要差),则以一定的概率接受移动,而且这个概率随着时间推移逐渐降低(逐渐降低才能趋向稳定)
这里的“一定的概率”的计算参考了金属冶炼的退火过程,这也是模拟退火算法名称的由来。
根据热力学的原理,在温度为T时,出现能量差为dE的降温的概率为P(dE),表示为:
P(dE) = exp( dE/(kT) )
其中k是一个常数,exp表示自然指数,且dE<0。这条公式说白了就是:温度越高,出现一次能量差为dE的降温的概率就越大;温度越低,则出现降温的概率就越小。又由于dE总是小于0(否则就不叫退火了),因此dE/kT < 0 ,所以P(dE)的函数取值范围是(0,1) 。
随着温度T的降低,P(dE)会逐渐降低。
我们将一次向较差解的移动看做一次温度跳变过程,我们以概率P(dE)来接受这样的移动。
下面给出模拟退火的伪代码表示。
/*
* J(y):在状态y时的评价函数值
* Y(i):表示当前状态
* Y(i+1):表示新的状态
* r: 用于控制降温的快慢
* T: 系统的温度,系统初始应该要处于一个高温的状态
* T_min :温度的下限,若温度T达到T_min,则停止搜索
*/
while( T > T_min )
{
dE = J( Y(i+1) ) - J( Y(i) ) ;
if ( dE >=0 ) //表达移动后得到更优解,则总是接受移动
Y(i+1) = Y(i) ; //接受从Y(i)到Y(i+1)的移动
else{ // 函数exp( dE/T )的取值范围是(0,1) ,dE/T越大,则exp( dE/T )也
if ( exp( dE/T ) > random( 0 , 1 ) )
Y(i+1) = Y(i) ; //接受从Y(i)到Y(i+1)的移动
}
T = r * T ; //降温退火 ,0<r<1 。r越大,降温越慢;r越小,降温越快
/*
* 若r过大,则搜索到全局最优解的可能会较高,但搜索的过程也就较长。若r过小,则搜索的过程会很快,但最终可能会达到一个局部最优值
*/
i ++ ;
}
注意事项:
- 温度T的初始值设置问题。 初始温度高,则搜索到全局最优解的可能性大,但因此要花费大量的计算时间;反之,则可节约计算时间,但全局搜索性能可能受到影响。
- 退火速度r问题。 模拟退火算法的全局搜索性能也与退火速度密切相关。同一温度下的“充分”搜索(退火)是相当必要的,但这需要计算时间。
- 温度管理问题 降温系数应为正的略小于1.00的常数。
有几题裸题:
BZOJ 3680
这道题,一切自然变化进行的方向都是使能量降低,因为能量较低的状态比较稳定。
因为物重一定,绳子越短,重物越低,势能越小,势能与物重成正比,所以,只要使得也就是总的重力势能最小,就可以使系统平衡
#include<cstdio>
#include<cctype>
#include<cmath>
#include<cstdlib>
struct Node
{
int x,y,weight;
}node[10005];
int n;
double potential_energy(double nowx,double nowy)
{
double sum=0;
for(int i=1;i<=n;i++)
{
double delx=nowx-node[i].x;
double dely=nowy-node[i].y;
sum+=(sqrt(delx*delx+dely*dely))*node[i].weight;
//物重一定,绳子越短,重物越低,势能越小
//势能又与物重成正比
}
return sum;//在(nowx,nowy)处的总势能
}
double xans,yans;//最终答案
double ans=1e18+7,t;//势能与温度
const double delta=0.993;//降温系数
void simulate_anneal()
{
double xx=xans;//钦定一个初始位置
double yy=yans;
t=1926;//t是温度
while(t>1e-14)
{
double xtemp=xans+(rand()*2-RAND_MAX)*t;
double ytemp=yans+(rand()*2-RAND_MAX)*t;
//随机一个新的坐标,变化幅度为t
//这里要注意rand()和rand()*2-RAND_MAX的区别
//rand()的范围是0~RAND_MAX-1
//rand()*2-RAND_MAX的范围是-RAND_MAX到RAND_MAX-1
double new_ans=potential_energy(xtemp,ytemp);//计算当前解的势能
double DE=new_ans-ans;
if(DE<0)//如果是一个更优解
{
xx=xtemp;
yy=ytemp;//就接受
xans=xx;
yans=yy;
ans=new_ans;
}
else if(exp(-DE/t)*RAND_MAX>rand())//能否接受这个差
{
//更新坐标
xx=xtemp;
yy=ytemp;
}
t*=delta;//降温
}
}
void SA()//洗把脸就AC了
{
simulate_anneal();
simulate_anneal();
simulate_anneal();
}
int main()
{
scanf("%d",&n);
for(int i=1;i<=n;i++)
{
scanf("%lf%lf%lf",&node[i].x,&node[i].y,&node[i].weight);
node[i].x=Read();
node[i].y=Read();
node[i].weight=Read();
}
SA();
printf("%.3lf %.3lf",xans,yans);
return 0;
}
POJ2680
这题算是简化版的模拟退火,并没有用到:以概率P(dE)来接受这样的移动。
所以一开始用普通的模拟退火来做卡了半天,最后发现其实是简单的移动。
我们首先任选一个点作为球心,并找到点集中与它距离最远的点,我们让球心靠近最远的点,不断重复此过程,就可以让球心达到稳定态了!此时我们就找到了最小球。
代码:
#include <iostream>
#include<cstdio>
#include<algorithm>
#include<cstring>
#include<cmath>
using namespace std;
const double eps=1e-7;
struct point3D
{
double x,y,z;
} data[35];
int n;
double dis(point3D a,point3D b)
{
return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y)+(a.z-b.z)*(a.z-b.z));
}
double solve()
{
double step=100,ans=1e30,mt;
point3D z;
z.x=z.y=z.z=0;
int s=0;
while(step>eps)
{
for(int i=0; i<n; i++)
if(dis(z,data[s])<dis(z,data[i])) s=i;
mt=dis(z,data[s]);
ans=min(ans,mt);
z.x+=(data[s].x-z.x)/mt*step;
z.y+=(data[s].y-z.y)/mt*step;
z.z+=(data[s].z-z.z)/mt*step;
step*=0.98;
}
return ans;
}
int main()
{ // freopen("t.txt","r",stdin);
double ans;
while(~scanf("%d",&n),n)
{
for(int i=0; i<n; i++)
scanf("%lf%lf%lf",&data[i].x,&data[i].y,&data[i].z);
ans=solve();
printf("%.5f\n",ans);
}
return 0;
}
本文系作者 @rinbn 原创发布在 噓だ。未经许可,禁止转载。