首页 > 编程语言 >玄学乱搞算法——模拟退火，SA

玄学乱搞算法——模拟退火，SA

时间：2024-08-27 21:26:27浏览次数：13  
标签：rand space int double 乱搞 模拟退火 SA 最优

\(\texttt{0x00:}\) 前言

在此之前只对模拟退火的大名有所耳闻，但并未在我的认知上激起太大的风浪，直到……

在外培的一场模拟赛上，队内大佬 yyc 在丝毫没有思路的情况下用 SA 骗了 70pts，赛后使得给我们上课的清华姚班老师惊掉下巴。

至此，在感叹 SA 的神力的同时，它也进入了我的学习计划中。

\(\texttt{0x01:}\) 介绍

退火原本是一个物理概念，是指一种热处理工艺，将固体材料加热到一定温度后，让其缓慢冷却的过程。

模拟退火算法（Simulate Anneal，SA）是一种通用概率演算法，用来在一个大的搜寻空间内找寻命题的最优解。模拟退火是由 S.Kirkpatrick，C.D.Gelatt 和 M.P.Vecchi 在1983年所发明的。V.Černý 在 1985 年也独立发明此演算法。模拟退火算法是解决 TSP 问题的有效方法之一。——百度百科

模拟退火是一种随机化算法。当一个问题的方案数量极大（甚至是无穷的）而且不是一个单峰函数时，我们常使用模拟退火求解。

简单来讲，模拟退火是一种启发式搜索，它通过模拟物理降温的过程，先随机取一个起始点，然后随机地搜索其他点，若更优就跳过去，否则以一定的概率跳过去。由于模拟退火的过程中有概率接受较劣的解，所以它有一定概率找到全局最优解而不是受困于某一局部最优解。

\(\texttt{0x02:}\) 算法过程

首先设定几个参数：初温 \(T_0\)，末温 \(T_E\)，降温系数 \(k\)。

在算法开始时钦定温度 \(T = T_0\)，然后进行最优解的寻找，重复找最优解的过程并不断将温度乘上 \(k\) 直到 \(T = T_E\)。

具体来讲，在找最优解时的准则：如果新状态的解更优则修改答案，否则以一定概率接受新状态。

我们定义当前温度为 \(T\)，新状态 \(S'\) 与已知状态 \(S\)（新状态由已知状态通过随机的方式得到）之间的能量（值）差为 \(\Delta_E\)（\(\Delta_E\geqslant 0\)），则发生状态转移（修改最优解）的概率为：

\[\large P(\Delta_E) = \begin{cases} 1, \space\space\space\space\space\space\space\space\space\space\space\space\text{S' is better than S,}\\ \frac{-\Delta_E}{T},\space\space\space\space\space\space \text{otherwise.} \end{cases}\]

注意为了使得解更为精确，我们通常不直接取当前解作为答案，而是在退火过程中维护遇到的所有解的最优值。

注意：模拟退火只能用来求解具有一定连续性的函数的最优解，这个连续性可以具体看题目中一个微小的策略改变是否导致了函数值的剧变，若导致了，则没多少连续性，否则可以认为有一定连续性。

模拟退火的过程如图所示：

\(\texttt{Code:}\)

void SA() {
    double T = 1e3;
    while(T > 1e-7) {
        //do something...
        ll now = calc(); //求当前的函数值
        ll dt = now - ans; //ans 是全局最优解
        if(dt < 0) ans = now; //更优则直接接受
        else if(exp(-(double)dt / T) * RAND_MAX < (double)rand()) 
			//不是更优则以一定概率接受
        T *= k;
    }
}

\(\texttt{0x03:}\) 一些例题

P1337 [JSOI2004] 平衡点 / 吊打XXX

题目大意

给定 \(n\) 的点表示每根绳子的固定点，每条绳子都系着一个重物，求出平衡点。

思路

这道题正解应该是三分套三分，但是不会，考虑乱搞。

根据基本的物理知识，平衡点会偏向重物所在的地方，事实上，假如钦定一个点 \(p\)，定义每条绳子的势能为点 \(p\) 到固定点的距离乘上所挂重物的重量，那么平衡点的势能应该是最小的，这个可以通过移动平衡点来证明确实很有道理。

同时在点 \(p\) 移动的距离很小的情况下，势能的变化也较小，所以函数具有连续性，可以用模拟退火来做。

这里我们顺便用当前温度代表当前的步长（随机的范围），便能轻松写出如下代码：

#include <cmath>
#include <ctime>
#include <iostream>

#define x first
#define y second

using namespace std;

const int N = 1010;
const double eps = 1e-9;
typedef pair<double, double> PDD;
int n;
struct node{
    int xx, yy, w;
}q[N];
double ans = 1e18;
int sumx, sumy;
PDD anss;

inline double rand(double l, double r) {
    return (double)rand() / RAND_MAX * (r - l) + l;
} 

inline double get_dist(PDD a, PDD b) {
    return sqrt((a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y));
}

inline double calc(PDD now) {
    double res = 0;
    for(int i = 1; i <= n; i++)
        res += get_dist(now, {q[i].xx, q[i].yy}) * q[i].w;
    return res;
}

void SA() {
    PDD now(rand(-10000.0, 10000.0), rand(-10000.0, 10000.0));
    double T = 1e4;
    while(T > 1e-5) {
        PDD newp(now.x + (double)((rand() << 1) - RAND_MAX) * T, now.y + (double)((rand() << 1) - RAND_MAX) * T);
        double fut = calc(newp);
        double dt = fut - ans;
        if(dt < 0) ans = fut, anss = newp, now = newp;
        else if(exp(-dt / T) * RAND_MAX > (double)rand()) now = newp;
        T *= 0.95;
    }
}

int main() {
    srand((unsigned)time(0));
    scanf("%d", &n);
    for(int i = 1; i <= n; i++) {
        scanf("%d%d%d", &q[i].xx, &q[i].yy, &q[i].w);
    }
    SA();
    printf("%.3lf %.3lf", anss.x, anss.y);
    return 0;
}

虽然函数都写对了，但是这样连样例都过不了，怎么办？

其实有模拟退火的三大调整策略：

  1. 改进随机数的生成或是用贪心等策略优化找最优解的过程；

  2. 在时间允许的条件下，可以适当调高降温系数 \(k\)，调整初温或末温，还可以多跑几遍，甚至卡时；（调参是模拟退火最好玩的一环，不爽不要玩！）

  3. 采用精度更高的数据类型，比如 long double。（其实这条没什么用，主要看前两条）

由于温度和步长有关，所以这道题可以将末温调低一点，同时时间很充裕，所以可以调大降温系数和多跑几遍。

这样我们就 AC 了这道题。

其实还可以将起点设为固定点中的平均值，这样更容易跑出最优解。

P2503 [HAOI2006] 均分数据

题目大意

给定一个长为 \(n\) 的序列 \(\{a_n\}\)，先要将其分成 \(m\) 段，使得各组数字之和的标准差最小。

思路

因为数据范围比较小，所以其实各种乱搞都能过（

考虑模拟退火。

比较朴素的一种想法是，每次找最优解时就随机将这 \(n\) 个数分到 \(m\) 个组里去，然后再计算标准差，更新最优解。

但是这样不容易搜到最优解，因为可能将许多大数丢到同一个集合去，所以很多情况都没有用。

综上，我们可以加一个贪心优化，具体地，先随机打乱原序列 \(a\)，然后每次将一个数加到总和最小的那个组里去，这样就能更快搜到最优解。

\(\texttt{Code:}\)

#include <cmath>
#include <ctime>
#include <cstring>
#include <iostream>
#include <algorithm>

using namespace std;

const int N = 25, M = 10;

int n, m;
int a[N], s[M];
double ans = 1e9;

double calc() {
    memset(s, 0, sizeof s);
    for(int i = 1; i <= n; i++) {
        int pos = 1;
        for(int j = 2; j <= m; j++)
            if(s[j] < s[pos]) pos = j;
        s[pos] += a[i];
    }

    double avg = 0;
    for(int i = 1; i <= m; i++) avg += 1.0 * s[i] / m;
    double res = 0;
    for(int i = 1;  i <= m; i++) res += (s[i] - avg) * (s[i] - avg);
    res = sqrt(res / m);
    ans = min(ans, res);
    return res;
}

void SA() {
    random_shuffle(a + 1, a + n + 1);
    double T = 1e4;
    while(T > 1e-6) {
        int x = rand() % n + 1, y = rand() % n + 1;
        if(x == y) continue;
        // printf("%d %d\n", x, y);
        double now = calc();
        swap(a[x], a[y]);
        double fut = calc();
        double dt = fut - now;
        if(exp(-dt / T) * RAND_MAX < (double)rand())
            swap(a[x], a[y]);
        T *= 0.99;
    }
}

int main() {
    scanf("%d%d", &n, &m);
    for(int i = 1; i <= n; i++) scanf("%d", &a[i]);
    for(int i = 0; i < 10; i++) SA();
    printf("%.2lf", ans);
    return 0;
}

P4044 [AHOI2014/JSOI2014] 保龄球

题目大意

一场保龄球比赛共有 \(N\) 轮，每一轮都会有 \(10\) 个木瓶放置在木板道的另一端。

每一轮中，选手都有两次投球的机会来尝试击倒这 \(10\) 个木瓶。

对于每一次投球机会，选手投球的得分等于这一次投球所击倒的木瓶数量。

选手每一轮的得分是他两次机会击倒全部木瓶的数量。

对于每一个轮次，有如下三种情况：

“全中”：若选手第一次尝试就击倒了全部 \(10\) 个木瓶，则“全中”。在一个“全中”轮中，不需要再进行第二次投球尝试。作为奖励，选手的下一轮的得分将会被乘 \(2\) 计入总分。

“补中”：若选手两次尝试击倒了 \(10\) 个木瓶，则“补中”。作为奖励，选手的下一轮中的第一次尝试的得分将会被乘以 \(2\) 计入总分。

“失误”：若选手未能通过两次尝试击倒全部的木瓶，则“失误”。失误没有奖励。

此外，如果第 \(N\) 轮是“全中”，那么选手将一共进行 \(N+1\) 轮比赛。

附加轮的规则只执行一次。

最后，选手的总得分就是附加轮规则执行过，并且分数按上述规则加倍后的每一轮分数之和。

JYY 刚刚进行了一场 \(N\) 个轮次的保龄球比赛，但是，JYY 非常不满意他的得分。

JYY 希望通过更改每一轮成绩的顺序来最大化他的得分。

当然了，JYY 不希望做的太假，他希望保证重新排列之后，所需要进行的轮数和重排前所进行的轮数是一致的：

比如如果重排前 JYY 在第 \(N\) 轮打出了“全中”，那么重排之后，第 \(N\) 轮还得是“全中”以保证比赛一共进行 \(N+1\) 轮；同样的，如果 JYY 第 \(N\) 轮没有打出“全中”，那么重排过后第 \(N\) 轮也不能是全中。

求 JYY 可以得到的最高的分数。

（题目太长了，已经尽力简化了）

思路

这道题的数据范围同样较小，所以可以乱搞（

考虑模拟退火。

先看一下有没有 \(n + 1\) 轮，然后每次随机生成一个排列，退火过程中随机交换此排列中的两个位置，再按规则计算分数，更新最大值即可。

\(\texttt{Code:}\)

#include <cmath>
#include <ctime>
#include <iostream>

#define x first
#define y second

using namespace std;

const int N = 55;
typedef pair<int, int> PII;
int n, nn;
PII q[N];
int ans;

inline void swap(PII &x, PII &y) {
    PII t = y;
    y = x;
    x = t;
}

inline int calc() {
    int res = 0;
    for(int i = 1; i <= nn; i++) {
        res += q[i].x + q[i].y;
        if(i <= n) {
            if(q[i].x == 10) res += q[i + 1].x + q[i + 1].y;
            else if(q[i].x + q[i].y == 10) res += q[i + 1].x;
        }
    }
    ans = max(ans, res);
    return res;
}

void SA() {
    double T = 1e6;
    while(T > 1e-6) {
        int x = rand() % nn + 1, y = rand() % nn + 1;
        int now = calc();
        swap(q[x], q[y]);
        if(n + (q[n].x == 10) == nn) {
            int fut = calc();
            int dt = fut - now;
            if(exp(dt / T) < (double)rand() / RAND_MAX)
                swap(q[x], q[y]);
        }
        else swap(q[x], q[y]);
        T *= 0.996;
    }
}

int main() {
    srand((unsigned)(time(0)));
    srand((unsigned)rand() * rand()); //玄学种子
    scanf("%d", &n);
    for(int i = 1; i <= n; i++) scanf("%d%d", &q[i].x, &q[i].y);
    if(q[n].x == 10) nn = n + 1, scanf("%d%d", &q[nn].x, &q[nn].y);
    else nn = n;
    for(int i = 1; i < 3; i++) SA();
    printf("%d\n", ans);
    return 0;
}

标签：rand,space,int,double,乱搞,模拟退火,SA,最优
From： https://www.cnblogs.com/Brilliant11001/p/18383586

相关文章

  • NSSCTF [NISACTF 2022]babyserialize
    <?phpinclude"waf.php";classNISA{public$fun="show_me_flag";public$txw4ever;//1shellpublicfunction__wakeup(){if($this->fun=="show_me_flag"){hint();}}......
  • Vue 3 + GSAP：解锁高性能动画的神奇组合
    在一个偶然的场景使用了GSAP（GreenSockAnimationPlatform），感觉挺好玩的，在此浅浅记录一下。GSAP是一个功能强大的JS动画库，常用于创建高性能、流畅的动画。它不仅支持基本的动画，还提供了时间轴控制、缓动效果等高级功能。将其与Vue3结合，可以轻松给Vue组件添加复杂的动......
  • 如何发明 SAM
    如何发明SAM我们想做一个结构，接受全部的子串，我们发现，如果考虑增量构造，每次加1个字符，增多的字串就是原来的后缀加上这个。那么我们就这样做。所以我们一直需要一个集合，这些点接受全部后缀。（以上由红线组成的到根的路径就是我们每个时刻维护的集合，即所有后缀）那么当在abc上......
  • C# yield keyword relieve congest and consume at the same time with produce
    usingSystem.Threading;namespaceConsoleApp57{internalclassProgram{staticvoidMain(string[]args){PrintNumers();Console.WriteLine("Hello,World!");}staticvoidPrintN......
  • esp-idf vscode debug command 'espIdf.getXtensaGdb' not found
    esp32idfvscodedebug错误vscode中配置文件采用的是正点原子的，调用gdb的时候，提示报错，找不到相应的命令launch.json文件中gdb的配置如下{"version":"0.2.0","configurations":[ { "name":"GDB", "type":"cppdbg", &......
  • ReconSage一个快速的在线子域名扫描工具
    ReconSage一个快速的在线子域名扫描工具https://www.reconsage.com/tools/scan-subdomain如下示例：子域名扫描工具的主要作用是帮助用户发现和枚举目标域名下的所有子域名。这些工具通过不同的技术手段，如字典枚举、API调用、搜索引擎查询等，来识别和展示目标域名的子域名信息，包括子......
  • 【图像分割】复合粒子群算法CPSOGSA图像多级阈值分割【含Matlab源码 7349期】
    ✅博主简介：热爱科研的Matlab仿真开发者，修心和技术同步精进，Matlab项目合作可私信或扫描文章底部QQ二维码。......
  • rasadhlp.dll揭秘：远程访问服务作用与修复丢失的实用手册
    rasadhlp.dll是一个与Windows操作系统相关的动态链接库（DLL）文件，通常与远程访问服务（RemoteAccessService,RAS）的功能实现有关。这个DLL文件可能包含了处理RAS连接管理、资源管理和与其他系统组件交互等功能所需的函数和资源，对于确保远程访问服务的正常运行非常重要。当rasad......
  • CISAW信息安全保障人员考试合格分数是多少？
    在网络安全领域，众多认证证书中，CISAW信息安全保障人员认证备受瞩目。随着信息安全专家的需求激增，国内面临着专业人才的短缺，越来越多人投身于CISAW认证考试。那么，要通过CISAW认证需多少分数呢？答案是84分（含）以上。如果考生作弊，其笔试成绩将被记为0分。考试总时长为150分钟，满分......
  • CF1554C Mikasa
    本蒟蒻的第一篇题解，如果有什么不足，请dalao轻喷。。。题意题目传送门给出两个数\(n\)和\(m\)，求出最小的整数\(p\)，使得\(p\)\(\notin\){\(n\oplus1,n\oplus2,n\oplus3,···,n\oplusm\)}。思路异或先思考这道题之前，还不如先思考这道题目所要用的符号......

赞助商

阅读排行

网站主页关于我们联系我们网站地图

本网站内容转载自其他媒体，侵权联系[admin##ips99.com]。

Copyright © 2020-2023 IPS99 版权所有 IPS99