【路径规划】基于改进差分算法实现三维多无人机协同航迹规划

一、差分进化算法的介绍

   差分进化算法(Differential Evolution, DE)是一种基于群体差异的启发式随机搜索算法，该算法是由R.Storn和K.Price为求解Chebyshev多项式而提出的。DE算法也属于智能优化算法，与前面的启发式算法，如ABC，PSO等类似，都属于启发式的优化算法。DE算法是我在一篇求解盒子覆盖问题论文中使用的一种优化算法。

二、差分进化算法的流程

1. 初始化种群
2. 变异
3. 交叉
4. 选择

(DE流程)

三、差分进化的具体步骤

   对于无约束优化问题

利用差分进化求解这样的优化问题，主要分为初始化、变异、交叉和选择等几项操作。

1、初始化

   如前面的的群智能优化算法一样，差分进化也需要初始化种群：

其中，是第个个体，表示第维。

其中，和分别为第维的下界和上界，表示在区间上的随机数。

2、变异

       DE算法通过差分策略实现个体变异，常见的差分策略是随机选取种群中两个不同的个体，将其向量差缩放后与待变异个体进行向量合成。

其中，，和是三个随机数，区间为，称为缩放因子，为一个确定的常数。表示第代。

3、交叉

   交叉操作的目的是随机选择个体，因为差分进化也是一种随机算法，交叉操作的方法是：

其中，称为交叉概率。通过概率的方式随机生成新的个体。

4、选择

   在DE中采用的是贪婪选择的策略，即选择较优的个体作为新的个体。

登录后复制

 
 
%清空环境
clear;
 
%定义地形对象
cTerrain = CTerrain3D;
%定义威胁区对象
cTread = CTread;
 
%调整为统一的图形环境
figure('Renderer','opengl');
 
%绘制地形和威胁区
cTerrain.draw();
cTread.draw();
 
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%数据形式说明
%xs 初始点位置集合，三维点集
%yt 目标点位置集合，三维点集
%MxsEna 能力矩阵，集合中的二元组表示最大最小飞行速度
%MxsDIs 最大航程矩阵，集合中的值表示UAV最大可飞的距离
%MytTOrd 任务点间的时序坐标，要求尽量在前面的先执行
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%调用目标分配算法
%AssignType表示不同的模型，1, N=M; 2, N>M; 3, N<M
%Key表示算法执行次数, 0，执行一次，1，执行多次
AssignType =1;
Key = 0;
 
%时间测试
Time = cputime;
 
%%    第一种情况的实验
if AssignType == 1
    %%%%%%%%%%%%%%%%%%%%%%%数据区%%%%%%%%%%%%
    xs = [%11 2 3;2 17 5; 34 26 7;5 2 2;6 33 3;27 34 6;21 45 5;25 12 3;22 17 11;13 31 12;
        %45 21 11; 36 22 11;13 21 16;41 25 10;25 25 11;53 12 14;22 12 12;18 32 11; 67 34 12;26 16 15;
        78 20 15;93 31 12;31 20 13;112 32 15; 134 26 17;45 52 12;36 63 13;67 34 16;21 85 15;32 62 13;];
    yt = [%40 210 12;170 90 13 ; 80 197 21 ;172 120 13;160 56 13;160 143 21;170 200 21;113 200 12;97 134 16;100 145 11;
        %81 101 21; 72 152 12; 133 60 11;84 164 21; 150 90 21;146 121 12; 167 54 12; 108 165 12;99 120 21;110 143 11;
        70 181 12;15 133 13 ; 19 151 12 ;160 192 13;101 120 23;160 113 12;82 101 11;173 140 12;182 65 14; 128 156 12;];
    
    MxsEna=[%0.2 0.3; 0.2 0.4;0.4 0.75;0.3 0.6;0.2 0.3;0.35 0.45;0.3 0.5;0.3 0.6;0.2 0.3; 0.2 0.4;
        % 0.4 0.75;0.3 0.6;0.2 0.3;0.35 0.45;0.3 0.5;0.3 0.6;0.2 0.3;0.35 0.45;0.3 0.5;0.3 0.6;
        0.2 0.3; 0.2 0.4;0.4 0.75;0.3 0.6;0.2 0.3;0.35 0.45;0.3 0.5;0.3 0.6;0.3 0.5;0.3 0.6;];
    
    MxsDis=[%400 700 650 500 700 900 450 610 400 700
        %650 500 700 900 450 610 700 900 450 610
        500 700 300 350 700 900 450 610 450 610];
    
    
    ytW = [%1 3 4 2 1 1 3 2 1 2
        %3 2 1 3 2 1 2 3 2 1
        1 1 1 1 1 1 1 1 1 1];
    
    MytTOrd=[3 4; 5 2; 6 8; 7 4;];
    
    %引入新的协同矩阵，最大起始时间矩阵
    Twait = [];
    Twindow=[];
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    if Key == 0
        
        
        %%%%%%%%%%%%%%%%%%%%%%算法执行区%%%%%%%%%%%%%%%%%%
        %调用目标分配算法，并取得解和相关数据,计算算法耗费时间
        %定义目标分配算法对象
        DDE1 = DMDEAssignTarget(cTerrain,xs,yt,MxsEna,MxsDis,ytW,MytTOrd,Twait,Twindow,1);
        % MC = Cal1.GetMcost();
        %Cal1 = CallocationTargetOld1(cTerrain);
        %标签显示文字
        xlabel('X/km')
        ylabel('Y/km')
        zlabel('Z/km')
        %DDEAssign1=DDEAssignTarget1(cTerrain,xs,yt,MxsEna,MxsDis,ytW,MytTOrd);
        %执行run并返回结果集
        t1 = cputime;
        [solU solT solC solF ] = DDE1.run();
        Time = cputime -t1;
        
        %         %结果
        solU
        solT
        solC
        solF
        Time
        
        
    else if Key == 1
            for i=1:10
                DDE1 = DMDEAssignTarget(cTerrain,xs,yt,MxsEna,MxsDis,ytW,MytTOrd,Twait,Twindow,i);
                % MC = Cal1.GetMcost();
                %Cal1 = CallocationTargetOld1(cTerrain);
                %标签显示文字
                xlabel('X/km')
                ylabel('Y/km')
                zlabel('Z/km')
                %DDEAssign1=DDEAssignTarget1(cTerrain,xs,yt,MxsEna,MxsDis,ytW,MytTOrd);
                %执行run并返回结果集
                t1 = cputime;
                [solU solT solC solF ] = DDE1.run();
                solu(i,:) = solU;
                solt(i,:) = solT;
                solc(i,:) = solC;
                solf(i,:) = solF;
                Time(i,:) = cputime -t1;
            end %i
            %         %结果
            solu
            solt
            solc
            solf
            Time
            
        end % key=1
    end % key=0
    
    
    
else if AssignType == 2
        %% 第二种实验的情况
        %%%%%%%%%%%%%%%%%%%%%%%数据区%%%%%%%%%%%%
        xs = [%11 2 3;2 17 5; 34 26 7;5 2 2;6 33 3;27 34 6;21 45 5;25 12 3;22 17 11;13 31 12;
            %45 21 11; 36 22 11;13 21 16;41 25 10;25 25 11;53 12 14;22 12 12;18 32 11; 67 34 12;26 16 15;
            78 20 15;93 31 12;31 20 13;112 32 15; 134 26 17;45 52 12;36 63 13;67 34 16;21 85 15;32 62 13;];
        yt = [%40 210 12;170 90 13 ; 80 197 21 ;172 120 13;160 56 13;160 143 21;170 200 21;113 200 12;97 134 16;100 145 11;
            %81 101 21; 72 152 12; 133 60 11;84 164 21; 150 90 21;146 121 12; 167 54 12; 108 165 12;99 120 21;110 143 11;
            70 181 12;15 133 13 ; 19 151 12 ;160 192 13 ];
        
        MxsEna=[%0.2 0.3; 0.2 0.4;0.4 0.75;0.3 0.6;0.2 0.3;0.35 0.45;0.3 0.5;0.3 0.6;0.2 0.3; 0.2 0.4;
            % 0.4 0.75;0.3 0.6;0.2 0.3;0.35 0.45;0.3 0.5;0.3 0.6;0.2 0.3;0.35 0.45;0.3 0.5;0.3 0.6;
            0.2 0.3; 0.2 0.4;0.4 0.75;0.3 0.6;0.2 0.3;0.35 0.45;0.3 0.5;0.3 0.6;0.3 0.5;0.3 0.6;];
        
        MxsDis=[%400 700 650 500 700 900 450 610 400 700
            %650 500 700 900 450 610 700 900 450 610
            500 700 300 350 700 900 450 610 450 610];
        
        
        
        ytW = [%1 3 4 2 1 1 3 2 1 2
            %3 2 1 3 2 1 2 3 2 1
            1 1 1 1];
        
        %这个数据里都代表的是目标点，是目标点与目标点的关系
        MytTOrd=[3 2];
        
        %引入新的协同矩阵，最大起始时间矩阵
        Twait = [30 40 20 15 50 70 30 80 90 65];
        Twindow=[];
        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
        if Key == 0
            
       
        else if Key == 1
                for i=1:20
                    DDE1 = DMDEAssignTarget(cTerrain,xs,yt,MxsEna,MxsDis,ytW,MytTOrd,Twait,Twindow,i);
                    % MC = Cal1.GetMcost();
                    %Cal1 = CallocationTargetOld1(cTerrain);
                    %标签显示文字
                    xlabel('X/km')
                    ylabel('Y/km')
                    zlabel('Z/km')
                    %DDEAssign1=DDEAssignTarget1(cTerrain,xs,yt,MxsEna,MxsDis,ytW,MytTOrd);
                    %执行run并返回结果集
                    t1 = cputime;
                    [solU solT solC solF ] = DDE1.run();
                    solu(i,:) = solU;
                    solt(i,:) = solT;
                    solc(i,:) = solC;
                    solf(i,:) = solF;
                    Time(i,:) = cputime -t1;
                end %i
                %         %结果
                solu
                solt
                solc
                solf
                Time
                
            end % key=1
        end % key=0
        
        
    else if AssignType == 3
            %% 第三种实验的情况
            %%%%%%%%%%%%%%%%%%%%%%%数据区%%%%%%%%%%%%
            xs = [%11 2 3;2 17 5; 34 26 7;5 2 2;6 33 3;27 34 6;21 45 5;25 12 3;22 17 11;13 31 12;
                %45 21 11; 36 22 11;13 21 16;41 25 10;25 25 11;53 12 14;22 12 12;18 32 11; 67 34 12;26 16 15;
                78 20 15;93 31 12;31 20 13;112 32 15];
            yt = [%40 210 12;170 90 13 ; 80 197 21 ;172 120 13;160 56 13;160 143 21;170 200 21;113 200 12;97 134 16;100 145 11;
                %81 101 21; 72 152 12; 133 60 11;84 164 21; 150 90 21;146 121 12; 167 54 12; 108 165 12;99 120 21;110 143 11;
                70 131 12;15 133 13 ; 19 151 12 ;160 192 13;101 120 23;160 113 12;82 101 11;173 140 12;182 65 14; 45 176 12;];
            
            MxsEna=[%0.2 0.3; 0.2 0.4;0.4 0.75;0.3 0.6;0.2 0.3;0.35 0.45;0.3 0.5;0.3 0.6;0.2 0.3; 0.2 0.4;
                % 0.4 0.75;0.3 0.6;0.2 0.3;0.35 0.45;0.3 0.5;0.3 0.6;0.2 0.3;0.35 0.45;0.3 0.5;0.3 0.6;
                0.2 0.5; 0.3 0.4;0.4 0.75;0.3 0.45];
            
            MxsDis=[%400 700 650 500 700 900 450 610 400 700
                %650 500 700 900 450 610 700 900 450 610
                700 700 900 650];
            
            ytW = [%1 3 4 2 1 1 3 2 1 2
                %3 2 1 3 2 1 2 3 2 1
                1 1 1 1 1 1 1 1 1 1];
            
            MytTOrd=[3 4; 5 2];
            
            %引入新的协同矩阵，最大起始时间矩阵
            Twait = [30 40 20 15];
            Twindow=[1300,6000];
            %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
            if Key == 0
                
                
                %%%%%%%%%%%%%%%%%%%%%%算法执行区%%%%%%%%%%%%%%%%%%
                %调用目标分配算法，并取得解和相关数据,计算算法耗费时间
                %定义目标分配算法对象
                DDE1 = DMDEAssignTarget(cTerrain,xs,yt,MxsEna,MxsDis,ytW,MytTOrd,Twait,Twindow,1);
                % MC = Cal1.GetMcost();
                %Cal1 = CallocationTargetOld1(cTerrain);
                %标签显示文字
                xlabel('X/km')
                ylabel('Y/km')
                zlabel('Z/km')
                %DDEAssign1=DDEAssignTarget1(cTerrain,xs,yt,MxsEna,MxsDis,ytW,MytTOrd);
                %执行run并返回结果集
                t1 = cputime;
                [solU solT solC solF]  = DDE1.run();
                Time = cputime -t1;
                
                %         %结果
                solU
                solT
                solC
                solF
                Time
                
                
            else if Key == 1
                    
                    for i=1:20
                        DDE1 = DMDEAssignTarget(cTerrain,xs,yt,MxsEna,MxsDis,ytW,MytTOrd,Twait,Twindow,i);
                        % MC = Cal1.GetMcost();
                        %Cal1 = CallocationTargetOld1(cTerrain);
                        %标签显示文字
                        xlabel('X/km')
                        ylabel('Y/km')
                        zlabel('Z/km')
                        
                        t1 = cputime;
                        [solU solT solC solF ] = DDE1.run();
                        solu(i,:) = solU;
                        solt(i,:) = solT;
                        solc(i,:) = solC;
                        solf(i,:) = solF;
                        Time(i,:) = cputime -t1;
                    end %i
                    %         %结果
                    solu
                    solt
                    solc
                    solf
                    Time
                end % key=1
            end % key=0
            
        end % type =3
        
    end % type =2
    
end % type=1

• 1.
• 2.
• 3.
• 4.
• 5.
• 6.
• 7.
• 8.
• 9.
• 10.
• 11.
• 12.
• 13.
• 14.
• 15.
• 16.
• 17.
• 18.
• 19.
• 20.
• 21.
• 22.
• 23.
• 24.
• 25.
• 26.
• 27.
• 28.
• 29.
• 30.
• 31.
• 32.
• 33.
• 34.
• 35.
• 36.
• 37.
• 38.
• 39.
• 40.
• 41.
• 42.
• 43.
• 44.
• 45.
• 46.
• 47.
• 48.
• 49.
• 50.
• 51.
• 52.
• 53.
• 54.
• 55.
• 56.
• 57.
• 58.
• 59.
• 60.
• 61.
• 62.
• 63.
• 64.
• 65.
• 66.
• 67.
• 68.
• 69.
• 70.
• 71.
• 72.
• 73.
• 74.
• 75.
• 76.
• 77.
• 78.
• 79.
• 80.
• 81.
• 82.
• 83.
• 84.
• 85.
• 86.
• 87.
• 88.
• 89.
• 90.
• 91.
• 92.
• 93.
• 94.
• 95.
• 96.
• 97.
• 98.
• 99.
• 100.
• 101.
• 102.
• 103.
• 104.
• 105.
• 106.
• 107.
• 108.
• 109.
• 110.
• 111.
• 112.
• 113.
• 114.
• 115.
• 116.
• 117.
• 118.
• 119.
• 120.
• 121.
• 122.
• 123.
• 124.
• 125.
• 126.
• 127.
• 128.
• 129.
• 130.
• 131.
• 132.
• 133.
• 134.
• 135.
• 136.
• 137.
• 138.
• 139.
• 140.
• 141.
• 142.
• 143.
• 144.
• 145.
• 146.
• 147.
• 148.
• 149.
• 150.
• 151.
• 152.
• 153.
• 154.
• 155.
• 156.
• 157.
• 158.
• 159.
• 160.
• 161.
• 162.
• 163.
• 164.
• 165.
• 166.
• 167.
• 168.
• 169.
• 170.
• 171.
• 172.
• 173.
• 174.
• 175.
• 176.
• 177.
• 178.
• 179.
• 180.
• 181.
• 182.
• 183.
• 184.
• 185.
• 186.
• 187.
• 188.
• 189.
• 190.
• 191.
• 192.
• 193.
• 194.
• 195.
• 196.
• 197.
• 198.
• 199.
• 200.
• 201.
• 202.
• 203.
• 204.
• 205.
• 206.
• 207.
• 208.
• 209.
• 210.
• 211.
• 212.
• 213.
• 214.
• 215.
• 216.
• 217.
• 218.
• 219.
• 220.
• 221.
• 222.
• 223.
• 224.
• 225.
• 226.
• 227.
• 228.
• 229.
• 230.
• 231.
• 232.
• 233.
• 234.
• 235.
• 236.
• 237.
• 238.
• 239.
• 240.
• 241.
• 242.
• 243.
• 244.
• 245.
• 246.
• 247.
• 248.
• 249.
• 250.
• 251.
• 252.
• 253.
• 254.
• 255.
• 256.
• 257.
• 258.
• 259.
• 260.
• 261.
• 262.
• 263.
• 264.
• 265.
• 266.
• 267.
• 268.
• 269.
• 270.
• 271.
• 272.
• 273.
• 274.
• 275.
• 276.