模拟退火算法

云上西

于 2024-08-01 12:45:00 发布

阅读量986

点赞数 22

分类专栏：清风建模数学建模文章标签：模拟退火算法算法机器学习回归开发语言数学建模 matlab

本文链接：https://blog.csdn.net/huahua121_/article/details/140834805

版权

数学建模同时被 2 个专栏收录

40 篇文章 4 订阅

订阅专栏

清风建模

28 篇文章 0 订阅

订阅专栏

该博客为个人学习清风建模的学习笔记，代码全部摘自清风老师，部分课程可以在B站：【强烈推荐】清风：数学建模算法、编程和写作培训的视频课程以及Matlab等软件教学_哔哩哔哩_bilibili

这个算法为公开课，地址在：数学建模清风第二次直播：模拟退火算法_哔哩哔哩_bilibili

这个笔记结合了清风老师的课堂讲义和课件。

1爬山法

1.1爬山法找最大值

怎么找到这个一元函数的最大值？（只有一个上下界约束，即函数的定义域）

1.2爬山法存在的问题

假如我们现在站在上帝视角来看这个图，很明显我们要找的最大值在左侧的蓝色点处，但是我们刚刚找到的是位于黄点处的局部最大值点。

为什么爬山法的搜索在黄色点处就停止了呢？是因为我们观察到该点左右两侧邻域的目标函数值都小于该点的目标函数值。若我们想要搜索到最大值处，就还需要往左走一段，因此爬山法实际上是一种“眼光狭隘”的搜索算法( 贪心 ) 。

模拟退火算法就能克服这个问题。它在到达黄色点处后，会以一定的概率接受一个比当前解要差的解（即黄色点处左右两侧的邻域），此时就有几率跳出黄色点，从而我们就有了找到蓝色点的可能性。

2模拟退火的介绍

2.1从物理角度看模拟退火

2.2从优化角度看模拟退火

3模拟退火的有关概念

3.1Metropolis准则

3.2冷却进度表

如何设置模拟退火的温度：

3.3新解的产生规则

4模拟退火求最值

4.1求解一元函数最大值

%% SA 模拟退火: 求解函数y = 11*sin(x) + 7*cos(5*x)在[-3,3]内的最大值(动画演示）
tic
clear; clc

%% 绘制函数的图形
x = -3:0.1:3;
y = 11*sin(x) + 7*cos(5*x);
figure
plot(x,y,'b-')
hold on  % 不关闭图形，继续在上面画图

%% 参数初始化
narvs = 1; % 变量个数
T0 = 100;   % 初始温度
T = T0; % 迭代中温度会发生改变，第一次迭代时温度就是T0
maxgen = 200;  % 最大迭代次数
Lk = 100;  % 每个温度下的迭代次数
alfa = 0.95;  % 温度衰减系数
x_lb = -3; % x的下界
x_ub = 3; % x的上界

%%  随机生成一个初始解
x0 = zeros(1,narvs);
for i = 1: narvs
    x0(i) = x_lb(i) + (x_ub(i)-x_lb(i))*rand(1);    
end
y0 = Obj_fun1(x0); % 计算当前解的函数值
h = scatter(x0,y0,'*r');  % scatter是绘制二维散点图的函数（这里返回h是为了得到图形的句柄，未来我们对其位置进行更新）

%% 定义一些保存中间过程的量，方便输出结果和画图
max_y = y0;     % 初始化找到的最佳的解对应的函数值为y0
MAXY = zeros(maxgen,1); % 记录每一次外层循环结束后找到的max_y (方便画图）

%% 模拟退火过程
for iter = 1 : maxgen  % 外循环, 我这里采用的是指定最大迭代次数
    for i = 1 : Lk  % 内循环，在每个温度下开始迭代
        y = randn(1,narvs);  % 生成1行narvs列的N(0,1)随机数
        z = y / sqrt(sum(y.^2)); % 根据新解的产生规则计算z
        x_new = x0 + z*T; % 根据新解的产生规则计算x_new的值
        % 如果这个新解的位置超出了定义域，就对其进行调整
        for j = 1: narvs
            if x_new(j) < x_lb(j)
                r = rand(1);
                x_new(j) = r*x_lb(j)+(1-r)*x0(j);
            elseif x_new(j) > x_ub(j)
                r = rand(1);
                x_new(j) = r*x_ub(j)+(1-r)*x0(j);
            end
        end
        x1 = x_new;    % 将调整后的x_new赋值给新解x1
        y1 = Obj_fun1(x1);  % 计算新解的函数值
        if y1 > y0    % 如果新解函数值大于当前解的函数值
            x0 = x1; % 更新当前解为新解
            y0 = y1;
        else
            p = exp(-(y0 - y1)/T); % 根据Metropolis准则计算一个概率
            if rand(1) < p   % 生成一个随机数和这个概率比较，如果该随机数小于这个概率
                x0 = x1; % 更新当前解为新解
                y0 = y1;
            end
        end
        % 判断是否要更新找到的最佳的解
        if y0 > max_y  % 如果当前解更好，则对其进行更新
            max_y = y0;  % 更新最大的y
            best_x = x0;  % 更新找到的最好的x
        end
    end
    MAXY(iter) = max_y; % 保存本轮外循环结束后找到的最大的y
    T = alfa*T;   % 温度下降
    pause(0.01)  % 暂停一段时间(单位：秒)后再接着画图
    h.XData = x0;  % 更新散点图句柄的x轴的数据（此时解的位置在图上发生了变化）
    h.YData = Obj_fun1(x0); % 更新散点图句柄的y轴的数据（此时解的位置在图上发生了变化）
end

disp('最佳的位置是：'); disp(best_x)
disp('此时最优值是：'); disp(max_y)

pause(0.5)
h.XData = [];  h.YData = [];  % 将原来的散点删除
scatter(best_x,max_y,'*r');  % 在最大值处重新标上散点
title(['模拟退火找到的最大值为', num2str(max_y)])  % 加上图的标题

%% 画出每次迭代后找到的最大y的图形
figure
plot(1:maxgen,MAXY,'b-');
xlabel('迭代次数');
ylabel('y的值');
toc

4.2求解二元函数最小值

%% SA 模拟退火: 求解函数y = x1^2+x2^2-x1*x2-10*x1-4*x2+60在[-15,15]内的最小值(动画演示）
tic

clear; clc
%% 绘制函数的图形
figure 
x1 = -15:1:15;
x2 = -15:1:15;
[x1,x2] = meshgrid(x1,x2);
y = x1.^2 + x2.^2 - x1.*x2 - 10*x1 - 4*x2 + 60;
mesh(x1,x2,y)
xlabel('x1');  ylabel('x2');  zlabel('y');  % 加上坐标轴的标签
axis vis3d % 冻结屏幕高宽比，使得一个三维对象的旋转不会改变坐标轴的刻度显示
hold on  % 不关闭图形，继续在上面画图

%% 参数初始化
narvs = 2; % 变量个数
T0 = 100;   % 初始温度
T = T0; % 迭代中温度会发生改变，第一次迭代时温度就是T0
maxgen = 200;  % 最大迭代次数
Lk = 100;  % 每个温度下的迭代次数
alfa = 0.95;  % 温度衰减系数
x_lb = [-15 -15]; % x的下界
x_ub = [15 15]; % x的上界

%%  随机生成一个初始解
x0 = zeros(1,narvs);
for i = 1: narvs
    x0(i) = x_lb(i) + (x_ub(i)-x_lb(i))*rand(1);    
end
y0 = Obj_fun2(x0); % 计算当前解的函数值
h = scatter3(x0(1),x0(2),y0,'*r');  % scatter3是绘制三维散点图的函数（这里返回h是为了得到图形的句柄，未来我们对其位置进行更新）

%% 定义一些保存中间过程的量，方便输出结果和画图
min_y = y0;     % 初始化找到的最佳的解对应的函数值为y0
MINY = zeros(maxgen,1); % 记录每一次外层循环结束后找到的min_y (方便画图）

%% 模拟退火过程
for iter = 1 : maxgen  % 外循环, 我这里采用的是指定最大迭代次数
    for i = 1 : Lk  % 内循环，在每个温度下开始迭代
        y = randn(1,narvs);  % 生成1行narvs列的N(0,1)随机数
        z = y / sqrt(sum(y.^2)); % 根据新解的产生规则计算z
        x_new = x0 + z*T; % 根据新解的产生规则计算x_new的值
        % 如果这个新解的位置超出了定义域，就对其进行调整
        for j = 1: narvs
            if x_new(j) < x_lb(j)
                r = rand(1);
                x_new(j) = r*x_lb(j)+(1-r)*x0(j);
            elseif x_new(j) > x_ub(j)
                r = rand(1);
                x_new(j) = r*x_ub(j)+(1-r)*x0(j);
            end
        end
        x1 = x_new;    % 将调整后的x_new赋值给新解x1
        y1 = Obj_fun2(x1);  % 计算新解的函数值
        if y1 < y0    % 如果新解函数值小于当前解的函数值
            x0 = x1; % 更新当前解为新解
            y0 = y1;
        else
            p = exp(-(y1 - y0)/T); % 根据Metropolis准则计算一个概率
            if rand(1) < p   % 生成一个随机数和这个概率比较，如果该随机数小于这个概率
                x0 = x1;  % 更新当前解为新解
                y0 = y1;
            end
        end
        % 判断是否要更新找到的最佳的解
        if y0 < min_y  % 如果当前解更好，则对其进行更新
            min_y = y0;  % 更新最小的y
            best_x = x0;  % 更新找到的最好的x
        end
    end
    MINY(iter) = min_y; % 保存本轮外循环结束后找到的最小的y
    T = alfa*T;   % 温度下降
    pause(0.02)  % 暂停一段时间后(单位：秒)再接着画图
    h.XData = x0(1);  % 更新散点图句柄的x轴的数据（此时解的位置在图上发生了变化）
    h.YData = x0(2); % 更新散点图句柄的y轴的数据（此时解的位置在图上发生了变化）
    h.ZData = Obj_fun2(x0);  % 更新散点图句柄的z轴的数据（此时粒子的位置在图上发生了变化）
end

disp('最佳的位置是：'); disp(best_x)
disp('此时最优值是：'); disp(min_y)

pause(0.5) 
h.XData = [];  h.YData = [];  h.ZData = [];  % 将原来的散点删除
scatter3(best_x(1), best_x(2), min_y,'*r');  % 在最小值处重新标上散点
title(['模拟退火找到的最小值为', num2str(min_y)])  % 加上图的标题

%% 画出每次迭代后找到的最小y的图形
figure
plot(1:maxgen,MINY,'b-');
xlabel('迭代次数');
ylabel('y的值');
toc

4.3Matlab自带的模拟退火函数

4.3.1打开工具箱

4.3.2设置目标函数

4.3.3停止准则

4.3.4退火设置

4.3.5接受准则

4.3.6后续使用其他函数求解

4.3.7模拟退火过程可视化

4.3.8展示迭代过程

4.3.9计算结果并生成函数代码

5旅行商问题

6书店买书问题

7总结

（ 1 ）模拟退火可以求解带有复杂约束条件的规划问题吗？

可以的，和粒子群算法一样，我们可以在目标函数上引入惩罚项来进行求解，也可以在产生新解后判断新解是否满足约束条件，若不满足则直接拒绝。

（ 2 ）模拟退火的衍生算法有哪些？

加温退火法、有记忆的模拟退火、带返回搜索的模拟退火、多次寻优法、回火退火法、并行模拟退火算法等，关于这些算法的具体介绍可参考康立山等1994 年编写的《非数值并行算法（第一册）模拟退火算法》，如果你不是专门研究模拟退火理论的，这些方法的细节我们可以不需要知道。

（ 3 ）模拟退火现在的研究还多吗，主要是研究什么？

模拟退火的研究可分为两类，一类基于马尔可夫链的有关理论，探究不同条件下模拟退火求解的收敛性（即纯理论研究）；另一类针对具体问题，给出了模拟退火的若干成功应用（即应用层面）。目前模拟退火相关的论文数量较十多年前有很大的下滑，知网搜索可以发现，智能算法中使用最多的是我们下一节要介绍的遗传算法。

（ 4 ）站在数学建模比赛的角度，模拟退火有什么用？

（ 1 ）求解无约束 ( 或只有上下界约束 ) 函数最值；

（ 2 ）求解简单的组合优化问题，例如旅行商问题和书店买书问题；

（ 3 ）未来再遇到非线性整数规划问题时，就可以想想能否用模拟退火求解。

注：第 (1) 点粒子群算法也可以求解哦，第 (2) 点和第 (3) 点我们下一节要介绍的遗传算法也可以求解，另外对于常规的规划问题，大家一定要优先考虑我们更新12 里面介绍的各种求解函数，智能算法只有在那些函数都无能为力的时候再想着使用。

云上西

关注

22
点赞
踩
24

收藏

觉得还不错? 一键收藏
0
评论
模拟退火算法

2.2从优化角度看模拟退火1）模拟退火可以求解带有复杂约束条件的规划问题吗？可以的，和粒子群算法一样，我们可以在目标函数上引入惩罚项来进行求解，也可以在产生新解后判断新解是否满足约束条件，若不满足则直接拒绝。2）模拟退火的衍生算法有哪些？加温退火法、有记忆的模拟退火、带返回搜索的模拟退火、多次寻优法、回火退火法、并行模拟退火算法等，关于这些算法的具体介绍可参考康立山等1994年编写的《非数值并行算法（第一册）模拟退火算法》，如果你不是专门研究模拟退火理论的，这些方法的细节我们可以不需要知道。3。
复制链接

扫一扫