先开个坑，期中考完填。
模拟退火的概念
退火原先是一个物理学上的概念，当固体的温度很高时，内能比较大，固体内的粒子处于快速无序运动状态，当温度慢慢降低，固体的内能减小，粒子逐渐趋于有序，最终固体处于常温状态，内能达到最小，此时粒子最为稳定。
那么模拟退火算法 $(\text{Simulated Annealing})$ 的思想就是指从某一较高初温出发，伴随温度参数的不断下降,结合概率突跳特性在解空间中随机寻找目标函数的全局最优解，即在局部最优解能概率性地跳出并最终趋于全局最优。
简单来说，模拟退火算法就是一个聪明的爬山算法，因为即使有时新搜到的解劣于已经得到的最优解，但我们还是有一定概率将最优解变为这个更劣解。而爬山算法的缺点很明显——容易被限制在一个次优的峰值上，模拟退火即能帮助我们跳到其它的峰值。
过程
模拟退火的过程就是如上的过程，温度较高的时候，红线抖动的范围会更广，更有可能找到全局最优峰值，温度较低时，红线向峰值收敛，最终得到最优解。
那么现在我们只需要解决一个问题：如何设计更新状态的概率？
由上述的概念可知，我们需要接受以一定概率去更换到最优解的情况，但这个概率具体是多少呢？
这个公式也有一个名字，叫Metropolis接受准则，可自行前往百度百科查看，我们可以总结成如下公式：
- 当新解优于当前最优解时，$p=1$，即一定将其更新为最优解。
- 当新解劣于当前最优解时，$p=\exp(\dfrac{-\Delta E}{T})$，其中的 $\exp$ 表示自然对数。
那么这里的 $\Delta E$ 和 $T$ 怎么解释？
$T$ 即为上文所提到的温度，$\Delta E$ 代表新解与旧解之间的变化量。
然后我们就解决了问题，每次随机一个情况，计算此时的代价，修改答案，直到退火结束。
参数
模拟退火算法的一个麻烦的地方。
一个适合的参数对题目的分数产生巨大的影响，我们要给模拟退火算法制定哪些参数呢？
初始温度 $T_0$
初始温度即为物体一开始时的温度，物体内的粒子最活泼的时刻，一般情况下，一个越大的初始温度搜索到最优解的可能性越大。
降温系数 $d$
降温系数顾名思义，就是物体降温速度的快慢，可以粗略的得知 $d\in (0,1)$，但一般情况下，$d$ 要取 $0.999$ 附近，注意 $d$ 对 $T$ 的影响是指数的，所以任何 $10^{-3}$ 量级下的改变也会对降温速度的快慢产生很大的影响，一般按照具体题目分析。
终止温度 $T_k$
物体停止降温的时刻，此时物体内的粒子趋于稳定，也是我们找到最优解的时刻，一般取一个接近 $0$ 的数。
理论上只有这 $3$ 个函数，但是还有一个玄学的因素：
随机种子 $seed$！
由于中途需要随机情况，所以一个好的种子也是得到答案的关键，只不过会大大增加模拟退火算法的不确定性，一般以时间作为种子，当然有比较好的种子也行。
当然，具体的参数还是要在具体的题目分析，被卡除了算法的问题就是参数问题，此时要有足够的耐心。
技巧
模拟退火算法还有着几个帮助解题的技巧。
分值域模拟退火
当函数的峰点很多时，模拟退火难以直接跑出最优解。
此时可以把整个值域分成几段，每段跑一遍模拟退火，然后再取最优解。
卡时
在题目的时限要求内尽可能多跑几次模拟退火，增加答案最优的可能性。
具体可写为：
lf TimeLimit = 0.9; // 单位为秒，略小于题目时限。
while((lf)clock() / CLOCKS_PER_SEC < TimeLimit) {
SimulateAnneal();
}
题目
P1337 [JSOI2004] 平衡点 / 吊打XXX
我们通过物理学上的定义，计算能量总和，能量总和越小越稳定。
lf getEnergy(lf px, lf py) {
lf dist = 0;
for(int i = 1;i <= n;i++) {
lf dx = x[i] - px, dy = y[i] - py;
dist += sqrt(dx * dx + dy * dy) * w[i];
}
if(res > dist) {
res = dist, ansx = px, ansy = py;
}
return dist;
}
因为我们是以坐标 $(x,y)$ 向周围随机拓展，根据温度越高随机范围越大，我们可以用如下方法随机新一轮的坐标。
lf nxtx = nowx + Rand(-1e4, 1e4) * t0, nxty = nowy + Rand(-1e4, 1e4) * t0;
当温度降到终止温度后，我们可以再随机 $T$ 次，求得最终的解。
代码
#include <bits/stdc++.h>
using namespace std;
typedef long double lf;
int n, x[1005], y[1005], w[1005];
lf ansx, ansy, res;
mt19937 engine(chrono::steady_clock::now().time_since_epoch().count());
int Rand(int l, int r) {
uniform_int_distribution<int> dist(l, r);
return dist(engine);
}
lf getEnergy(lf px, lf py) {
lf dist = 0;
for(int i = 1;i <= n;i++) {
lf dx = x[i] - px, dy = y[i] - py;
dist += sqrt(dx * dx + dy * dy) * w[i];
}
if(res > dist) {
res = dist, ansx = px, ansy = py;
}
return dist;
}
void SimulateAnneal() {
lf t0 = 10000, tk = 1e-8, d = 0.999;
lf nowx = ansx, nowy = ansy;
while(t0 > tk) {
lf nxtx = nowx + Rand(-1e4, 1e4) * t0, nxty = nowy + Rand(-1e4, 1e4) * t0;
lf delta = getEnergy(nxtx, nxty) - getEnergy(nowx, nowy);
if(exp(-delta / t0) > (lf)Rand(0, 1e4) / 1e4) {
nowx = nxtx, nowy = nxty;
}
t0 *= d;
}
int T = 5000;
for(int test = 1;test <= T;test++) {
lf nxtx = ansx + Rand(-1e4, 1e4) * t0, nxty = ansy + Rand(-1e4, 1e4) * t0;
getEnergy(nxtx, nxty);
}
}
int main() {
ios::sync_with_stdio(0);
cin.tie(0), cout.tie(0);
cin >> n;
for(int i = 1;i <= n;i++) {
cin >> x[i] >> y[i] >> w[i];
ansx += x[i], ansy += y[i];
}
ansx /= n, ansy /= n, res = getEnergy(ansx, ansy);
SimulateAnneal();
SimulateAnneal();
cout << fixed << setprecision(3) << ansx << ' ' << ansy << '\n';
cout.flush();
return 0;
}
开始的时候跑一遍 SimulateAnneal() 错了，改成了跑两遍就过了。
P2210 Haywire
我们对于答案的最优排列做模拟退火，本题的能量计算就是对于一个排列的答案。
随机的时候就没有能用温度来控制变化范围的需要，我们随机两个点 $x,y$，交换 $now_x,now_y$ 作为新一轮的排列。
代码
#include <bits/stdc++.h>
using namespace std;
typedef long double lf;
int n, fri[15][3];
int ans[15], res;
mt19937 engine(chrono::steady_clock::now().time_since_epoch().count());
int Rand(int l, int r) {
uniform_int_distribution<int> dist(l, r);
return dist(engine);
}
int getEnergy(int p[]) {
int sum = 0;
for(int i = 1;i <= n;i++) {
for(int j = 0;j < 3;j++) {
sum += abs(p[i] - p[fri[i][j]]);
}
}
if(res > sum) {
res = sum;
memcpy(ans, p, sizeof(ans));
}
return sum;
}
void SimulateAnneal() {
lf t0 = 10000, tk = 1e-8, d = 0.999;
int now[15], nxt[15], x, y;
memcpy(now, ans, sizeof(now));
while(t0 > tk) {
memcpy(nxt, now, sizeof(nxt));
do {
x = Rand(1, n), y = Rand(1, n);
} while(x == y);
swap(nxt[x], nxt[y]);
double delta = getEnergy(nxt) - getEnergy(now);
if(exp(-delta / t0) > (lf)Rand(0, 1e4) / 1e4) {
memcpy(now, nxt, sizeof(now));
}
t0 *= d;
}
int T = 10000;
for(int i = 1;i <= T;i++) {
random_shuffle(now + 1, now + n + 1);
getEnergy(now);
}
}
int main() {
ios::sync_with_stdio(0);
cin.tie(0), cout.tie(0);
cin >> n;
for(int i = 1;i <= n;i++) {
cin >> fri[i][0] >> fri[i][1] >> fri[i][2];
ans[i] = i;
}
res = getEnergy(ans);
SimulateAnneal();
SimulateAnneal();
cout << res / 2 << '\n';
cout.flush();
return 0;
}
P3878 [TJOI2010]分金币
首先将前 $\dfrac{n}{2}$ 个分为一组，剩下 $(n-\dfrac{n}{2})$ 个为一组。
每次在两组间随机挑选一对换一下跑模拟退火。
思路十分简单，但参数的调整十分困难。
由于本题是多测，所以不好控制每个询问的时间，直接对于每组询问跑了 $450$ 遍。
代码
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef long double lf;
int T, n;
ll a[30];
ll seq[30];
ll ans = 1e18;
mt19937 engine(chrono::steady_clock::now().time_since_epoch().count());
int Rand(int l, int r) {
uniform_int_distribution<int> dist(l, r);
return dist(engine);
}
inline ll Abs(ll x) {
return x > 0 ? x : (-x);
}
ll getEnergy(ll now[]) {
ll sum = 0;
for(int i = 0;i < n / 2;i++) {
sum += now[i];
}
for(int i = n / 2;i < n;i++) {
sum -= now[i];
}
sum = Abs(sum);
if(ans > sum) {
ans = sum;
}
return sum;
}
void SimulateAnneal() {
lf t0 = 2000, tk = 1e-10, d = 0.95;
ll now[30], nxt[30];
memcpy(now, a, sizeof(now));
while(t0 > tk) {
int x = Rand(1, n / 2) - 1, y = Rand(n / 2 + 1, n) - 1;
memcpy(nxt, now, sizeof(nxt));
swap(nxt[x], nxt[y]);
ll delta = getEnergy(nxt) - getEnergy(now);
if(exp(-delta / t0) > (lf)Rand(0, 1e4) / 1e4) {
memcpy(now, nxt, sizeof(now));
}
t0 *= d;
}
}
int main() {
ios::sync_with_stdio(false);
cin.tie(nullptr), cout.tie(nullptr);
cin >> T;
while(T--) {
cin >> n;
for(int i = 0;i < n;i++) {
cin >> a[i];
}
if(n == 1) {
cout << a[0] << '\n';
continue;
}
ans = getEnergy(a);
for(int test = 1;test <= 450;test++) {
random_shuffle(a, a + n);
SimulateAnneal();
}
cout << ans << '\n';
}
cout.flush();
return 0;
}
P7962 [NOIP2021] 方差
这道题就是我想在 $\text{NOIP}$ 前学模拟退火的原因，基于一个正常的考场情况，我们肯定是希望在 $1h$ 内拿到这题的一个比较高的分数然后花时间冲 $\text{T4}$ 暴力，至于使用模拟退火通过那大可不必，因为：
- 在数据稍微加强一些的情况下，模拟退火是大概率过不去后面的点的。
- 这篇题解所提到的模拟退火正解已经用到了此题几乎全部的性质，那么我相信了解了所有性质后会更倾向于写正确性有保证的动态规划。
所以最后的程序十分简单：观察出交换差分的性质后，随机找两个位置交换跑模拟退火，洛谷上可以得 $72$ 分，相比于观察差分单谷的性质，这个做法会更暴力无脑一些，适合一些一般同学性价比较高的做法。
代码
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef long double lf;
int n, a[405];
int D[405];
ll ans;
mt19937 engine(chrono::steady_clock::now().time_since_epoch().count());
int Rand(int l, int r) {
uniform_int_distribution<int> dist(l, r);
return dist(engine);
}
ll getEnergy(int now[]) {
ll res = 0, sum = 0, s = 0;
for(int i = 1;i <= n;i++) {
sum += now[i];
res += 1ll * n * sum * sum;
s += sum;
}
res -= s * s;
if(ans > res) {
ans = res;
}
return res;
}
void SimulateAnneal() {
lf t0 = 3000, tk = 1e-8, d = 0.99;
int now[405], nxt[405], x, y;
memcpy(now, D, sizeof(now));
while(t0 > tk) {
do {
x = Rand(2, n), y = Rand(2, n);
} while(x == y);
memcpy(nxt, now, sizeof(nxt));
swap(nxt[x], nxt[y]);
lf delta = getEnergy(nxt) - getEnergy(now);
if(exp(-delta / t0) >= (lf)Rand(0, 1e4) / 1e4) {
memcpy(now, nxt, sizeof(now));
}
t0 *= d;
}
}
int main() {
ios::sync_with_stdio(false);
cin.tie(nullptr), cout.tie(nullptr);
cin >> n;
for(int i = 1;i <= n;i++) {
cin >> a[i];
D[i] = a[i] - a[i - 1];
}
ans = getEnergy(D);
lf TimeLimit = 0.96;
while((lf)clock() / CLOCKS_PER_SEC < TimeLimit) {
SimulateAnneal();
}
cout << ans << '\n';
cout.flush();
return 0;
}
总结
模拟退火终究是乱搞，在赛时争取高分时可以视情况使用，但在算法学习中不建议乱搞，所以大概率不会再往里加题，除非有什么模拟退火的优化技巧或者场上有用模拟退火拿高分的题目会更新到这里。