航迹规划——模拟退火算法

╰半夏微凉° 2022-05-29 11:05 225阅读 0赞

# 模拟退火算法简介 #

模拟退火是一种通用概率算法，可来在固定时间内寻求在一个大的搜寻空间内找到的最优解，也可以用来求解函数最优解。模拟退火是S. Kirkpatrick, C. D. Gelatt和M. P. Vecchi在1983年所发明。而V. Černý在1985年也独立发明此算法。  
  在这里不讲解模拟退火算法的基础知识，重点对模拟退火算法的两种应用程序做解析，需要了解基础知识的可以参考文章最后的参考资料。

--------------------

# 模拟退火算法流程图 #

模拟退火算法流程图

![这里写图片描述][SouthEast]

--------------------

# 模拟退火算法Matlab程序 #

## 求解函数最优解 ##

变量声明初始化

clear  
    clc  
    
    %生成初始解  
    sol_new2=1;%（1）解空间（初始解）  
    sol_new1=2-sol_new2^2;  
    sol_current1 = sol_new1;   
    sol_best1 = sol_new1;  
    sol_current2 = sol_new2;   
    sol_best2 = sol_new2;  
    E_current = inf;  
    E_best = inf;

退火算法参数设置

rand('state',sum(clock)); %初始化随机数发生器  
    t=90; %初始温度  
    tf=89.9; %结束温度  
    a = 0.99; %温度下降比例

主循环

while t>=tf%（7）结束条件 
        for r=1:100 %退火次数 
    
            %产生随机扰动（3）新解的产生 
            sol_new2=sol_new2+rand*0.2;  
            sol_new1=2-sol_new2^2;  
    
            %检查是否满足约束 
            if sol_new1^2-sol_new2>=0 && -sol_new1-sol_new2^2+2==0 && sol_new1>=0 &&sol_new2>=0  
            else  
                sol_new2=rand*2;  
                sol_new1=2-sol_new2^2;  
                continue;  
            end  
    
            %退火过程 
            E_new=sol_new1^2+sol_new2^2+8;%（2）目标函数 
            if E_new<E_current%（5）接受准则 
                    E_current=E_new;  
                    sol_current1=sol_new1;  
                    sol_current2=sol_new2;  
                    if E_new<E_best  
                        %把冷却过程中最好的解保存下来 
                        E_best=E_new;  
                        sol_best1=sol_new1;  
                        sol_best2=sol_new2;  
                    end  
            else  
                    if rand<exp(-(E_new-E_current)/t)%（4）代价函数差 
                        E_current=E_new;  
                        sol_current1=sol_new1;  
                        sol_current2=sol_new2;  
                    else  
                        sol_new1=sol_current1;  
                        sol_new2=sol_current2;  
                    end  
            end  
            plot(r,E_best,'.')  
            hold on  
        end  
        t=t*a;%（6）降温 
    end

显示结果

disp('最优解为：')  
    disp(sol_best1)  
    disp(sol_best2)  
    disp('目标表达式的最小值等于：')  
    disp(E_best)

### 运行结果如图 ###

![这里写图片描述][SouthEast 1]

通过运行结果可以直观的看出最后程序的结果趋近于最优解，在这里最多退火100次，多次迭代后结果会更好。

## 求解二维空间最优路径 ##

程序初始化，设置50个需要访问的点

clc, clear
    % 设置处理点的数量
    N = 50;
    M = zeros(N, 2);
    % 产生100个位置不同的点
    M(:, 1) = randperm(N);
    M(:, 2) = randperm(N);

计算两点之间的距离

% 初始化距离矩阵
    d = zeros(N);
    % 统计每两个点之间的距离
    for i=1:N-1
        for j=i+1:N
            d(i, j) = sqrt((M(i,1)-M(j,1))*(M(i,1)-M(j,1)) + (M(i,2)-M(j,2))*(M(i,2)-M(j,2)));
        end
    end

for循环生成的两点之间的距离思路为

1.  先计算（1，2），（1，3），（1，4）……（1，N）的距离存在二维数组d中。
2.  计算（2，3），（2，4）（2，5）……（2，N）的距离存在二维数组d中。
3.  计算（3，4）（3，5），（3，6）……（3，N）的距离存在二维数组d中。
4.  ……
5.  计算（N-1，N）的距离存在二维数组d中，结束。

此时d中存的数据为上三角矩阵，因为（1，2）和（2，1）的距离是一样的，因此将d中的数据扩充。

d = d + d';

通过打乱50个坐标点的链接顺序，多次计算其总路线长度，得出近似最优解。目的是为了和退火算法比较。

% 初始化路径及路线长度
    path = zeros(1, N);
    length = inf;
    % 先确定一个比较好的初始解
    for i=1:1000
        path0 = randperm(N); % randperm 1-N 乱序排列
        temp = 0;
        for j=1:N-1
            temp = temp + d(path0(j) , path0(j+1));
        end
        if temp < length
           length = temp;
           path = path0;
        end
    end
    length
    subplot(1, 2, 1)
    plot(M(path, 1), M(path, 2), '-o');
    title("Initial Solution");

退火算法求解50个坐标点的最优路径  
初始化参数

% 设置模拟退火的参数
    e = 0.1^30; % 结束温度
    L = 20000; % 退火次数
    at = 0.999; % 温度下降比例
    T = 1; % 初始温度

主循环

for k = 1:L
        % 产生新解
        c = floor(rand(1, 2)*(N-2)) + 2; % floor 比它小的最大整数 
        c = sort(c);
        c1 = c(1);
        c2 = c(2);
        change = d(path(c1-1), path(c2)) + d(path(c1), path(c2+1)) - d(path(c1-1), path(c1)) - d(path(c2), path(c2+1));
        if change < 0
            path = [path(1:c1-1), path(c2:-1:c1), path(c2+1: N)];
            length = length + change;
        elseif exp(-change/T) >= rand % 接受条件
            path = [path(1:c1-1), path(c2:-1:c1), path(c2+1: N)];
            length = length + change;
        end
        T = T * at;
        if T < e
            break;
        end
    end

输出结果对比

length
    subplot(1, 2, 2)
    plot(M(path, 1), M(path, 2), '-*');
    title("Simulte Anneal Solution");
    hold off

产生新解的思路

1.  随机生成两个数，A，B。
2.  将A，B之间的序列调换。  
    ![这里写图片描述][SouthEast 2]
3.  计算变化量。
4.  总路长比之前解短的话保留结果，总路长比之前解长的话按接受条件`exp(-change/T) >= rand`接收。
5.  满足终止条件后推出循环。

### 运行结果 ###

![这里写图片描述][SouthEast 3]

从上图可以很明显的看出模拟退火算法比初始的最优解效果好，通过最后输出的总路径长度也可以看出模拟退火算法比初始的最优解效果好N倍~

--------------------

# Marlab 完整代码 #

两个代码在一起放着哦，大家注意看下哦  
[模拟退火算法Matlab程序二合一下载链接][Matlab]

--------------------

# 参考资料 #

[优化算法——模拟退火算法 ][Link 1]

[SouthEast]: /images/20220529/79dde8498b804d6184afa143a3eca60d.png
[SouthEast 1]: /images/20220529/d1d7396c4cad455b804febce790fec65.png
[SouthEast 2]: /images/20220529/d3056d8430124eafadc20a6beddb34f5.png
[SouthEast 3]: /images/20220529/86bf427cc0ae48a8941b18a278a05353.png
[Matlab]: https://download.csdn.net/download/qq_16775293/10302457
[Link 1]: http://blog.csdn.net/google19890102/article/details/45395257