多目标规划_matlab程序

如何使用Matlab进行最优化和多目标优化问题求解Matlab是一种强大的数学计算工具，广泛应用于各个领域的科学研究和工程实践中。

其中，最优化和多目标优化问题的求解是Matlab的一项重要功能。

本文将介绍如何使用Matlab进行最优化和多目标优化问题的求解，并提供一些实际应用案例。

一、最优化问题求解最优化问题求解是指在给定的约束条件下，寻找一个使得目标函数取得最大（或最小）值的变量组合。

Matlab提供了多种最优化算法，如线性规划、二次规划、非线性规划等。

下面以非线性规划为例，介绍如何使用Matlab进行最优化问题的求解。

1. 准备工作在使用Matlab进行最优化问题求解之前，需要先定义目标函数和约束条件。

目标函数是最优化问题的核心，可以是线性的或非线性的。

约束条件可以是等式约束或不等式约束。

同时，还需要确定变量的取值范围和初值。

2. 选择合适的算法Matlab提供了多个最优化算法，根据问题的特点选择合适的算法是非常重要的。

常用的算法有fmincon、fminunc、fminsearch等。

例如，fmincon函数适用于求解具有约束条件的非线性规划问题，而fminunc函数适用于求解无约束或有约束的非线性规划问题。

3. 调用相应的函数根据选择的算法，调用相应的函数进行求解。

以fmincon函数为例，其调用方式为：```[x, fval] = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options)```其中，fun为目标函数，x0为变量的初值，A、b为不等式约束矩阵和向量，Aeq、beq为等式约束矩阵和向量，lb、ub为变量的下界和上界，nonlcon为非线性约束函数，options为求解选项。

4. 解析结果求解完成后，可以通过解析结果来评估求解器的性能。

Matlab提供了fval和exitflag两个输出参数，其中fval表示最优解的目标函数值，exitflag表示求解器的退出标志。

Gurobi多目标问题在Matlab中的解决一、Gurobi简介Gurobi是一款强大的商业数学建模工具，广泛应用于优化领域。

它提供了多种优化算法，能够高效地解决线性规划、整数规划、二次规划等各种优化问题。

在实际工程和科学研究中，经常遇到多目标优化问题，即需要同时优化多个目标函数。

本文将介绍如何使用Gurobi在Matlab中解决多目标优化问题。

二、多目标优化问题的定义在多目标优化问题中，我们需要最小化或最大化多个目标函数，而且这些目标函数之间往往存在相互矛盾的关系。

在生产计划中，一个目标函数可能是最大化产量，另一个目标函数可能是最小化成本。

在实际应用中，我们需要找到一组可行的解，使得所有目标函数都达到一个较好的平衡。

三、Gurobi在Matlab中的调用在Matlab中调用Gurobi需要先安装Gurobi的Matlab接口。

安装完成后，我们可以在Matlab命令窗口中输入命令"gurobi"来验证是否成功安装。

接下来，我们需要在Matlab中编写代码，定义优化问题的目标函数、约束条件和变量类型。

在定义目标函数时，我们需要考虑多个目标函数之间的相关性，以及它们之间的权重关系。

在定义约束条件和变量类型时，我们需要考虑多目标函数之间可能存在的约束条件和变量之间的相互制约关系。

四、多目标优化问题的解决方法Gurobi提供了多种解决多目标优化问题的方法，包括加权法、约束法和Pareto最优解法等。

在加权法中，我们将多个目标函数进行线性组合，并引入权重因子来平衡各个目标函数之间的重要性。

在约束法中，我们将多个目标函数作为多个约束条件，通过逐步添加约束条件来找到最优解。

在Pareto最优解法中，我们寻找一组可行解，使得没有其他可行解能比它在所有目标函数上都更好。

五、案例分析以生产计划为例，假设我们需要同时考虑最大化产量和最小化成本两个目标。

我们可以先使用加权法，通过调整权重因子来平衡这两个目标的重要性，找到一个较好的解。

使用Matlab进行多目标优化和约束优化引言：多目标优化和约束优化是现代科学和工程领域中的重要问题。

在很多实际应用中，我们常常面对的是多个目标参数之间存在冲突的情况，同时还需要满足一定的约束条件。

这就需要我们采用适当的方法和工具进行多目标优化和约束优化。

本文将介绍如何使用Matlab进行多目标优化和约束优化。

一、多目标优化多目标优化是指在优化问题中存在多个目标函数，我们的目标是同时优化这些目标函数。

在Matlab中，可以使用多种方法进行多目标优化，其中常用的方法包括遗传算法、粒子群算法和模拟退火等。

1.1 遗传算法遗传算法是一种模拟生物进化过程的优化算法。

它模拟了遗传的过程，通过交叉、变异和选择等操作，利用群体中不断进化的个体来搜索最优解。

在多目标优化中，遗传算法常用于生成一组非支配解，即没有解能同时优于其他解的情况。

Matlab中提供了相关的工具箱，如Global Optimization Toolbox和Multiobjective Optimization Toolbox，可以方便地进行多目标优化。

1.2 粒子群算法粒子群算法是一种基于群体行为的优化算法。

它通过模拟鸟群或鱼群等群体的行为，寻找最优解。

在多目标优化中，粒子群算法也可以生成一组非支配解。

Matlab中的Particle Swarm Optimization Toolbox提供了相关函数和工具，可以实现多目标优化。

1.3 模拟退火模拟退火是一种模拟金属冶炼过程的优化算法。

它通过模拟金属在高温下退火的过程，通过温度控制来逃离局部最优解，最终达到全局最优解。

在多目标优化中，模拟退火算法可以通过调整温度参数来生成一组非支配解。

Matlab中也提供了相关的函数和工具，可以进行多目标优化。

二、约束优化约束优化是指在优化问题中存在一定的约束条件，我们的目标是在满足这些约束条件的前提下，使目标函数达到最优。

在Matlab中，也有多种方法可以进行约束优化，其中常用的方法包括罚函数法、惩罚函数法和内点法等。

如何在MATLAB中进行多目标优化多目标优化问题是指在存在多个冲突目标的情况下，求解一个能够同时最小化或最大化多个目标函数的问题。

在实际应用中，多目标优化问题被广泛应用于工程优化、金融投资、交通规划等领域。

在MATLAB中，有多种方法可以用来解决多目标优化问题，本文将介绍其中的几种常用方法。

一、多目标优化问题的定义在开始使用MATLAB进行多目标优化之前，首先需要明确多目标优化问题的数学定义。

一般而言，多目标优化问题可以表示为：```minimize f(x) = [f1(x), f2(x), ..., fm(x)]subject to g(x) ≤ 0, h(x) = 0lb ≤ x ≤ ub```其中，f(x)为多个目标函数，g(x)和h(x)为约束条件，lb和ub分别为决策变量的下界和上界。

问题的目标是找到一组决策变量x，使得目标函数f(x)取得最小值。

二、多目标优化问题的解法在MATLAB中，有多种方法可以用来解决多目标优化问题。

下面将介绍其中的几种常见方法。

1. 非支配排序遗传算法（Non-dominated Sorting Genetic Algorithm，NSGA）NSGA是一种经典的多目标优化算法，它将候选解集划分为多个等级或层次，从而使得每个解在候选解集内具备非劣势性。

在MATLAB中，可以使用多目标遗传算法工具箱（Multi-Objective Optimization Toolbox）中的`gamultiobj`函数来实现NSGA算法。

该函数可以通过指定目标函数、约束条件和决策变量范围等参数来求解多目标优化问题。

2. 多目标粒子群优化算法（Multi-objective Particle Swarm Optimization，MOPSO）MOPSO是一种基于群体智能的多目标优化算法，它模拟了粒子的行为，通过不断迭代寻找最优解。

在MATLAB中，可以使用多目标粒子群优化工具箱（Multi-Objective Particle Swarm Optimization Toolbox）中的`mopso`函数来实现MOPSO算法。

如何使用Matlab进行多目标优化使用Matlab进行多目标优化概述：多目标优化是在现实问题中常见的一种优化方法，即需要优化多个目标函数，而非只有一个目标函数。

这篇文章将介绍如何使用Matlab进行多目标优化，包括问题建模、求解方法和实例分析。

1. 问题建模在进行多目标优化之前，需要将实际问题建模为数学模型。

首先，明确问题的决策变量和目标函数。

决策变量是需要优化的参数或变量，而目标函数是需要最小化或最大化的指标。

例如，我们要优化一个生产系统的成本和产量，可以将成本设为一个目标函数，产量设为另一个目标函数。

2. 目标权重设定由于多目标优化存在矛盾或折衷的情况，需要设定目标函数的权重。

权重反映了各个目标函数的重要性，较高的权重意味着对应的目标更重要。

例如，在上述生产系统的例子中，如果成本比产量更重要，可以给成本赋予较高的权重。

3. 多目标优化求解方法Matlab提供了多种多目标优化求解方法，常用的有基于进化算法的优化方法，例如遗传算法、粒子群优化算法等。

这些方法通过不断迭代搜索解空间，逐步找到最优解。

以下是使用Matlab进行多目标优化的一般步骤：a) 定义优化问题的问题函数，包括目标函数和约束条件。

b) 设定优化问题的求解选项，例如优化算法、迭代次数和收敛准则等。

c) 运行优化求解器，获得最优解或近似最优解。

d) 对求解结果进行分析和评价。

4. 多目标优化实例分析为了更好地理解如何使用Matlab进行多目标优化，我们以一个简单的例子进行分析。

假设有一个三维空间内的旅行商问题，即找到一条路径，使得旅行距离最短、花费最少以及时间最短。

我们可以将问题建模为一个三目标优化问题：目标一：最小化旅行距离。

目标二：最小化旅行花费。

目标三：最小化旅行时间。

通过定义目标函数和约束条件，我们可以使用Matlab的多目标优化求解器，如gamultiobj函数，来获得近似最优解。

在求解过程中，可以通过设置收敛准则、种群大小等选项来调节求解参数。

优化与决策——多目标线性规划的若干解法及MATLAB 实现摘要：求解多目标线性规划的基本思想大都是将多目标问题转化为单目标规划，本文介绍了理想点法、线性加权和法、最大最小法、目标规划法，然后给出多目标线性规划的模糊数学解法，最后举例进行说明，并用Matlab 软件加以实现。

关键词：多目标线性规划 Matlab 模糊数学。

注：本文仅供参考，如有疑问，还望指正。

一．引言多目标线性规划是多目标最优化理论的重要组成部分，由于多个目标之间的矛盾性和不可公度性,要求使所有目标均达到最优解是不可能的,因此多目标规划问题往往只是求其有效解(非劣解)。

目前求解多目标线性规划问题有效解的方法,有理想点法、线性加权和法、最大最小法、目标规划法。

本文也给出多目标线性规划的模糊数学解法。

二．多目标线性规划模型多目标线性规划有着两个和两个以上的目标函数,且目标函数和约束条件全是线性函数,其数学模型表示为：11111221221122221122max n n n nr r r rn nz c x c x c x z c x c x c x z c x c x c x =+++⎧⎪=+++⎪⎨ ⎪⎪=+++⎩ （1）约束条件为：1111221121122222112212,,,0n n n n m m mn n mn a x a x a x b a x a x a x b a x a x a x bx x x +++≤⎧⎪+++≤⎪⎪ ⎨⎪+++≤⎪≥⎪⎩ （2） 若（1）式中只有一个1122i i i in n z c x c x c x =+++ ，则该问题为典型的单目标线性规划。

我们记:()ij m n A a ⨯=,()ij r n C c ⨯=,12(,,,)T m b b b b = ,12(,,,)T n x x x x = ，12(,,,)T r Z Z Z Z = .则上述多目标线性规划可用矩阵形式表示为:max Z Cx =约束条件：0Ax bx ≤⎧⎨≥⎩（3）三．MATLAB 优化工具箱常用函数[3]在MA TLAB 软件中，有几个专门求解最优化问题的函数，如求线性规划问题的linprog 、求有约束非线性函数的fmincon 、求最大最小化问题的fminimax 、求多目标达到问题的fgoalattain 等,它们的调用形式分别为：①.[x,fval]=linprog(f,A,b,Aeq,beq,lb,ub)f 为目标函数系数,A,b 为不等式约束的系数, Aeq,beq 为等式约束系数, lb,ub 为x 的下限和上限, fval 求解的x 所对应的值。

多目标非线性规划程序M a t l a bDocument serial number【NL89WT-NY98YT-NC8CB-NNUUT-NUT108】f u n c t i o n[e r r m s g,Z,X,t,c,f a i l]=BNB18(fun,x0,xstat,xl,xu,A,B,Aeq,Beq,nonlcon,setts,options1,options2,maxSQPit,varargin );%·Dêy1￡Díóa·§¨μü′ú·¨￡úDê1ó￡DèOptimization toolbox §3% Minimize F(x)%subject to: xlb <= x <=xub% A*x <= B% Aeq*x=Beq% C(x)<=0% Ceq(x)=0%% x(i)éaáD±á￡êy￡ò1ì¨μ% ê1óê%[errmsg,Z,X]=BNB18('fun',x0,xstat,xl,xu,A,B,Aeq,Beq,'nonlcon',setts)%fun￡o Mt￡±íê×Dˉ±êoˉêyf=fun(x)%x0: áDòá￡±íê±á3μ%xstat￡o áDòá￡xstat(i)=0±íêx(i)aáD±á￡1±íêêy￡2±íê1ì¨μ%xl￡o áDòá￡±íê±á%xu: áDòá￡±íê±áé%A: ó, ±íêD2μèêêμêy%B: áDòá, ±íêD2μèêêé%Aeq: ó, ±íêDμèêêμêy%Beg: áDòá, ±íêD2μèêêóòμ%nonlcon: Mt￡±íê·Dêoˉêy[C,Ceq]=nonlin(x),DC(x)a2μèêê,% Ceq(x)aμèêê%setts: ·¨éè%errmsq: ·μ′íóìáê%Z: ·μ±êoˉêy×Dμ%X: ·μ×óa%%àyìa% max x1*x2*x3% -x1+2*x2+2*x3>=0% x1+2*x2+2*x3<=72% 10<=x2<=20% x1-x2=10% èD′ Moˉêy% function f=discfun(x)% f=-x(1)*x(2)*x(3);%óa% clear;x0=[25,15,10]';xstat=[1 1 1]';% xl=[20 10 -10]';xu=[30 20 20]';% A=[1 -2 -2;1 2 2];B=[0 72]';Aeq=[1 -1 0];Beq=10;% [err,Z,X]=BNB18('discfun',x0,xstat,xl,xu,A,B,Aeq,Beq);% XMAX=X',ZMAX=-Z%% BNB18 Finds the constrained minimum of a function of several possibly integer variables.% Usage: [errmsg,Z,X,t,c,fail] =%BNB18(fun,x0,xstatus,xlb,xub,A,B,Aeq,Beq,nonlcon,settings,options1,options2,maxSQPiter ,P1,P2,...)%% BNB solves problems of the form:% Minimize F(x) subject to: xlb <= x0 <=xub% A*x <= B Aeq*x=Beq% C(x)<=0 Ceq(x)=0% x(i) is continuous for xstatus(i)=0% x(i) integer for xstatus(i)= 1% x(i) fixed for xstatus(i)=2%% BNB uses:% Optimization Toolbox Version (R11) 09-Oct-1998% From this toolbox is called. For more info type help fmincon.%% fun is the function to be minimized and should return a scalar. F(x)=feval(fun,x).% x0 is the starting point for x. x0 should be a column vector.% xstatus is a column vector describing the status of every variable x(i).% xlb and xub are column vectors with lower and upper bounds for x.% A and Aeq are matrices for the linear constrains.% B and Beq are column vectors for the linear constrains.% nonlcon is the function for the nonlinear constrains.% [C(x);Ceq(x)]=feval(nonlcon,x). Both C(x) and Ceq(x) should be column vectors.%% errmsg is a string containing an error message if BNB found an error in the input.% Z is the scalar result of the minimization, X the values of the accompanying variables.% t is the time elapsed while the algorithm BNB has run, c is the number of BNB cycles and% fail is the number of unsolved leaf sub-problems.%% settings is a row vector with settings for BNB:% settings(1) (standard 0) if 1: use phase 1 by relaxation. This sometimes makes the algorithm% faster, because phase 1 means the algorithm first checks if there is a feasible solution% for a sub-problem before trying to find a best solution. If there is no feasible solution BNB% will not try to find a best solution.% settings(2) (standard 0) if 1: if the sub-problem did not converge do not branch. If a sub-% problem did not converge this means BNB did not find a solution for it. Normally BNB will% branch the problem so it can try again to find a solution.% A sub-problem that is a leaf of the branch-and-bound-three can not be branched. If such% a problem does not converge it will be considered unfeasible and the parameter fail will be% raised by one.% settings(3) (standard 0) if 1: if 1 a sub-problem that did not converge but did return a feasible% point will be considered convergent. This might be useful if fmincon is having a hard time with% a certain problem but you do want some results.% options1 and options2 are options structures for phase 1 and phase 2.% For details about the options structure type help optimset.% maxSQPiter is a global variable used by fmincon (if modified as described in .% maxSQPiter is 1000 by default.% P1,P2,... are parameters to be passed to fun and nonlcon.% F(x)=feval(fun,x,P1,P2,...). [C(x);Ceq(x)]=feval(nonlcon,x,P1,P2,...).% Type edit BNB18 for more info.% . Kuipers% e-mail% FI-Lab% Applied Physics% Rijksuniversiteit Groningen% To get rid of bugs and to stop fmincon from hanging make the following chances:%% In optim/private/ ($Revision: $ $Date: 1998/08/24 13:46:15 $):% Get EXITFLAG independent of verbosity.% After the lines: disp(' less than 2* but constraints are not satisfied.')% end% EXITFLAG = -1;% end% end% status=1;% add the line: if (strncmp(howqp, 'i',1) & mg > 0), EXITFLAG = -1; end;%% In optim/private/ ($Revision: $ $Date: 1998/09/01 21:37:56 $):% Stop qpsub from hanging.% After the line: % Andy Grace 7-9-90. Mary Ann Branch 9-30-96.% add the line: global maxSQPiter;% and changed the line: maxSQPiters = Inf;% to the line: if exist('maxSQPiter','var'), maxSQPiters = maxSQPiter; else maxSQPiters=inf; end;% I guess there was a reason to put maxSQPiters at infinity, but this works fine for me.global maxSQPiter;% STEP 0 CHECKING INPUTZ=[]; X=[]; t=0; c=0; fail=0;if nargin<2, errmsg='BNB needs at least 2 input arguments.'; return; end;if isempty(fun), errmsg='No fun found.'; return; end;if isempty(x0), errmsg='No x0 found.'; return;elseif size(x0,2)>1, errmsg='x0 must be a column vector.'; return; end;xstatus=zeros(size(x0));if nargin>2 & ~isempty(xstat)if all(size(xstat)<=size(x0))xstatus(1:size(xstat))=xstat;else errmsg='xstatus must be a column vector the same size as x0.'; return;end;if any(xstatus~=round(xstatus) | xstatus<0 | 2<xstatus)errmsg='xstatus must consist of the integers 0,1 en 2.'; return;end;end;xlb=zeros(size(x0));xlb(find(xstatus==0))=-inf;if nargin>3 & ~isempty(xl)if all(size(xl)<=size(x0))xlb(1:size(xl,1))=xl;else errmsg='xlb must be a column vector the same size as x0.'; return;end;end;if any(x0<xlb)errmsg='x0 must be in the range xlb <= x0.'; return;elseif any(xstatus==1 & (~isfinite(xlb) | xlb~=round(xlb)))errmsg='xlb(i) must be an integer if x(i) is an integer variabele.'; return;end;xlb(find(xstatus==2))=x0(find(xstatus==2));xub=ones(size(x0));xub(find(xstatus==0))=inf;if nargin>4 & ~isempty(xu)if all(size(xu)<=size(x0))xub(1:size(xu,1))=xu;else errmsg='xub must be a column vector the same size as x0.'; return;end;end;if any(x0>xub)errmsg='x0 must be in the range x0 <=xub.'; return;elseif any(xstatus==1 & (~isfinite(xub) | xub~=round(xub)))errmsg='xub(i) must be an integer if x(i) is an integer variabale.'; return;end;xub(find(xstatus==2))=x0(find(xstatus==2));if nargin>5if ~isempty(A) & size(A,2)~=size(x0,1), errmsg='Matrix A not correct.'; return; end; else A=[]; end;if nargin>6if ~isempty(B) & any(size(B)~=[size(A,1) 1]), errmsg='Column vector B not correct.'; return; end;else B=[]; end;if isempty(A) & ~isempty(B), errmsg='A and B should only be nonempty together.'; return; end;if isempty(B) & ~isempty(A), B=zeros(size(A,1),1); end;if nargin>7 & ~isempty(Aeq)if size(Aeq,2)~=size(x0,1), errmsg='Matrix Aeq not correct.'; return; end;else Aeq=[]; end;if nargin>8if ~isempty(Beq) & any(size(Beq)~=[size(Aeq,1) 1]), errmsg='Column vector Beq not correct.'; return; end;else Beq=[]; end;if isempty(Aeq) & ~isempty(Beq), errmsg='Aeq and Beq should only be nonempty together'; return; end;if isempty(Beq) & ~isempty(Aeq), Beq=zeros(size(Aeq,1),1); end;if nargin<10, nonlcon=''; end;settings = [0 0 0];if nargin>10 & ~isempty(setts)if all(size(setts)<=size(settings))settings(setts~=0)=setts(setts~=0);else errmsg='settings should be a row vector of length 3.'; return; end;end;if nargin<12, options1=[]; end;options1=optimset(optimset('fmincon'),options1);if nargin<13, options2=[]; end;options2=optimset(optimset('fmincon'),options2);if nargin<14, maxSQPiter=1000;elseif isnumeric(maxSQPit) & all(size(maxSQPit))==1 & maxSQPit>0 &round(maxSQPit)==maxSQPitmaxSQPiter=maxSQPit;else errmsg='maxSQPiter must be an integer >0'; return; end;eval(['z=',fun,'(x0,varargin{:});'],'errmsg=''fun caused error.''; return;');if ~isempty(nonlcon)eval(['[C, Ceq]=',nonlcon,'(x0,varargin{:});'],'errmsg=''nonlcon caused error.''; return;');if size(C,2)>1 | size(Ceq,2)>1, errmsg='C en Ceq must be column vectors.'; return; end;end;% STEP 1 INITIALISATIONcurrentwarningstate=warning;warning off;tic;lx = size(x0,1);z_incumbent=inf;x_incumbent=inf*ones(size(x0));I = ceil(sum(log2(xub(find(xstatus==1))-xlb(find(xstatus==1))+1))+size(find(xstatus==1),1)+1);stackx0=zeros(lx,I);stackx0(:,1)=x0;stackxlb=zeros(lx,I);stackxlb(:,1)=xlb;stackxub=zeros(lx,I);stackxub(:,1)=xub;stacksize=1;xchoice=zeros(size(x0));if ~isempty(Aeq)j=0;for i=1:size(Aeq,1)if Beq(i)==1 & all(Aeq(i,:)==0 | Aeq(i,:)==1)J=find(Aeq(i,:)==1);if all(xstatus(J)~=0 & xchoice(J)==0 & xlb(J)==0 & xub(J)==1) if all(xstatus(J)~=2) | all(x0(J(find(xstatus(J)==2)))==0)j=j+1;xchoice(J)=j;if sum(x0(J))==0, errmsg='x0 not correct.'; return; end;end;end;end;end;end;errx=optimget(options2,'TolX');errcon=optimget(options2,'TolCon');fail=0;c=0;% STEP 2 TERMINIATIONwhile stacksize>0c=c+1;% STEP 3 LOADING OF CSPx0=stackx0(:,stacksize);xlb=stackxlb(:,stacksize);xub=stackxub(:,stacksize);x0(find(x0<xlb))=xlb(find(x0<xlb));x0(find(x0>xub))=xub(find(x0>xub));stacksize=stacksize-1;% STEP 4 RELAXATION% PHASE 1con=BNBCON(x0,A,B,Aeq,Beq,xlb,xub,nonlcon,varargin{:});if abs(con)>errcon & settings(1)~=0[x1 dummyfeasflag]=fmincon('0',x0,A,B,Aeq,Beq,xlb,xub,nonlcon,options1,varargin{:});if settings(3) & feasflag==0con=BNBCON(x1,A,B,Aeq,Beq,xlb,xub,nonlcon,varargin{:});if con<errcon, feasflag=1; end;end;else x1=x0; feasflag=1; end;% PHASE 2if feasflag>0[x z convflag]=fmincon(fun,x1,A,B,Aeq,Beq,xlb,xub,nonlcon,options2,varargin{:});if settings(3) & convflag==0con=BNBCON(x,A,B,Aeq,Beq,xlb,xub,nonlcon,varargin{:});if con<errcon, convflag=1; end;end;else convflag=feasflag; end;% STEP 5 FATHOMINGK = find(xstatus==1 & xlb~=xub);separation=1;if convflag<0 | (convflag==0 & settings(2))% FC 1separation=0;elseif z>=z_incumbent & convflag>0% FC 2separation=0;elseif all(abs(round(x(K))-x(K))<errx) & convflag>0% FC 3z_incumbent = z;x_incumbent = x;separation = 0;end;% STEP 6 SELECTIONif separation == 1 & ~isempty(K)dzsep=-1;for i=1:size(K,1)dxsepc = abs(round(x(K(i)))-x(K(i)));if dxsepc>=errx | convflag==0xsepc = x; xsepc(K(i))=round(x(K(i)));dzsepc = abs(feval(fun,xsepc,varargin{:})-z);if dzsepc>dzsepdzsep=dzsepc;ixsep=K(i);end;end;end;% STEP 7 SEPARATIONif xchoice(ixsep)==0% XCHOICE==0branch=1;domain=[xlb(ixsep) xub(ixsep)];while branch==1xboundary=(domain(1)+domain(2))/2;if x(ixsep)<xboundarydomainA=[domain(1) floor(xboundary)];domainB=[floor(xboundary+1) domain(2)];elsedomainA=[floor(xboundary+1) domain(2)];domainB=[domain(1) floor(xboundary)];end;stacksize=stacksize+1;stackx0(:,stacksize)=x;stackxlb(:,stacksize)=xlb;stackxlb(ixsep,stacksize)=domainB(1);stackxub(:,stacksize)=xub;stackxub(ixsep,stacksize)=domainB(2);if domainA(1)==domainA(2)stacksize=stacksize+1;stackx0(:,stacksize)=x;stackxlb(:,stacksize)=xlb;stackxlb(ixsep,stacksize)=domainA(1);stackxub(:,stacksize)=xub;stackxub(ixsep,stacksize)=domainA(2);branch=0;elsedomain=domainA;branch=1;end;end;else% XCHOICE~=0L=find(xchoice==xchoice(ixsep));M=intersect(K,L);[dummy,N]=sort(x(M));part1=M(N(1:floor(size(N)/2))); part2=M(N(floor(size(N)/2)+1:size(N))); stacksize=stacksize+1;stackx0(:,stacksize)=x;O = (1-sum(stackx0(part1,stacksize)))/size(part1,1);stackx0(part1,stacksize)=stackx0(part1,stacksize)+O;stackxlb(:,stacksize)=xlb;stackxub(:,stacksize)=xub;stackxub(part2,stacksize)=0;stacksize=stacksize+1;stackx0(:,stacksize)=x;O = (1-sum(stackx0(part2,stacksize)))/size(part2,1);stackx0(part2,stacksize)=stackx0(part2,stacksize)+O;stackxlb(:,stacksize)=xlb;stackxub(:,stacksize)=xub;stackxub(part1,stacksize)=0;if size(part2,1)==1, stackxlb(part2,stacksize)=1; end;end;elseif separation==1 & isempty(K)fail=fail+1;end;end;% STEP 8 OUTPUTt=toc;Z = z_incumbent;X = x_incumbent;errmsg='';eval(['warning ',currentwarningstate]);function CON=BNBCON(x,A,B,Aeq,Beq,xlb,xub,nonlcon,varargin); if isempty(A), CON1=[]; else CON1 = max(A*x-B,0); end;if isempty(Aeq), CON2=[]; else CON2 = abs(Aeq*x-Beq); end; CON3 = max(xlb-x,0);CON4 = max(x-xub,0);if isempty(nonlcon)CON5=[]; CON6=[];else[C Ceq]=feval(nonlcon,x,varargin{:});CON5 = max(C,0);CON6 = abs(Ceq);end;CON = max([CON1; CON2; CON3; CON4; CON5; CON6]);。

NSGA-II 算法实例目前的多目标优化算法有很多, Kalyanmoy Deb 的带精英策略的快速非支配排序遗传算法(NSGA-II) 无疑是其中应用最为广泛也是最为成功的一种。

本文用的算法是MATLAB 自带的函数gamultiobj ，该函数是基于NSGA-II 改进的一种多目标优化算法。

一、 数值例子多目标优化问题424221*********422421221211212min (,)10min (,)55..55f x x x x x x x x x f x x x x x x x x x s t x =-++-=-++-≤≤⎧⎨-≤≤⎩二、 Matlab 文件1． 适应值函数m 文件： function y=f(x) y(1)=x(1)^4-10*x(1)^2+x(1)*x(2)+x(2)^4-x(1)^2*x(2)^2;y(2)=x(2)^4-x(1)^2*x(2)^2+x(1)^4+x(1)*x(2);2． 调用gamultiobj 函数，及参数设置：clearclcfitnessfcn=@f; %适应度函数句柄nvars=2; %变量个数lb=[-5,-5]; %下限ub=[5,5]; %上限A=[];b=[]; %线性不等式约束Aeq=[];beq=[]; %线性等式约束options=gaoptimset('paretoFraction',,'populationsize',100,'ge nerations',200,'stallGenLimit',200,'TolFun',1e-100,'PlotFc ns',@gaplotpareto);% 最优个体系数paretoFraction 为；种群大小populationsize 为100，最大进化代数generations 为200，% 停止代数stallGenLimit 为200， 适应度函数偏差TolFun 设为1e-100，函数gaplotpareto ：绘制Pareto 前端[x,fval]=gamultiobj(fitnessfcn,nvars,A,b,Aeq,beq,lb,ub,options)3. 计算结果-40-35-30-25-20-15-10-5-505101520253035Objective 1O b j e c t i v e 2Pareto front图1. 实例1对应的Pareto 前沿图从图1可以看出Pareto 前分布较均匀，多样性较好。