模拟退火

本篇将介绍模拟退火算法

---

引入

聪明的读者啊，在你绞尽脑汁寻找完成题目的算法的时候，是否想过大力出奇迹，直接通过随机的方法来猜出答案。模拟退火是一种随机化算法，当我们遇到一个最优解问题，并且答案不是一个单峰函数（单峰可以使用爬山算法，直接靠近每个比当前值大的即可），即可使用模拟退火求解。（如果方案数量过大，有可能会受困于时间复杂度）

简述

我们参考爬山算法，当遇到比当前答案更优的解则修改答案。与其不同的是，由于我们不是单峰函数，遇到非更优的解的时候，我们选择概率修改答案，来跳出当前局部最优解，尝试得到更好的全局最优解。但是大多数情况一次模拟退火得不到最优解，我们会多次进行模拟退火来尝试得到最优解。

针对跳转的概率，我们定义当前温度为 $T$ ，新状态 $S'$ 和目前最优解状态 $S$ 以及他们的差 $\Delta E$ 。我们可以将修改答案的概率描述为：

$p=\left\{\begin{matrix} 1 &S'优于S \\ e^{\frac{-\Delta E}{T} } &否则 \end{matrix}\right.$

构造退火函数

我们定义初始温度 $T_0$ ，降温系数 $d$ ，终止温度 $T_k$ 。其中 $T_0$ 是一个比较大的数，$d$ 是一个小于 $1$ 但是非常接近 $1$ 的数， $T_k$ 是一个接近 $0$ 的正数。

让温度 $T=T_0$ ，然后进行一次修改尝试，再让 $T=d\times T$ ，当 $T<T_k$ 的时候模拟退火结束。

void sa() {
int oldval = calc(); //通过计算获取默认值
for(double t = 1e4;t >= 1e-4;t *= 0.99) {
     ... //进行随机交换，来尝试别的可能
     int newval = calc(); //通过计算获取当前值
     ans = max(ans, newval); //更新最优解答案
     int dt = newval - oldval;
     if (exp(dt / t) < double(rand()) / RAND_MAX) swap(a[x], a[y]); //取消更新，回溯
     else oldval = newval; //保存新答案
 }
}

实现

我们以 洛谷P1559 为例，$ans$ 为答案，$n$ 为人数，$p$ 和 $q$ 为优势值，$a$ 为当前的女生和哪位男生进行配对

1. 制作一个calc函数去计算，在当前匹配情况下优势的总和
2. 在sa函数中随机交换的位置，随机选择数组中的两位进行交换
3. 在主函数中通过double(clock()) / CLOCKS_PER_SEC < 0.8的时间来尽量多的跑模拟退火

代码

#include <bits/stdc++.h>
using namespace std;
const int N = 21;

int ans, n, p[N][N], q[N][N];
vector<int> a;

int calc() {
    int res = 0;
    for(int i = 0;i < n;i++) {
        int x = a[i];
        res += p[x][i+1] * q[i+1][x];
    }
    return res;
}

void sa() {
    int oldval = calc();
    for(double t = 1e4;t >= 1e-4;t *= 0.99) {
        int x = rand() % n, y = rand() % n;
        swap(a[x], a[y]);
        int newval = calc();
        ans = max(ans, newval);
        int dt = newval - oldval;
        if (exp(dt / t) < double(rand()) / RAND_MAX) swap(a[x], a[y]);
        else oldval = newval;
    }
}

signed main(){
    cin >> n;
    for(int i = 1;i <= n;i++) a.push_back(i);
    for(int i = 1;i <= n;i++) for(int j = 1;j <= n;j++) cin >> p[i][j];
    for(int i = 1;i <= n;i++) for(int j = 1;j <= n;j++) cin >> q[i][j];
    while(double(clock()) / CLOCKS_PER_SEC < 0.8) sa();
    cout << ans << endl;
    return 0;
}