【Coel.学习笔记】随机化算法：模拟退火与爬山法

简介

模拟退火（ $\text{Simulate Anneal}$ ）和爬山法是随机化算法，二者的原理都在于通过随机生成答案并检查，把答案逐步缩小在一个可行的区间，尽可能地靠近正确答案。
在考场中，如果某道题目的正解比较难想（例如性质复杂的贪心）或者对应的函数取值可能数极大且难以二分/三分求解时，模拟退火和爬山法可以一定程度上得到正确答案并获得可观的分数。
顺带一提，现在的人工智能其实也是基于随机化算法的，例如遗传算法就是一种典型的随机化算法。
~~高情商：利用人工智能算法获取近似解~~
~~低情商：随机骗分~~

原理

下面主要讲模拟退火，爬山法只提一点。

爬山法

爬山法通常用来求解一个函数的最值，当然得是一个凸函数。

可能有人要问了，既然是凸函数，为什么不用三分？

一方面，某些时候看不出是凸函数，爬山法硬上有时能过；另一方面，如果函数有多维，就要写嵌套三分，三维套两个，四维套三个……如果有 $k$ 维（比如下面这一题），即使耐得住性子写完，时间复杂度也是 $O(n\log ^k n)$ ，效率很低。

回到爬山法。爬山法的原理可以总结成“从一个胜利走向另一个胜利”~~（经典英语造句）~~，即每次寻找一个最优解，按照当前给定的参数控制爬山距离，并根据当前局面选定爬山的方向。当然有可能跳过最优解或者陷在某一个局部最优解的位置，所以爬山法的局限性比较大。

模拟退火

顾名思义，模拟退火就是模拟物理学上的退火过程。类比分子在高温时能量更高、碰撞次数更大的原理（参见《化学反应原理》P28），我们定义几个量：

温度 $T$ ：相当于答案区间，控制每次随机区域。一开始的 $T$ 为较大值，接下来逐步减小。
$T$ 有初始态 $T_0$ 、终止态 $T_k$ 和降温系数 $d$ 三个关联量。一般来说 $T$ 取值较大， $T_k$ 接近 $0$ ， $d$ 略小于 $1$ 。这样每次让 $T$ 乘以 $d$ 就可以实现“降温”，让答案稳定下来。
能量差 $\Delta E$ ，表示新得到答案和原本答案的差值。如果 $\Delta E<0$ ，说明答案更优，应取该答案；反之，以一定概率取这个答案，从而减小落入局部最优解的可能性。
通常取概率值为 $e^{\frac{-\Delta E}{T}}$ ，这样可以保证能量差越高取答案的可能越大。

下面这个图反映了落入局部最优解的坏处（来自 OI-Wiki）。

例题讲解

[JSOI2008]球形空间产生器 - 爬山法浅谈

洛谷传送门
给定 $n+1$ 个 $n$ 维空间上的点，求经过这些点的球的球心坐标。

解析：这道题正解是用高斯消元求解 $n$ 元一次方程组（球的方程可以通过作差消掉高次项），效率确实很高，但使用爬山法求解写起来更直观。

具体操作为：把球心初始化成这些点的重心（即坐标均值），然后每次爬山时，通过尽可能让距离平均来寻找球心，远了就拉近一点，近了就推远一点。

当然这道题也可以用模拟退火做，这里只是简单介绍一下爬山法。事实上爬山法比较鸡肋，一般不会用quq

#include <iostream>
#include <iomanip>
#include <cmath>

using namespace std;

const int maxn = 20;

int n;
double dis[maxn][maxn], cans[maxn], ans[maxn], delta[maxn]; 

void hillClimbing() {
	double sum = 0;
	for (int i = 1; i <= n + 1; i++) {
		cans[i] = delta[i] = 0;
		for (int j = 1; j <= n; j++)
			cans[i] += (dis[i][j] - ans[j]) * (dis[i][j] - ans[j]);
		cans[i] = sqrt(cans[i]);
		sum += cans[i] / (n + 1); //求出距离均值
	}
	for (int i = 1; i <= n + 1; i++)
		for (int j = 1; j <= n; j++)
			delta[j] += (cans[i] - sum) * (dis[i][j] - ans[j]) / sum;
}

int main(void) {
	ios::sync_with_stdio(false);
	cin.tie(nullptr);
	cin >> n;
	for (int i = 1; i <= n + 1; i++)
		for (int j = 1; j <= n; j++) {
			cin >> dis[i][j];
			ans[j] += dis[i][j] / (n - 1); //预处理重心
		}
	for (double step = 1e4; step > 1e-6; step *= 0.99995) {
		hillClimbing();
		for (int i = 1; i <= n; i++)
			ans[i] += delta[i] * step;
	}
	for (int i = 1; i <= n; i++)
		cout << fixed << setprecision(3) << ans[i] << ' ';
	return 0;
}

[POJ2420]A Star not a Tree? - 模拟退火入门

洛谷传送门
给定平面直角坐标系上的 $n$ 个点，找到一个点使得到这些点的距离最小（即费马点），输出距离之和（保留整数）。

解析：这道题的答案是一个凸函数，可以用三分法求解，这里演示一下模拟退火。
模拟退火写起来其实非常无脑xwx

#include <iostream>
#include <algorithm>
#include <cmath>
#include <random>
#include <iomanip>

using namespace std;

const int maxn = 150;
const double maxTime = 0.95;

random_device rd;
mt19937 Rand(rd());
uniform_real_distribution<double> Dis(0, 1); //随机三件套

struct Point {
    double x, y;
} q[maxn];

int n;
double ans;

inline double rand(double l, double r) { //得到 [l, r) 的数字
    return (double) Dis(Rand) * (r - l) + l;
}

inline double dis(Point a, Point b) {
    double dx = a.x - b.x, dy = a.y - b.y;
    return sqrt(dx * dx + dy * dy);
}

double calcDistance(Point p) { //计算距离和，顺便更新答案
    double res = 0;
    for (int i = 1; i <= n; i++) res += dis(p, q[i]);
    ans = min(ans, res);
    return res;
}

void simulateAnneal() {
    Point cur = {rand(0, 1e4), rand(0, 1e4)};
    for (double T = 1e4; T > 1e-4; T *= 0.99) {
        Point p = {rand(cur.x - T, cur.x + T), rand(cur.y - T, cur.y + T)}; //在范围内随机取一个新点
        double delta = calcDistance(p) - calcDistance(cur);
        if (exp(-delta / T) > rand(0, 1)) cur = p; //概率取解
    }
}

int main(void) {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    int T;
    cin >> T;
    for (int j = 1; j <= T; j++) {
        ans = 1e8;
        cin >> n;
        for (int i = 1; i <= n; i++) cin >> q[i].x >> q[i].y;
        for (int i = 1; i <= 100; i++) simulateAnneal();
        /*如不是多组数据则可以利用“卡时”技巧，下一题会详细讲解*/
        cout << fixed << setprecision(0) << ans << '\n';
        if (j != T) cout << '\n';
    }
    return 0;
}

[JSOI/AHOI2014]保龄球 - 非数学的模拟退火

洛谷传送门
进行若干轮保龄球，对于每一轮，全倒时下一轮有双倍奖励分，补中时下一轮第一次有双倍奖励分，失误无奖励分。特别地，如果第 $n$ 轮全中则可以额外进行一个奖励轮且有双倍奖励分，奖励轮只能进行一次。

某人可以随意调整这 $n$ 轮的得分情况，但必须保证调整前后奖励轮的情况一致。求出在该操作下可以得到的最高分数。

解析：这道题看起来不像上一题那么“数学”，不过本质上都是求最优解；如果用动态规划或者搜索，会发现状态太多不好储存，且速度很慢，考虑使用模拟退火。事实上，之前提到的很多最优解问题（动态规划、二分、网络流……）都可以用模拟退火，只不过答案不够精确罢了。

不难发现，这道题的调整为保龄球情况的排列，对于排列上的模拟退火，我们可以随机选取两个数据交换一下，把它作为新的状态。如果状态合法，就按照解的优劣，以一定概率更新答案；否则把位置交换回去，还原现场。此时温度 $T$ 的“区间长度”意义被淡化，仅仅用来控制退火次数以及取答案的概率。

由于模拟退火次数比较难估计，这里采用卡时技巧。C++11 的chrono 头文件提供了一个叫 chrono::system_clock::now().time_since_epoch().count() 的函数（好长的名字……），可以得到以纳秒为单位的时间，除以 $10^9$ 就是秒为单位了。

时间上限可以设置为一个略小于时间限制的值（不能等于，否则还没输出就 TLE 了），实际发现数据很水， $0.05$ 秒也可以过。

#include <iostream>
#include <cmath>
#include <chrono>
#include <random>
#include <vector>

using namespace std;

const int maxn = 80;

int n, ans;
bool bonus;

random_device rd; 
mt19937 Rand(rd());
uniform_real_distribution<double> check(0, 1);

pair<int, int> a[maxn];

inline auto getTime() { //得到秒为单位的时间（函数名有点长……）
    return chrono::system_clock::now().time_since_epoch().count() / 1e9;
}

int newRecord() { // 计算当前得分
    int res = 0;
    for (int i = 1; i <= n + bonus; i++) {
        res += a[i].first + a[i].second;
        if (i <= n) {
            if (a[i].first == 10) res += a[i + 1].second;
            if (a[i].first + a[i].second == 10) res += a[i + 1].first;
        }
    }
    ans = max(ans, res);
    return res;
}

void simulateAnneal() {
    uniform_int_distribution<int> dis(1, n + bonus + 1);
   for (double T = 1e4; T > 1e-4; T *= 0.99) {
       int pos1 = dis(Rand), pos2 = dis(Rand), now = newRecord();
       swap(a[pos1], a[pos2]);
       if (!((a[n].first == 10) ^ bonus)) {
           int delta = newRecord() - now;
           if (exp(delta / T) < check(Rand)) swap(a[pos1], a[pos2]); //答案较劣，一定概率还原
		/*这里 delta 不带负号是因为答案劣才还原，而不是答案优还原*/
       } else swap(a[pos1], a[pos2]); // 答案不合法，直接还原
   }
}

int main(void) {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    auto start = getTime(); //起始时间
    cin >> n;
    for (int i = 1; i <= n; i++) cin >> a[i].first >> a[i].second;
    if (a[n].first == 10) { //判断是否有额外轮
        bonus = true;
        cin >> a[n + 1].first >> a[n + 1].second;
    }
    while (getTime() - start < 0.05) simulateAnneal();
    cout << ans;
    return 0;
}

[HAOI2006] 均分数据 - 模拟退火配合贪心

洛谷传送门
已知 $n$ 个正整数 $a_1,a_2 ... a_n$ 。今要将它们分成 $m$ 组，使得各组数据的数值和最平均，即各组的均方差最小。均方差公式如下：

σ = \sqrt{\frac{1}{n} \sum_{i = 1}^{n} (\bar{x} - x_{i})^{2}}, \bar{x} = \frac{1}{n} \sum_{i = 1}^{n} x_{i}

$\sigma = \sqrt{\frac 1n \sum\limits_{i=1}^n(\overline x - x_i)^2},\overline x = \frac 1n \sum\limits_{i=1}^n x_i$

其中 $\sigma$ 为均方差， $\bar{x}$ 为各组数据和的平均值， $x_i$ 为第 $i$ 组数据的数值和。

解析：暴力做法为枚举每个数据所在的位置，方案可以（？）用 $n$ 位 $m$ 进制储存。当然时间复杂度为 $O(m^n)$ ，跑得过就有鬼了ovo

由于更换某几个数据后答案差别不会太大，所以对应的“函数”具有一定的连续性，可以使用模拟退火。事实上，这是模拟退火可用的一个重要条件。

但是如果直接用模拟退火做的话答案区间过大，很难得到正确答案，考虑进行优化。使用贪心思想，从前往后找要选择的数据，每次把数据插入到数值和最小的数组中，可以证明这样做可以让均方差尽可能小。当然数据插入顺序要用模拟退火随机修改，从而提高贪心的正确率。

当然本题也可以不用模拟退火，而是 random_shuffle 随机生成数据顺序，不过正确率没有模拟退火高。

#include <iostream>
#include <cmath>
#include <chrono>
#include <iomanip>
#include <algorithm>
#include <random>
#include <cstring>

const int maxn = 50;

using namespace std;

int n, m;
int a[maxn], s[maxn];
double ans = 1e8;

random_device rd;
mt19937 Rand(rd());
uniform_real_distribution<double> check(0, 1);

inline auto getTime() {
	return chrono::system_clock::now().time_since_epoch().count() / 1e9;
}

double calcVariance() {
	memset(s, 0, sizeof(s));
	double sum = 0, res = 0;
	for (int i = 1; i <= n; i++) {
		int k = 1;
		for (int j = 2; j <= m; j++)
			if (s[j] < s[k]) k = j; // 贪心寻找数据和最小的数组
		s[k] += a[i];
	}
	for (int i = 1; i <= m; i++) sum += (double) s[i] / m;
	for (int i = 1; i <= m; i++)
		res += (s[i] - sum) * (s[i] - sum);
	res = sqrt(res / m);
	ans = min(ans, res);
	return res;
}

void simulateAnneal() {
	static uniform_int_distribution<int> dis(1, n);
	random_shuffle(a + 1, a + n + 1); //初始化一个排列（其实也可以不初始化）
	for (double T = 1e4; T >= 1e-4; T *= 0.99) {
		int x = dis(Rand), y = dis(Rand);
		double now = calcVariance(), delta;
		swap(a[x], a[y]);
		delta = calcVariance() - now;
		if (exp(-delta / T) < (double) check(Rand)) swap(a[x], a[y]);
	}
}

int main(void) {
	ios::sync_with_stdio(false);
	cin.tie(nullptr);
	auto start = getTime();
	cin >> n >> m;
	for (int i = 1; i <= n; i++) cin >> a[i];
	simulateAnneal(); //数据很水，所以跑一遍就够了
	cout << fixed << setprecision(2) << ans;
	return 0;
}

本文作者：Coel's Blog

本文链接：https://www.cnblogs.com/Coel-Flannette/p/16712447.html

posted @ 2022-09-20 20:25 秋泉こあい阅读(262) 评论(0) 编辑收藏举报

刷新页面返回顶部

登录后才能查看或发表评论，立即登录或者逛逛博客园首页

2025年3月

日

一

二

三

四

五

六

Coel's Blog

【Coel.学习笔记】随机化算法：模拟退火与爬山法

简介

原理

爬山法

模拟退火

例题讲解

[JSOI2008]球形空间产生器 - 爬山法浅谈

[POJ2420]A Star not a Tree? - 模拟退火入门

[JSOI/AHOI2014]保龄球 - 非数学的模拟退火

[HAOI2006] 均分数据 - 模拟退火配合贪心

搜索

常用链接

随笔分类

随笔档案

阅读排行榜

评论排行榜

推荐排行榜

最新评论