膜模魔 模拟退火

 

一. 爬山算法 ( 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 ++ ;
}

注意事项:

  1. 温度T的初始值设置问题。 初始温度高,则搜索到全局最优解的可能性大,但因此要花费大量的计算时间;反之,则可节约计算时间,但全局搜索性能可能受到影响。
  2. 退火速度r问题。 模拟退火算法的全局搜索性能也与退火速度密切相关。同一温度下的“充分”搜索(退火)是相当必要的,但这需要计算时间。
  3. 温度管理问题 降温系数应为正的略小于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;  
}

 

本文系作者 @ 原创发布在 噓だ。未经许可,禁止转载。

喜欢()
评论 (0)
热门搜索
37 文章
3 评论
11 喜欢
Top