差分进化算法(DE)求函数最小值

差分进化算法求函数 Z = 3 * cos(X .* Y) + X + Y , -4 <= X <= 4, -4 <= Y <= 4。

函数图片

计算目标函数值

计算目标函数值的函数:

function z = calobj (pop)
% 计算目标函数值
% pop       input  种群
% z         output 目标函数值
z = 3 * cos(pop(:,1) .* pop(:,2)) + pop(:,1) + pop(:,2);
end

初始化种群

目标函数有两个参数,生成每个个体有两个基因的种群:

function pop = initpop(popsize, chromlength, xl, xu)
% 生成初始种群
% popsize           input  种群规模
% chromlengt        input  染色体长度
% xl                input  x下限
% xu                input  x上限
% pop               output 种群
pop = rand(popsize, chromlength) * (xu - xl) + xl;
end

变异

变异函数如下:

function mutationpop = mutation (pop, F)
% 变异操作
% pop           input  种群
% F             input  缩放因子
% mutationpop   output 变异后种群
[popsize, chromlength] = size(pop);
mutationpop = zeros(popsize, chromlength);
for i = 1:popsize
    % 取3个互异的索引 r0 r1 r2
    r = randperm(popsize);
    index = find (r ~= i);
    rn = r(index(1:3));
    r0 = rn(1); r1 = rn(2); r2 = rn(3);
%     fprintf('i = %d, r0 = %d, r1 = %d, r2 = %d\n', i, r0, r1, r2);

    mutationpop(i,:) = pop(r0,:) + F .* (pop(r1,:) - pop(r2,:));
end
end

交叉

交叉函数如下:

function crossoverpop = crossover(pop, mpop, cr)
% 交叉
% pop           input  种群
% mpop          input  变异后的种群
% cr            input  交叉概率
% crossoverpop  output 交叉后的种群
[popsize, chromlength] = size(pop);
crossoverpop = mpop;
r = rand(popsize, chromlength);
index = find (r > cr);
crossoverpop(index) = pop(index);
jrand = randi(chromlength, 1, popsize);
crossoverpop(sub2ind(size(crossoverpop), [1:popsize], jrand)) ...
    = mpop(sub2ind(size(mpop), [1:popsize], jrand));
end

在交叉操作之后,应约束边界:

function newpop = constrictboundary(pop, xl, xu)
% 约束边界(边界吸收)
% pop       input  种群
% xl        input  自变量最小值(包含)
% xu        input  自变量最大值(包含)
% newpop    output 约束边界后的种群
newpop = pop;
newpop(newpop < xl) = xl;
newpop(newpop > xu) = xu;
end

选择

function newpop = selection(pop, npop)
% 选择(小值优化)
% pop           input  种群1(原始种群)
% pop           input  种群2(变异-交叉种群)
% newpop        output 选择后的种群
newpop = pop;
index = find(calobj(npop) <= calobj(pop));
newpop(index, :) = npop(index, :);
end

主程序

主程序如下:

clc;
clear;

NP  = 20;       % 种群规模
D   = 2;        % 参数个数
G   = 30;       % 最大进化代数
F   = 0.5;      % 缩放因子
Cr  = 0.8;      % 交叉因子

xl  = -4;       % x下限(也是y下限)
xu  = 4;        % x上限(也是y上限)

bestvalue = zeros(3, G);

% 优化
gen = 0;
pop = initpop(NP, D, xl, xu);
objvalue = calobj(pop);
while gen < G
    mpop = mutation(pop, F);                    % 变异
    cpop = crossover(pop, mpop, Cr);            % 交叉
    cpop = constrictboundary(cpop, xl, xu);     % 约束边界
    pop = selection(pop, cpop);                 % 选择
    objvalue = calobj(pop);
    gen = gen + 1;

    % 记录最优
    [~, index] = min(objvalue);
    bestvalue(1:2, gen) = pop(index,:)';
    bestvalue(3,gen) = objvalue(index);
end

fprintf('bestX = %f, bestY = %f, bestZ = %f\n', ...
    bestvalue(1,end), bestvalue(2,end), bestvalue(3,end));

% 绘图
figure(1);
x = [-4:0.1:4]; y = [-4:0.1:4];
[X, Y] = meshgrid(x, y);
Z = 3 * cos(X .* Y) + X + Y;
surf(X, Y, Z);
hold on;
scatter3(bestvalue(1,:), bestvalue(2,:), bestvalue(3,:), ...
    'MarkerEdgeColor','k', 'MarkerFaceColor',[0 .75 .75]);
xlabel('x'); ylabel('y'); zlabel('z'); title('函数图');
hold off;

figure(2);
plot(bestvalue(3,:));
xlabel('进化代数'); ylabel('最优目标函数值'); title('目标函数值变化图');

执行结果

bestX = -3.947841, bestY = -4.000000, bestZ = -10.937414

最优点图

目标值变化图

加入对话

6条评论

  1. 你好,最近在用差分进化优化vmd分解,如果两个自变量x与y的取值区间分别是x为3:7,y
    y为50:50:500
    源程序中的xl = -4; % x下限(也是y下限) xu = 4; % x上限(也是y上限)该怎么修改呢

    1. 你好,今天才看到你的留言,不好意思回复晚了。
      我不知道vmd分解的问题,也没有研究过,所以很抱歉无法就该问题提供更多的想法了。
      差分进化算法多用于解决连续问题,如果自变量是在连续域变化的,那么是直接修改代码中自变量的上下界就可以。
      从您说的x=3:7和y=50:50:500中,我猜测x是{3 4 5 6 7}中的一个元素,y是{50 100 150 … 500}中的一个元素,这并不是连续域问题,是离散域中的问题,所以肯定无法通过简单的更改上下界处理。
      离散域可以是多变量无关的(x的取值和y的取值互不影响)和多变量相关的(x的取值和y的取值相互影响)。我有2个简单的想法,只是想法,不能保证算法性能。
      1. 对于多变量无关的情况,可以把变量取成连续域,然后通过解码变成离散域,带入目标函数求值。比如x的取值为x={x|2http://mwangblog.com/?p=1712
      http://mwangblog.com/?p=1714
      http://mwangblog.com/?p=1721
      http://mwangblog.com/?p=1723
      希望这些对你有所帮助。

      1. 还有一个问题,我这个程序运行起来需要20分钟左右,怎么能缩短优化进行的时间,需要改哪些差分参数呢

        1. 你好,
          博客上的程序(http://mwangblog.com/?p=1660),在我的电脑上大概只需要0.2秒。
          如果您指的是文章里程序,那么怎么也不会运行二十分钟的,运行这么长时间的可能是您代码写错了,比如参数设置错误等。
          如果您指的是您修改后的程序,那么我也不清楚是怎么回事儿。
          您可以把程序拷贝到MATLAB而不要做修改试一试。

      2. 您好!如果-4<=x<=4,-1<=y<=2(即x和y的取值范围不同)该怎么办呢??亦或目标函数是七个自变量的函数并且每个自变量的取值范围都不一样该怎么处理呢??种群初始化、边界条件等该怎么修改呢??期待您的回复,万分感谢。

        1. 你好,已在简书回复,回复内容如下:

          你好,能用DE处理,但是函数需要改写。
          在本文中,自变量的取值范围在初始化种群时用了一次,在交叉之后约束边界时用了一次。
          假设x和y的边界不一样,染色体的意义为[ x值 y值],以交叉后的约束边界函数为例:
          function newpop = constrictboundary(pop, xl, xu, yl, yu) % 增加两个参数yl yu
          newpop = pop;
          newpop中第1列(x列)中大于xu的限制为xu
          newpop中第1列(x列)中小于xl的限制为xl
          newpop中第2列(y列)中大于yu的限制为yu
          newpop中第2列(y列)中小于yl的限制为yl
          end

          这种方法对于两个参数是可以的,但是对于多个参数(比如你说的7个)是繁琐的。对于多个参数,可以设定一个取值范围的变量,比如叫limit, limit=[xl xu; yl yu; zl zu] ,这是3个变量的情况,每行两个值,分别是变量的下限和上限。染色体的意义是[x值 y值 z值],还是以约束函数为例:
          function newpop = constrictboundary(pop, limit) % 用limit变量表示取值边界
          newpop = pop;
          对于newpop中的每一列,假设是 i 列
          该列最小值为limit中 i 行的第1个值
          该列的最大值为limit中 i 行的第2个值
          超出边界的进行约束
          end
          end
          使用一个循环可以解决。

留下评论

电子邮件地址不会被公开。 必填项已用*标注