目錄
寫在前面
感謝 caq 的傾情講解 
模擬退火是個隨機化算法，正確性有一定保證，但如果你想我一樣臉黑的話......
實測模擬退火做多了 rp 會掉
正文
簡介
模擬退火是一種隨機化算法，當一個問題的方案數量極大（甚至是無窮的）而且不是一個單峰函數時，我們常使用模擬退火求解。- Oi-wiki
什么是退火？
退火是一種金屬熱處理工藝，指的是將金屬緩慢加熱到一定溫度，保持足夠時間，然后以適宜速度冷卻。目的是降低硬度，改善切削加工性；消除殘余應力，穩定尺寸，減少變形與裂紋傾向；細化晶粒，調整組織，消除組織缺陷。准確的說，退火是一種對材料的熱處理工藝，包括金屬材料、非金屬材料。而且新材料的退火目的也與傳統金屬退火存在異同。---百度百科
扯遠了。
這個算法就是在溫度不斷降低的過程中，不斷地從當前位置尋找別的位置進行計算，溫度越低，也就是它的動能越小時，位置就會變化的越小，最后逐漸停留在最優解（或者附近）
算法流程
每次隨機一個新的狀態，如果狀態更優就更新答案，否則以一定概率接受這個狀態。
Metropolis准則
以求最小值為例。
-
如果 \(\Delta E < 0\)，說明當前解更優，直接更新即可
-
否則，如果
\[e^{\frac{-\Delta E}{T}} > \frac{\text{rand()}}{\text{RAND_MAX}} \]
就接受這個狀態。
- 否則 直接跳過。
為什么？
第一步因為是最優解所以一定選擇更新答案
第二步后邊的是一個隨機值我們暫且不論。
考慮整個退火過程，
假設溫度 \(T\) 不變，新的解越劣，\(\Delta E\) 越大，左項的值越小，接受的概率也越小。
假設 \(\Delta E\) 不變，隨着溫度的下降，求解的范圍也趨於穩定，\(T\) 越小，左項得值也越小，接受的概率也越小
扔一張圖可能會更好理解：
聽上去很扯 ，但它還是有一定的正確性的。
SA 函數
通常降溫系數 \(d\) 是一個很接近 \(1\) 的數，終止溫度 \(T_0\) 是一個很接近 \(0\) 的數
這里給一個偽代碼：
const double lim = ... // 溫度最小值，通常為 1e-10 左右
const double d = ... // 變化系數，通常為 0.996 左右
void SA() {
double T = ... // 初始溫度，通常為 2021 左右
while(T > lim) {
... // 獲取一個隨機的位置
now = calc(); // 計算當前位置的答案
del = now - ans; // 計算 變化量
if(del < 0) { // 以最小值為例
ans = now; // 更新答案
... // 更新答案和中間量的狀態
} else if(exp(-del/T) > (double)rand()/RAND_MAX) {
... // 一定概率選擇當前當前狀態
}
T *= d; // 降溫
}
}
計算函數 calc
依據題目而定，這里不給出
一些技巧/思想
如果想要隨機一個無限大平面內的一個點，可以這樣：
double nowx = limx + ((rand() << 1) - RAND_MAX) * T;
double nowy = limy + ((rand() << 1) - RAND_MAX) * T;
其中 nowx,nowy 是我們隨機的位置， limx, limy 是我們一個中間狀態的位置（注意不是答案的位置），
后面的那一坨剛好對應着溫度越小變化越小的實際情況。
我們有時為了使得到的解更有質量，會在模擬退火結束后，以當前溫度在得到的解附近多次隨機狀態，嘗試得到更優的解（其過程與模擬退火相似）。
模擬退火是個隨機的算法，執行次數越多獲得的解越有可能更優，所以我們可以執行多遍 SA 函數。至於如何控制時間？
while((double)clock()/CLOCKS_PER_SEC < 0.90) SA();
上面這個代碼控制時間在 \(0.90s\) 左右，如果時間限制為 \(1s\)，而每次 SA 函數執行時間略長時，就要小心可能會 \(\text{T}\) 掉了。
如果一個代碼不行，就考慮換個種子吧。
srand(...);
為了獲得更精確的解，也可以把 \(d\) 和 \(T_0\) 調的更精准一點
const double d = 0.996 -> 0.99996;
const double lim = 1e-10 -> 1e-15;
還有，隨機亂搞一些 初溫，終溫，降溫系數 也是可以的。
隨機微調：
對於序列分組這類問題，每次可以隨機兩個位置交換計算新解的貢獻。
但要注意如果沒有接受這個狀態時要把兩個位置在交換回去。
打亂排序：
對於一些統計貢獻與序列順序無關的題可以考慮。
函數 random_shuffle() 可以在 \(O(n)\) 的時間內隨機打亂一個序列
包含在 ctime 中，使用方法和 sort 類似。
可以用下面那道 P2503 來練練手。
Tips
模擬退火通常用來解決一些數據范圍比較小，但對正解的計算不好算的問題。
通過隨機一些狀態來不斷逼近正解。
其計算過程通常是 簡單粗暴 的。
其正確性可能是 玄學 的。
這些例題中都有體現，注意體會。
例題
UVA10228 A Star not a Tree?
初始點的位置考慮選擇所有點的平均值（人類本質？）
按照上面說的每次隨機一個點即可。
如何計算一個點的貢獻？\(n\) 的數據范圍只有 \(100\)，暴力計算即可。
看一下代碼可能會更好理解。
/*
Work by: Suzt_ilymics
Problem: 不知名屑題
Knowledge: 垃圾算法
Time: O(能過)
CAQ yyds !!!!
*/
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<queue>
#include<ctime>
#define LL long long
#define orz cout<<"lkp AK IOI!"<<endl
using namespace std;
const int MAXN = 1e5+5;
const int INF = 1e9+7;
const int mod = 1e9+7;
const double d = 0.996;
const double lim = 1e-10;
int n;
double ax[MAXN], ay[MAXN];
double limx, limy; // 上一次選擇保留的位置
double ansx, ansy; // 存答案的位置;
double ans;
int read(){
int s = 0, f = 0;
char ch = getchar();
while(!isdigit(ch)) f |= (ch == '-'), ch = getchar();
while(isdigit(ch)) s = (s << 1) + (s << 3) + ch - '0' , ch = getchar();
return f ? -s : s;
}
double dis(double sx, double sy, double ex, double ey) {
return sqrt((sx - ex) * (sx - ex) + (sy - ey) * (sy - ey));
}
double calc(double sx, double sy) {
double sum = 0.0;
for(int i = 1; i <= n; ++i) {
sum += dis(sx, sy, ax[i], ay[i]);
}
return sum;
}
void SA() {
double T = 2021; // 設置初始溫度
while(T > lim) { // 當沒降完溫繼續降溫
double nowx = limx + ((rand() << 1) - RAND_MAX) * T; // 獨特的隨機一個點的公式
double nowy = limy + ((rand() << 1) - RAND_MAX) * T;
double now = calc(nowx, nowy);//計算貢獻
double del = now - ans; // 計算 改變量
if(del < 0) {
ansx = nowx, ansy = nowy;
limx = nowx, limy = nowy;
ans = now;
} else if(exp(-del / T) > (double)rand() / RAND_MAX) { // 概率的選擇接受這個解
limx = nowx, limy = nowy;
}
T *= d; // 降溫
}
}
void work() {
limx /= (1.0 * n); // 取一個平均值
limy /= (1.0 * n);
ansx = limx, ansy = limy;
ans = calc(limx, limy);
for(int i = 1; i <= 10; ++i) SA(); // 多執行幾遍 SA 函數
}
int main()
{
int t = read();
for(int i = 1; i < t; ++i) { // UVA 獨特的數據測試方式
limx = limy = ansx = ansy = 0.0;
ans = 0;
n = read();
for(int j = 1; j <= n; ++j) {
scanf("%lf%lf", &ax[j], &ay[j]);
limx += ax[j], limy += ay[j];
}
work();
printf("%.0lf\n\n", ans);
}
n = read();
for(int j = 1; j <= n; ++j) {
scanf("%lf%lf", &ax[j], &ay[j]);
limx += ax[j], limy += ay[j];
}
work();
printf("%.0lf\n", ans);
return 0;
}
P5544 [JSOI2016]炸彈攻擊1
肯定要先隨機一個點，考慮如何確定 \(R\)。
隨機嗎？
注意隨機一個量的時候只能有一個參照量，也就是說這里不能同時隨機點的位置和爆炸半徑。
因為我們無法簡單的確定一個決策依據。
人類本質？
能炸多少炸多少啊！
所以暴力找到最大爆炸范圍統計貢獻就行啊。
/*
Work by: Suzt_ilymics
Problem: 不知名屑題
Knowledge: 垃圾算法
Time: O(能過)
*/
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<queue>
#include<cstdlib>
#include<ctime>
#define LL long long
#define orz cout<<"lkp AK IOI!"<<endl
using namespace std;
const int MAXN = 1e5+5;
const int INF = 1e9+7;
const int mod = 1e9+7;
const double d = 0.996;
const double lim = 1e-10;
struct node {
int x, y, r;
}a[MAXN];
struct node2 {
int x, y;
}b[MAXN];
int n, m, R;
double limx, limy, ansx, ansy;
int ans;
int read(){
int s = 0, f = 0;
char ch = getchar();
while(!isdigit(ch)) f |= (ch == '-'), ch = getchar();
while(isdigit(ch)) s = (s << 1) + (s << 3) + ch - '0' , ch = getchar();
return f ? -s : s;
}
double GetR(double sx, double sy, int x) {
return sqrt((sx - a[x].x) * (sx - a[x].x) + (sy - a[x].y) * (sy - a[x].y)) - a[x].r;
}
double Dis(double sx, double sy, int x) {
return sqrt((sx - b[x].x) * (sx - b[x].x) + (sy - b[x].y) * (sy - b[x].y));
}
int calc(double sx, double sy) {
double maxR = 1.0 * R;
for(int i = 1; i <= n; ++i) {
maxR = min(maxR, GetR(sx, sy, i));
}
int cnt = 0;
for(int i = 1; i <= m; ++i) {
if(Dis(sx, sy, i) <= maxR) cnt++;
}
return cnt;
}
void SA() {
double T = 2021;
while(T > lim) {
double nowx = limx + ((rand() << 1) - RAND_MAX) * T;
double nowy = limy + ((rand() << 1) - RAND_MAX) * T;
int now = calc(nowx, nowy);
int del = now - ans;
if(del > 0) {
ansx = nowx, ansy = nowy;
limx = nowx, limy = nowy;
ans = now;
} else if(exp(-del/T) < (double) rand() / RAND_MAX) {
limx = nowx, limy = nowy;
}
T *= d;
}
}
void work() {
ansx = ansy = limx = limy = 0.0;
ans = calc(limx, limy);
while((double)clock() / CLOCKS_PER_SEC < 0.70) SA(); // CLOCKS_PER_SEC
}
int main()
{
n = read(), m = read(), R = read();
for(int i = 1; i <= n; ++i) a[i].x = read(), a[i].y = read(), a[i].r = read();
for(int i = 1; i <= m; ++i) b[i].x = read(), b[i].y = read();
work();
printf("%d", ans);
return 0;
}
P2503 [HAOI2006]均分數據
原題面中的式子可能有點錯誤，應該把 \(m\) 組當成 \(m\) 個元素來計算方差
介紹一個函數 random_shuffle()，可以在 \(O(n)\) 的時間內，將一段區間隨機打亂。
考慮最簡單的想法，
每次隨機打亂，貪心的分組。
每次加進一個元素，將它放進當前 \(m\) 組中最小的那一組中。
然后計算貢獻。
多隨機幾遍，總有一遍我們會碰到最優解的。
/*
Work by: Suzt_ilymics
Problem: 不知名屑題
Knowledge: 垃圾算法
Time: O(能過)
*/
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<queue>
#include<ctime>
#define LL long long
#define orz cout<<"lkp AK IOI!"<<endl
using namespace std;
const int MAXN = 1e5+5;
const int INF = 1e9+7;
const int mod = 1e9+7;
const double d = 0.996;
const double lim = 1e-10;
int n, m;
int a[22];
int b[7];
int sum;
double ave, ans = 10000000.0;
int read(){
int s = 0, f = 0;
char ch = getchar();
while(!isdigit(ch)) f |= (ch == '-'), ch = getchar();
while(isdigit(ch)) s = (s << 1) + (s << 3) + ch - '0' , ch = getchar();
return f ? -s : s;
}
double calc() {
memset(b, false, sizeof b);
for(int i = 1; i <= n; ++i) {
int wz = 1;
for(int j = 1; j <= m; ++j) if(b[wz] > b[j]) wz = j;
b[wz] += a[i];
}
double Sum = 0;
for(int i = 1; i <= m; ++i) {
Sum += (ave - 1.0 * b[i]) * (ave - 1.0 * b[i]);
}
Sum /= (1.0 * m);
return sqrt(Sum);
}
void SA() {
double T = 2021;
while(T > lim) {
random_shuffle(a + 1, a + n + 1);
double now = calc();
ans = min(ans, now);
// double del = now - ans;
//// cout<<now<<endl;
// if(del < 0) {
// ans = now;
// }
T *= d;
}
}
void work() {
ave = 1.0 * sum / m;
for(int i = 1; i <= 20; ++i) SA();
}
int main()
{
n = read(), m = read();
for(int i = 1; i <= n; ++i) a[i] = read(), sum += a[i];
work();
printf("%.2lf", ans);
return 0;
}
P3878 [TJOI2010]分金幣
前一半分一組，后一半分一組。
然后每次隨機兩個位置交換一下計算貢獻。
nowx = rand() % n + 1; // 整數位置的隨機直接 % n + 1 即可
這里與板子略有不同
如果更優就更新，
否則考慮概率接受，
否則你要把元素再換回去！
計算函數應該很好想，依舊是暴力，如果還不會就看代碼吧。
/*
Work by: Suzt_ilymics
Problem: 不知名屑題
Knowledge: 垃圾算法
Time: O(能過)
*/
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<queue>
#include<ctime>
#define int long long
#define LL long long
#define orz cout<<"lkp AK IOI!"<<endl
using namespace std;
const int MAXN = 1e5+5;
const int INF = 1e12+7;
const int mod = 1e9+7;
const double d = 0.996;
const double lim = 1e-10;
int t, n, ans = INF;
int a[40];
int read(){
int s = 0, f = 0;
char ch = getchar();
while(!isdigit(ch)) f |= (ch == '-'), ch = getchar();
while(isdigit(ch)) s = (s << 1) + (s << 3) + ch - '0' , ch = getchar();
return f ? -s : s;
}
int calc() {
int sum1 = 0, sum2 = 0;
for(int i = 1; i <= n / 2; ++i) sum1 += a[i];
for(int i = n / 2 + 1; i <= n; ++i) sum2 += a[i];
return abs(sum1 - sum2);
}
void SA() {
double T = 2021;
while(T > lim) {
int nowx = rand() % n + 1;
int nowy = rand() % n + 1;
swap(a[nowx], a[nowy]);
int now = calc();
int del = now - ans;
if(del < 0) {
ans = now;
} else if(exp(-del/T) > (double)rand() / RAND_MAX) {
} else swap(a[nowx], a[nowy]);
T *= d;
}
}
void work() {
for(int i = 1; i <= 100; ++i) SA();
}
signed main()
{
t = read();
while(t--) {
n = read();
ans = INF;
for(int i = 1; i <= n; ++i) a[i] = read();
random_shuffle(a + 1, a + n + 1);
work();
printf("%lld\n", ans);
}
return 0;
}
P3936 Coloring
選擇按順序染色作為初始狀態。
然后開始微調。
如何更快的計算貢獻？
在原來的基礎上減去兩個位置原來位置的貢獻，然后加上交換后的貢獻。
推薦參數：
d = 0.9999，lim = 1e-15，T = 500，srand(20041011)
這樣的話分數就可以在 \(97pts\) 以上了
有個疑問？
我計算 \(del\) 的時候是當前值與交換前的值作差，而不是和最優值作差。可前面幾道題都是和最優值作差，可這道題這樣做卻死活過不去。有沒有大佬可以解釋一下啊，也可以到我這個帖子看一看。
/*
Work by: Suzt_ilymics
Problem: 不知名屑題
Knowledge: 垃圾算法
Time: O(能過)
*/
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<queue>
#include<ctime>
#define LL long long
#define orz cout<<"lkp AK IOI!"<<endl
using namespace std;
const int MAXN = 1e5+5;
const int INF = 1e9+7;
const int mod = 1e9+7;
const double d = 0.99998;
const double lim = 1e-15;
int dx[] = {1, -1, 0, 0};
int dy[] = {0, 0, 1, -1};
int n, m, c, ans, now_;
int p[55], cnt[55], top = 1;
int col[22][22], acol[22][22];
int read(){
int s = 0, f = 0;
char ch = getchar();
while(!isdigit(ch)) f |= (ch == '-'), ch = getchar();
while(isdigit(ch)) s = (s << 1) + (s << 3) + ch - '0' , ch = getchar();
return f ? -s : s;
}
bool check(int x, int y) { return x < 1 || y < 1 || x > n || y > m; }
int Dev(int sx, int sy) {
int Cnt = 0;
for(int i = 0; i < 4; ++i) {
int x = sx + dx[i], y = sy + dy[i]; // 找相鄰的四個格子計算貢獻
if(check(x, y)) continue;
if(col[sx][sy] != col[x][y]) Cnt++;
}
return Cnt;
}
//int Dev(int x, int y) { return (col[x][y] != col[x-1][y] && x!=1) + (col[x][y] != col[x][y-1] && y!=1) + (col[x][y] != col[x+1][y] && x != n) + (col[x][y] != col[x][y+1] && y != m); }
int calc(int sx, int sy, int ex, int ey) {
int res = 0;
// orz;
res = res - Dev(sx, sy) - Dev(ex, ey); // 減去原來的貢獻
swap(col[sx][sy], col[ex][ey]);
res = res + Dev(sx, sy) + Dev(ex, ey); // 加上現在的貢獻
return res;
}
void SA() {
// orz;
// for(int i = 1; i <= n; ++i) for(int j = 1; j <= m; ++j) col[i][j] = acol[i][j];
// now_ = ans;
double T = 500;
while(T > lim) {
int nowsx = rand() % n + 1;
int nowsy = rand() % m + 1;
int nowex = rand() % n + 1;
int nowey = rand() % m + 1;
// int pre = Calc();
// swap(col[nowsx][nowsy], col[nowex][nowey]); // 先交換
// cout<<calc(nowsx, nowsy, nowex, nowey)<<"\n";
int now = now_ + calc(nowsx, nowsy, nowex, nowey); // 再計算
// cout<<ans<<" "<<now_<<" "<<now<<"\n";
int del = now - now_; // now_ 存的是當前 col 矩陣的答案
// int del = now - ans; // ans 存的是最優解
if(del < 0) { // 更新這個解
ans = now;
now_ = now;
for(int i = 1; i <= n; ++i) for(int j = 1; j <= m; ++j) acol[i][j] = col[i][j];
} else if(exp(-del/T) > (double) rand() / RAND_MAX) {
now_ = now; // 一定概率接受當前這個狀態
} else {
swap(col[nowsx][nowsy], col[nowex][nowey]); // 交換回去
}
T *= d;
}
}
void work() {
// for(int i = 1; i <= n; ++i) {
// for(int j = 1; j <= m; ++j) {
// ans += Dev(i, j);
// }
// }
// ans /= 2;
for(int i = 2; i <= n; ++i)
for(int j = 1;j <= m; ++j)
if(col[i][j] != col[i - 1][j]) ans++;
for(int i = 1; i <= n; ++i)
for(int j = 2;j <= m; ++j)
if(col[i][j] != col[i][j - 1]) ans++;
now_ = ans;
// cout<<ans<<" "<<now_<<endl;
// cout<<ans<<endl;
// for(int i = 1; i <= 30; ++i) SA();
while((double)clock() / CLOCKS_PER_SEC < 3.90) SA();//, cout<<ans<<" "<<now_<<endl;
}
int main()
{
srand(20041011);
n = read(), m = read(), c = read();
for(int i = 1; i <= c; ++i) p[i] = read();
for(int i = 1; i <= n; ++i) {
for(int j = 1; j <= m; ++j) {
if(cnt[top] == p[top]) ++top;
acol[i][j] = col[i][j] = top;
++cnt[top];
}
}
work();
// printf("%d\n", ans);
for(int i = 1; i <= n; ++i) {
for(int j = 1; j <= m; ++j) {
printf("%d ", acol[i][j]);
}puts("");
}
return 0;
}
/*
int dev(int sx, int sy, int ex, int ey) {
int cnt = 0;
for(int i = 0; i < 4; ++i) {
int x = sx + dx[i], y = sy + dy[i]; // 與起點相鄰的點
if(check(x, y)) continue;
if(col[ex][ey] != col[x][y]) cnt++;
}
for(int i = 0; i < 4; ++i) {
int x = ex + dx[i], y = ey + dy[i]; // 與終點相鄰的點
if(check(x, y)) continue;
if(col[sx][sy] != col[x][y]) cnt++;
}
// cout<<"cnt1 :"<< cnt<<" ";
return cnt;
}
int Calc() {
int cnt = 0;
for(int i = 2; i <= n; ++i)
for(int j = 1;j <= m; ++j)
if(col[i][j] != col[i - 1][j]) cnt++;
for(int i = 1; i <= n; ++i)
for(int j = 2;j <= m; ++j)
if(col[i][j] != col[i][j - 1]) cnt++;
return cnt;
}
*/