双层优化入门（3）-基于智能优化算法的解决方法（附 Matlab 代码）
前面两篇博客介绍了双层优化的基本原理和使用KKT条件求解双层优化的方法，以及使用yalmip工具箱求解双层优化的方法。
除了数学规划方法之外，还可采用智能优化算法求解双层优化问题，一般在上层优化中采用智能优化算法，下层优化使用数学规划方法；也可以在上下层优化中都采用智能优化算法，这篇博客将进行详细介绍。
算例依旧使用上面两篇博客中的线性双层优化问题，由于这个优化问题比较简单，我们采用最基础的粒子群算法进行求解。
1.粒子群算法
1995年，受到鸟群觅食行为的规律性启发，James Kennedy和Russell Eberhar建立了一个简化算法模型，经过多年改进最终形成了粒子群优化算法(Particle Swarm Optimization, PSO) ，也可称为粒子群算法。
关于算法的原理及步骤可以参考这篇博客：
粒子群优化算法(Particle Swarm Optimization, PSO)的详细解读 - 知乎 (zhihu.com)
我们假设目标函数是平方和最小，也就是：
假设种群数为200，问题维度为10，最大迭代次数为200，位置上下限分别为10，-10，速度上下限为5，-5，惯性权重取固定值0.9，自我学习因子以及群体学习因子都取1.5，Matlab代码如下：
%% 清除变量 clc clear close all warning off %% 设置种群参数 sizepop = 200; % 初始种群个数 dim = 10; % 空间维数 ger = 200; % 最大迭代次数 x_max = 10*ones(1,dim); % 位置上限 x_min = -10*ones(1,dim); % 位置下限 v_max = 5*ones(1,dim); % 速度上限 v_min = -5*ones(1,dim); % 速度下限 w=0.9; % 惯性权重 c_1 = 1.5; % 自我学习因子 c_2 = 1.5; % 群体学习因子 %% 种群初始化 pop=x_min+rand(sizepop,dim).*(x_max-x_min); % 初始化种群 pop_v=v_min+rand(sizepop,dim).*(v_max-v_min); % 初始化种群速度 pop_zbest=pop(1,:); % 初始化群体最优位置 pop_gbest=pop; % 初始化个体最优位置 fitness=zeros(1,sizepop); % 所有个体的适应度 fitness_zbest=inf; % 初始化群体最优适应度 fitness_gbest=inf*ones(1,sizepop); % 初始化个体最优适应度 % 初始的适应度 for k=1:sizepop % 计算适应度值 fitness(k)=sum(pop(k,:).^2); if fitness(k)<fitness_zbest fitness_zbest=fitness(k); pop_zbest=pop(k,:); end end history_pso=zeros(1,ger); % 粒子群历史最优适应度值 %% 迭代求最优解 iter=1; while iter <= ger for k=1:sizepop % 更新速度并对速度进行边界处理 pop_v(k,:)= w * pop_v(k,:) + c_1*rand*(pop_gbest(k,:)-pop(k,:))+c_2*rand*(pop_zbest-pop(k,:)); for kk=1:dim if pop_v(k,kk) > v_max(kk) pop_v(k,kk) = v_max(kk); end if pop_v(k,kk) < v_min(kk) pop_v(k,kk) = v_min(kk); end end % 更新位置并对位置进行边界处理 pop(k,:)=pop(k,:)+pop_v(k,:); for kk=1:dim if pop(k,kk) > x_max(kk) pop(k,kk) = x_max(kk); end if pop(k,kk) < x_min(kk) pop(k,kk) = x_min(kk); end end % 更新适应度值 fitness(k)=sum(pop(k,:).^2); if fitness(k)<fitness_zbest fitness_zbest=fitness(k); pop_zbest=pop(k,:); end if fitness(k)<fitness_gbest(k) fitness_gbest(k)=fitness(k); pop_gbest(k,:)=pop(k,:); end end history_pso(iter)=fitness_zbest; disp(['PSO第',num2str(iter),'次迭代最优适应度=',num2str(fitness_zbest)]) iter=iter+1; end disp(['最优解：x=',num2str(pop_zbest)]) disp(['最优函数值=',num2str(fitness_zbest)]) plot(history_pso,'linewidth',1) ylabel('最优适应度值') xlabel('迭代次数')
运行结果如下：
我们知道最优解是当所有变量都取0时，最优函数值为0，上面显示采用粒子群算法求出的显然并不是全局最优解，只是一个局部最优解。这也是智能优化算法无法避免的问题，即使是一个非常简单的目标函数，求出的结果也无法保证是全局最优，那么当目标函数变复杂时，情况将会更糟糕。现在对智能优化算法的研究非常多，各种动植物园算法、各种改进都层出不穷，但还是无法从根本上解决算法无法保证全局收敛的问题。
所以，只有在数学模型比较复杂，非线性条件很多，而且对结果的误差是可以接受的情况下，才建议使用智能优化算法进行求解。
下面以粒子群算法为例，介绍采用智能优化算法求解双层优化问题的方法。
2.智能优化算法中对约束条件的处理
上面提供的代码中，含有的约束条件是x∈[-10,10]，直接可以通过粒子位置的上下限进行约束，但其他形式的约束就没办法这样处理了，常用的处理方式有三类，以上面提到的双层优化问题的上层优化为例，分别介绍这几种方法：
1）通过加罚函数的方式将有约束的问题转为无约束问题：
首先定义罚函数：
然后将罚函数加入原来的目标函数中形成新的目标函数：
其中，p是惩罚因子。这样处理，可以使得违反约束的粒子目标函数值变得很大，从而促使粒子朝着遵守约束的区域聚集。这种方式在过去非常常用，但我自己不太喜欢这种方式，原因待会我会解释。假设惩罚因子为1，新的目标函数可以用matlab代码表示出来：
function fitness=fitness_fun_method1(pop) x=pop(1); y=pop(2); punish=0; p=1; if -2*x+3*y > 12 punish=punish+p*(-2*x+3*y-12); end if x+y > 14 punish=punish+p*(x+y-14); end fitness=-x-2*y+punish; end
运行结果如下：
上一篇：
谷歌翻译 谷歌翻译 API 教程 PHP