Simulate Anneal模擬退火演算法，是一種用於得到最優解的隨機化演算法。

如果可以打一手漂亮的隨機化搜尋，也許當你面對一籌莫展的神仙題時就有一把趁手的兵器了。

這篇題解將教你什麼？SA的基本思路，什麼時候能用SA。

標題是淺談，所以本篇部落格參雜了些許個人簡介，若有疑問或異議，歡迎提出指正。

我也很感謝你們給出的建議，它們真的能讓我變好、變強。

那麼我們進入本篇正題。

1. 什麼是模擬退火：

模擬退火是一種在廣大的搜尋空間尋找最優解的隨機化演算法。我們看名字就明白，這個演算法實在模擬物理中退火的過程。知識很多時候都是相通的，我們學習的大部分知識都是有用的。

為什麼通過模擬退火的過程可以得出最優解呢？在物理中，固體物質的退火過程和一般的組合優化問題有著很高的相似度。

2. 模擬退火基本要素：

滿足這些你就可以將SA演算法愉快地嵌入你的程式裡了。

第一個要素就是狀態空間和狀態的產生函式，我們需要有一個搜尋空間，而且範圍較大。

什麼是搜尋空間呢？我們學習搜尋的時候應該接觸過，我們的搜尋演算法本質就是在問題的求解空間中進行遍歷，尋找最優解。模擬退火的狀態空間也類似，它就是我們自定義出來的最優解的有限集合。

我們光有狀態空間還不夠啊，我們真正需要的是狀態。學函式的時候，老師就說過了，一個函式要有“域”。我們的狀態空間就是狀態產生函式的“域”。而我們如果想寫一個好的退火，我們的搜尋空間要足夠大，足夠讓函式生成不同的新解。

第二個就是候選解，我們在狀態空間內通過生成的隨機數在一定密度內隨機選擇我們的候選解。

其實還有一條是概率分佈，大部分採用均勻分佈，少數情況我們會用到指數分佈。

至於狀態轉移概率，我目前接觸到的SA演算法統統採用了Metropolis準則，我接下來會介紹它。

3. 模擬退火基本流程：

一. 由一個狀態產生函式從當前解產生一個新解，一般採用增量構造（即由原解加上/減去產生函式產生的值）來得到。

二. 計算新解的目標函式差Δt'。

三. 判斷新解是否被接受，

這裡介紹Metropolis準則：若Δt'<0，我們接受它，否則以exp(-Δt'/T)的多項式概率接受它。

**T是我們模擬退火過程中的溫度**

四. 當新解被接受時，用新解代替當前解，否則繼續下一輪試驗。

4. 模擬退火中的引數控制：

調參可以說是SA演算法最難的部分，也是決定你的演算法得分率的部分。

這裡分享一下神仙FlashHu（LCT導師Orz）寫SA時調參的經驗，我總結了一下就是這樣：

對於eps，可以根據資料範圍和精度要求粗略得到eps的大概大小，手動微調可以得到最終的eps。

對於T和ΔT（溫度的變化率，一般在0.95-0.99間），先開大一點跑出最優解，然後一邊退火一邊輸出當前的T、ans等資訊，大致感受解的下降速率，越均勻表示引數越好（神仙稱之為觀察法）

5. 注意事項：

SA可能會陷入當前解卡在區域性最小的情況，也就是這樣：

圖中，紅色是我們的當前解，藍色是我們的最優解，紅色的線代表SA演算法得到的新解，我們會發現，如果引數不佳，我們就有可能出現卡在區域性最優解的情況。本身演算法的設計就有避免這一種情況，還記得嗎？Metropolis準則，就算這不是最優解，我們也以一定概率接受它（答案不更新），就是為了防止這種情況，然而在引數不佳的情況下，它仍可能發生。所以我們要改進對溫度的控制方式。

這裡還是以那道經典到不能再經典的模擬退火模板做例題：平衡點/吊打xxx

這裡直接貼上程式碼，關於演算法中需要注意的地方我會打上註釋。

#include<bits/stdc++.h>
#define down 0.997//ΔT，模擬徐徐降溫 
using namespace std;
inline int read(){
    int data=0,w=1;char ch=0;
    while(ch!='-' && (ch<'0'||ch>'9'))ch=getchar();
    if(ch=='-')w=-1,ch=getchar();
    while(ch>='0' && ch<='9')data=data*10+ch-'0',ch=getchar();
    return data*w;
}
int n;
struct point{
    int x,y,w;
}object[2019];//物體資訊 
double ansx,ansy,answ;//答案 
double energy(double x,double y){//物理學知識：能量總和越小越穩定 
    double r=0,dx,dy;
    for(int i=1;i<=n;i++){
        dx=x-object[i].x;dy=y-object[i].y;
        r+=sqrt(dx*dx+dy*dy)*object[i].w;//力臂乘重力 
    }
    return r;
}
void SA(){
    double t=3e3+10;//初始溫度要高 
    while(t>1e-15){//exp略大於0 
        double ex=ansx+(rand()*2-RAND_MAX)*t;
        double ey=ansy+(rand()*2-RAND_MAX)*t;
        double ew=energy(ex,ey);
        double de=ew-answ;
        if(de<0){//此答案更優 
            ansx=ex;ansy=ey;answ=ew;
        }else if(exp(-de/t)*RAND_MAX>rand()){//Metropolis準則，以多項式概率接受 
            ansx=ex;ansy=ey;
        }t*=down;//逐步降溫 
    }
}
void solve(){
    SA();SA();SA();SA();SA();
}
int main(){
    n=read();
    for(int i=1;i<=n;i++){
        object[i].x=read();object[i].y=read();object[i].w=read();
        ansx+=object[i].x;ansy+=object[i].y;
    }
    ansx/=n;ansy/=n;//設平均值為初值 
    answ=energy(ansx,ansy);
    solve();
    printf("%.3lf %.3lf\n",ansx,ansy);
    return 0;
}

遊戲結