MATLAB编程：非线性规划

合集下载

matlab解决非线性规划问题（凸优化问题）

matlab解决⾮线性规划问题（凸优化问题）当⽬标函数含有⾮线性函数或者含有⾮线性约束的时候该规划问题变为⾮线性规划问题，⾮线性规划问题的最优解不⼀定在定义域的边界，可能在定义域内部，这点与线性规划不同；例如：编写⽬标函数，定义放在⼀个m⽂件中；编写⾮线性约束条件函数矩阵，放在另⼀个m⽂件中；function f = optf(x);f = sum(x.^2)+8;function [g, h] = limf(x);g = [-x(1)^2+x(2)-x(3)^2x(1)+x(2)^2+x(3)^3-20]; %⾮线性不等式约束h = [-x(1)-x(2)^2+2x(2)+2*x(3)^2-3]; %⾮线性等式约束options = optimset('largescale','off');[x y] = fmincon('optf',rand(3,1),[],[],[],[],zeros(3,1),[],...'limf',options)输出为：最速下降法（求最⼩值）：代码如下：function [f df] = detaf(x);f = x(1)^2+25*x(2)^2;df = [2*x(1)50*x(2)];clc,clear;x = [2;2];[f0 g] = detaf(x);while norm(g)>1e-6 %收敛条件为⼀阶导数趋近于0p = -g/norm(g);t = 1.0; %设置初始步长为1个单位f = detaf(x+t*p);while f>f0t = t/2;f = detaf(x+t*p);end %这⼀步很重要，为了保证最后收敛，保持f序列为⼀个单调递减的序列，否则很有可能在极值点两端来回震荡，最后⽆法收敛到最优值。

x = x+t*p;[f0,g] = detaf(x);endx,f0所得到的最优值为近似解。

基于matlab求解非线性规划问题

计算结果为：
x =[ 3.0000 5.0000 0.0000 7.0000 0.0000 1.0000 0.0000 0.0000 4.0000 0.0000 6.0000 10.0000]’ fval = 136.2275
即由料场 A、 B 向 6 个工地运料方案为： 1 料场 A 料场 B 3 0 2 5 0 3 0 4 4 7 0 5 0 6 6 1 10
（二）使用临时料场的情形
使用两个临时料场A(5,1)，B(2,7).求从料场j向工地i的运送量为Xij，在各工地用量必须满足和各料场运送量不超过日储量的条件下，使总的吨千米数最小，这是线性规划问题. 线性规划模型为：
min f
aa ( i , j ) X
j 1 i 1
2
6
ij
2、先建立M-文件 fun3.m: function f=fun3(x); f=-x(1)-2*x(2)+(1/2)*x(1)^2+(1/2)*x(2)^2
3、再建立主程序youh2.m： x0=[1;1]; A=[2 3 ;1 4]; b=[6;5]; Aeq=[];beq=[]; VLB=[0;0]; VUB=[]; [x,fval]=fmincon('fun3',x0,A,b,Aeq,beq,VLB,VUB) 4、运算结果为： x = 0.7647 1.0588 fval = -2.0294
工地位置（ a， b）及水泥日用量 d 1 a b d 1 .2 5 1 .2 5 3 2 8 .7 5 0 .7 5 5 3 0 .5 4 .7 5 4 4 5 .7 5 5 7 5 3 6 .5 6 6 7 .2 5 7 .2 5 11

非线性规划的MATLAB解法

特点
非线性规划问题通常具有多个局部最优解，解的稳定性与初始条件有关，需要使用特定的算法来找到全局最优解。
非线性规划的应用场景
数据拟合、模型选择、参数估计等。
生产计划、物流优化、设备布局等。
投资组合优化、风险管理、资本预算等。
金融
工业
科研
非线性规划的挑战与解决方法
挑战
非线性规划问题可能存在多个局部最优解，且解的稳定性与初始条件密切相关，需要使用特定的算法来找到全局最优解。
共轭梯度法
总结词
灵活、适用于大型问题、迭代方向交替
详细描述
共轭梯度法结合了梯度下降法和牛顿法的思想，通过迭代更新搜索方向，交替使用梯度和共轭方向进行搜索。该方法适用于大型非线性规划问题，具有较好的灵活性和收敛性。
04
非线性规划问题的约束处理
不等式约束处理
处理方式
在Matlab中，可以使用 `fmincon`函数来求解非线性规划问题，该函数可以处理不等式约束。
要点二
详细描述
这类问题需要同时考虑多个目标函数，每个目标函数可能有不同的优先级和权重。在Matlab中，可以使用 `gamultiobj`函数来求解这类问题。该函数可以处理具有多个目标函数的约束优化问题，并允许用户指定每个目标函数的权重和优先级。
谢谢观看
具体操作
将等式约束条件表示为线性方程组，并使用`Aeq`参数指定系数矩阵，使用`beq`参数指定常数向量。
注意事项
等式约束条件需要在可行域内满足，否则会导致求解失败。
边界约束处理
处理方式
边界约束可以通过在目标函数中添加惩罚项来处理，或者使用专门的优化算法来处理。
具体操作
在目标函数中添加惩罚项时，需要在目标函数中添加一个与边界约束相关的项，并调整其权重以控制边界约束的重要性。

MATLAB优化工具箱--线性规划-非线性规划

数学实验
linprog输入参数说明： f, A, b, Aeq, beq lb,ub 边界设置说明：如果x(i)无边界，则 lb(i) = -inf, ub(i) = inf
6
linprog 输出参数说明： x 决策变量取值 fval 目标函数最优值
exitflag > 0 成功找到最优解 0 达到最大迭代次数也没有找到最优解 < 0 该线性规划问题不可行或者linprog计
10
fmincon函数求解形如下面的有约束非线性规划模型
一般形式：
min f ( X ) s.t. AX b
Aeq X beq l X u c(X ) 0 ceq ( X ) 0
Matlab求解有约束非线性最小化 1.约束中可以有等式约束 2.可以含线性、非线性约束均可
数学实验
输入参数语法：
例子：某农场种植两种作物A、B，需要甲、乙两种化肥。种植每亩作物A和作物B分别需用的化肥数，可得利润及农场现有化
肥数量如下表所示：
问在现有条件下，如何安排种植，才能使利润最大？
作物
每亩所需化肥（百公斤）
现有化肥
AB （百公斤）
化肥
甲
23
100） 6 4
数学实验
例题建模
[x,fval,exitflag,output,lambda]=fmincon(fun,x0,...)
数学实验
输入参数的几点说明
模型中如果没有A,b,Aeq,beq,lb,ub的限制，则以空矩阵[ ]作为参数传入； nonlcon：如果包含非线性等式或不等式约束，则将这些函数
编写为一个Matlab函数， nonlcon就是定义这些函数的程序文件名；
3

多元非线性规划Matlab,非线性规划MATLAB代码

多元⾮线性规划Matlab,⾮线性规划MATLAB代码下⾯是三个⾮线性规划领域的算法。

课堂上给予了详细的讲解，在实践环节让学⽣编程实现，从⽽可以实验复杂⼀些的例⼦，加深对算法的理解。

下⾯共有四个程序grad，simplelinesearch，bfgs和phr，全部使⽤MATLAB语⾔编写。

这些代码远未完善，可修改余地很⼤，仅供教学之⽤。

function gradf=grad(hfun,x)%GRAD 数值法求函数在给定点处的导数值(⼀元函数)或梯度(多元函数)% gradf = grad(hfun,x0) hfun是函数句柄或内联函数，x0是⼀定点或⼀批点(按列)；返回% 值gradf是函数在该点处的导数或梯度。

% 要求函数能对成批的点求函数值。

⽐如：feval(hfun,X)返回⼀个与X同列的⾏向量，对应于以% X每⼀列作为函数⾃变量⽽求得的函数值。

%% Reference: 《最优化计算原理与算法程序设计》，粟塔⼭等编著，国防科技⼤学出版社% (湖南)， 2001.%% $Author: WBC $ $Date: 2003/10/25 $n = length(x);h=1e-3; % 数值法求梯度的步长w1=zeros(n,1);%h/2w2=w1; %-h/2w3=w1; %hw4=w1; %-hfor i=1:nx(i)=x(i)+h/2;w1(i)=hfun(x);x(i)=x(i)-h;w2(i)=hfun(x);x(i)=x(i)-h/2;w4(i)=hfun(x);x(i)=x(i)+2*h;w3(i)=hfun(x);x(i)=x(i)-h;endgradf=(8*(w1-w2)-(w3-w4))/(6*h);function [fv,x,lambda,exitflag]=simplelinesearch(hf,rho,l,u,lambda0,fv0,x0,g,d)%LINESEARCH 简单线搜索%输⼊参数：% hf -- 函数句柄，⽬标函数% rho -- 实标量，简单线搜索的参数，介于0到0.5之间的数% l -- 实标量，简单线搜索的参数，介于0到1之间的数% u -- 实标量，简单线搜索的参数，介于0到1之间的数，满⾜u>l% alpha0 -- 实标量，步长的初始值% fv0 -- 实标量，⼀维搜索的初始⽬标函数值，即 hf(x0+alpha_0*d)% x0 -- 实列向量，当前点% g -- 实列向量，函数hf在当前点处的梯度% d -- 实列向量，函数hf在当前点处的搜索⽅向%输出参数：% fv -- 实标量，⼀维搜索完成后的⽬标函数值，即 hf(x0+alpha*d)% x -- 实列向量，下⼀个点% lambda -- 实标量，可接受步长% exitflag -- 整型标量，等于0表⽰线搜索成功，等于-1表⽰线搜索失败(内部迭代次数⼤于iterMax) %参考⽂献：倪勤，最优化⽅法与程序设计，科学出版社% Date: 2009/12/20lambda = lambda0;x = x0;i = 0;imax = 30; % 最⼤迭代次数，⽤户可以修改% 主循环gd = dot(g, d);while i <= imaxfv = hf(x + lambda*d);i = i + 1;if fv < fv0 + lambda * rho * gd;x = x + lambda*d;exitflag = 0;return;endlambda_bar = -gd*lambda^2*0.5/(fv-fv0-lambda*gd);lambda = min(lambda_bar, u*lambda);endexitflag=0;if i >= imax && fv >= fv0fv = fv0;x = x0;lambda = 0;exitflag = -1;endfunction [fv, x, exitflag] = bfgs(hf, x0, epsi)%BFGS ⽆约束问题的BFGS算法%输⼊参数：% hf -- 函数句柄，⽬标函数% x0 -- 实列向量，初始点% epsi -- 实标量，终⽌误差%输出参数：% fv -- 实标量，⼀维搜索完成后的⽬标函数值，即 hf(x0+alpha*d)% x -- 实列向量，下⼀个点% exitflag -- 整型标量，等于0表⽰成功，等于-1表⽰失败(迭代次数⼤于iter_max) %参考⽂献：倪勤，最优化⽅法与程序设计，科学出版社% Date: 2009/12/20%%初始化%k = 0;rho = 0.01;l = 0.15;u = 0.85;x =x0;fv = hf(x);n = length(x);H = eye(n);iter_max = 100;%检查终⽌条件%g = grad(hf, x);while norm(g) > epsi && k<= iter_maxd = -H*g;lambda0 = 1.0;%%做线搜索，如果成功，则返回更新当前点%[fv, x, lambda, exitflag] = simplelinesearch(hf, rho, l, u, lambda0, fv, x, g, d); if exitflag == -1 %重开始H = eye(n);else%%更新H%g_old = g;g = grad(hf, x);p = lambda*d;q = g- g_old;Hq = H*q;pq = p'*q;qHq = q'*Hq;v = sqrt(qHq) * (p/pq-Hq/qHq);H = H + p*p'/pq - Hq*Hq'/qHq + v*v';endk = k+1;endexitflag = 0;if k > iter_maxexitflag = -1;endfunction [fv, x, exitflag] = phr(hf, cf, x0)%输⼊参数：% hf -- 函数句柄，⽬标函数% cf -- 函数句柄，约束条件，包含等式约束和不等式约束% x0 -- 实列向量，初始点%输出参数：% fv -- 实标量，⼀维搜索完成后的⽬标函数值，即 hf(x0+alpha*d)% x -- 实列向量，下⼀个点% exitflag -- 整型标量，等于0表⽰成功，等于-1表⽰失败(迭代次数⼤于iter_max) %参考⽂献：倪勤，最优化⽅法与程序设计，科学出版社% Date: 2009/12/20%%初始化%epsi = 1.0e-4;k = 0;sigma = 0.8;c = 1.5;theta = 0.8;x = x0;[ce, ci] = cf(x);l = length(ce);li = length(ci);lambda = ones(l+li, 1) * 0.1;iter_max = 100;phi = 0;if lphi = phi + ce'*ce;endif liphi = phi + sum(min(ci,lambda(l+1:end)/sigma).^2);endwhile phi > epsi && k <= iter_max%%hmf = @(x) mfun(x, hf, cf,lambda, sigma);[fv, x] = bfgs(hmf, x, epsi);[ce, ci] = cf(x);phi_old = phi;phi = 0;if lphi = phi + ce'*ce;endif liphi = phi + sum(min(ci,lambda(l+1:end)/sigma).^2); endif phi > epsi%%更新罚因⼦%if k >= 2 && phi/phi_old > thetasigma = c * sigma;end%%更新乘⼦%if llambda(1:l) = lambda(1:l) - sigma*ce;endif lilambda(l+1:end) = max(0, lambda(l+1:end) - sigma*ci); endendk = k+1;endexitflag = 0;fv = hf(x);exitflag = -1;end%%乘⼦罚函数%function fv = mfun(x, hf, cf, lambda, sigma)[ce, ci] = cf(x);l = length(ce);li = length(ci);fv = 0;fv = fv + hf(x);if lfv = fv - lambda(1:l)'*ce + 0.5*sigma*ce'*ce;endif lifv = fv + 0.5/sigma*sum(max(0,lambda(l+1:end) - sigma*ci).^2 - lambda(l+1:end).^2); end这⾥是演⽰代码：%% bfgs演⽰%教材P328.1-3hf1 = @(x) 100 * (x(2) - x(1).^2).^2 + (1 - x(1)).^2; %banana函数hf2 = @(x) (6 + x(1) + x(2)).^2 + (2 - 3*x(1) - 3*x(2) - x(1)*x(2)).^2;hf3 = @(x) x(1).^2 - 2*x(1)*x(2) + 4*x(2).^2 + x(1) - 3*x(2);[fv,x,exitflag]=bfgs(hf1,[0;0],0.001);fvxexitflag[fv,x,exitflag]=bfgs(hf2,[4;6],0.001);fvxexitflag[fv,x,exitflag]=bfgs(hf3,[1;1],0.001);fvexitflag%% phr算法%教材P414.5hf1 = @(x) x(1).^2 + x(2).^2;cf1 = @(x) deal([], x(1) - 1);hf2 = @(x) x(1) + (x(2) + 1).^2/3;cf2 = @(x) deal([],[x(1); x(2) - 1]);%教材P392.2hf3 = @(x) x(1).^2 + x(1).*x(2) + 2*x(2).^2 - 6*x(1) - 2*x(2) - 12*x(3);cf3 = @(x) deal(x(1)+x(2)+x(3)-2,...[x(1) - 2*x(2) + 3; x(1); x(2); x(3)]);%其它例⼦hf4 = @(x) 6*x(2)*x(5) + 7*x(1)*x(3) + 3*x(2)^2;cf4 = @(x) deal([3*x(2)^2*x(5) + 3*x(1)^2*x(3) - 20.875;x(1) - 0.3*x(2)],...[-x(1) + 0.2*x(2)*x(5) + 71-0.9*x(3) + x(4)^2 + 67x(3)x(5) - 1-x(3) + 20x(4) - 0.1*x(5)-x(4) + 0.5*x(5)x(3) - 0.9*x(5)]);hf5 = @(x) exp(x(1)) * (4*x(1)^2 + 2*x(2)^2 + 4*x(1)*x(2) + 2*x(2) + 1);cf5 = @(x) deal([], [x(1) + x(2) - x(1)*x(2) - 1.5; x(1)*x(2) + 10]);[fv, x, exitflag] = phr(hf1, cf1, [3;2]); fv x exitflag [fv, x, exitflag] = phr(hf2, cf2, [3;2]); fv x exitflag [fv, x, exitflag] =phr(hf3, cf3, [1;1;0]); fv x exitflag [fv, x, exitflag] = phr(hf4, cf4, [1; 4; 5; 2; 5]); fv x exitflag [fv, x, exitflag] = phr(hf5, cf5, [-1; 1]); fv x exitflag。

MATLAB线性规划非线性规划

解编写M文件xxgh1.m如下： c=[-0.4 -0.28 -0.32 -0.72 -0.64 -0.6]; A=[0.01 0.01 0.01 0.03 0.03 0.03;0.02 0 0 0.05 0 0;0 0.02 0 0 0.05 0;0 0 0.03 0 0 0.08]; b=[850;700;100;900];
返回
解答
线性规划模型的一般形式
目标函数和所有的约束条件都是设计变量的线性函数.
min u ci xi
i 1
n
n aik xk bi , i 1, 2,..., n. s.t. k 1 x 0, i 1, 2,..., n. i
矩阵形式： min u cx Ax b s.t. vlb x vub
Aeq=[]; beq=[];
vlb=[0;0;0;0;0;0]; vub=[];
To MATLAB (xxgh1)
[x,fval]=linprog(c,A,b,Aeq,beq,vlb,vub)
例2
min z 6x1 3x2 4x3 s.t. x1 x2 x3 120 x1 30 0 x2 50 x3 20
编写M文件xxgh3.m如下: f = [13 9 10 11 12 8]; A = [0.4 1.1 1 0 0 0 0 0 0 0.5 1.2 1.3]; b = [800; 900]; Aeq=[1 0 0 1 0 0 0 1 0 0 1 0 0 0 1 0 0 1]; To MATLAB (xxgh3) beq=[400 600 500]; vlb = zeros(6,1); vub=[]; [x,fval] = linprog(f,A,b,Aeq,beq,vlb,vub)

多目标非线性规划程序Matlab完整版

多目标非线性规划程序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]);。

非线性规划及matlab应用

非线性规划及matlab 应用目录1.概念 ............................................................................................................................................... 1 2.二次规划........................................................................................................................................ 1 3.一般非线性规划 ............................................................................................................................ 2 4. 案例：供应与选址 . (4)1.概念定义：如果目标函数或约束条件中至少有一个是非线性函数时的最优化问题就叫做非线性规划问题．其它情况: 求目标函数的最大值或约束条件为小于等于零的情况，都可通过取其相反数化为上述一般形式．2.二次规划用MATLAB 软件求解,其输入格式如下: 1. x=quadprog(H,C,A,b); 2. x=quadprog(H,C,A,b,Aeq,beq); 3. x=quadprog(H,C,A,b,Aeq,beq,VLB,VUB); 4. x=quadprog(H,C,A,b, Aeq,beq ,VLB,VUB,X0); 5. x=quadprog(H,C,A,b, Aeq,beq ,VLB,VUB,X0,options); 6. [x,fval]=quadprog(...); 7. [x,fval,exitflag]=quaprog(...); 实例1：2212121122121212min (,)2622..2220,0=--+-++≤-+≤≥≥f x x x x x x x x s tx x x x x x写成标准形式标准型为： Min Z= 21X T HX+c T Xs.t. AX<=b beq X Aeq =⋅ VLB ≤X ≤VUB111222 1 -12min (,) 1 26Tx x z x x x x ⎛⎫⎛⎫⎛⎫⎛⎫-=+ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪--⎝⎭⎝⎭⎝⎭⎝⎭1212 1 121 2200x x x x ⎛⎫⎛⎫⎛⎫≤ ⎪ ⎪ ⎪ ⎪-⎝⎭⎝⎭⎝⎭⎛⎫⎛⎫≤ ⎪ ⎪ ⎪⎝⎭⎝⎭Matlab 命令H=[1 -1; -1 2]; c=[-2 ;-6]; A=[1 1; -1 2]; b=[2;2]; Aeq=[]; beq=[]; VLB=[0;0]; VUB=[];[x,z]=quadprog(H,c,A,b,Aeq,beq,VLB,VUB)运算结果为：x =0.6667 1.3333 z = -8.22223.一般非线性规划标准型为：min ()..()0()0≤≤=≤≤F X s t AX bG X Ceq X VLB X VUB其中X 为n 维变元向量，G(X)与Ceq(X)均为非线性函数组成的向量，其它变量的含义与线性规划、二次规划中相同.非线性规划求解的函数是fmincon,命令的基本格式如下： ● x=fmincon(‘fun’,X0,A,b)● x=fmincon(‘fun’,X0,A,b,Aeq,beq)● x=fmincon(‘fun’,X0,A,b, Aeq,beq,VLB,VUB)● x=fmincon(‘fun’,X0,A,b,Aeq,beq,VLB,VUB,’nonlcon’)● x=fmincon(‘fun’,X0,A,b,Aeq,beq,VLB,VUB,’nonlcon’,options) ● [x,fval]= fmincon(...)● [x,fval,exitflag]= fmincon(...)● [x,fval,exitflag,output]= fmincon(...) 1.fun 为目标函数2.x0为初始值3.A 是不等式约束AX<=b 的系数矩阵4.b 是不等式约束AX<=b 的常数项4.Aeq 是等式约束AeqX=beq 的系数矩阵，5.beq 是等式约束AeqX=beq 的常数项，6.lb 是X 的下限，7.ub 是X 的上限，8.nonlcon 为非线性不等式约束 9.option 为设置fmincon 的参数注意：fmincon 函数提供了大型优化算法和中型优化算法。

多目标非线性规划程序(Matlab)

function [errmsg,Z,X,t,c,fail] =BNB18(fun,x0,xstat,xl,xu,A,B,Aeq,Beq,nonlcon,setts,options1,options2,ma xSQPit,varargin);%·ÇÏßÐÔÕûÊý¹æ»®Ä£ÐÍÇó½â·ÖÖ§¶¨½çµü´úËã·¨¡£ÔÚMATLAB5.3ÖÐÊ¹ÓÃ£¬ÐèOptimizat ion toolbox 2.0Ö§³Ö?% Minimize F(x)%subject to: xlb <= x <=xub% A*x <= B% Aeq*x=Beq% C(x)<=0% Ceq(x)=0%% x(i)¿ÉÎªÁ¬Ðø±äÁ¿£¬ÕûÊý£¬»ò¹Ì¶¨Öµ% Ê¹ÓÃ¸ñÊ½%[errmsg,Z,X]=BNB18('fun',x0,xstat,xl,xu,A,B,Aeq,Beq,'nonlcon',setts)%fun£º MÎÄ¼þÃû£¬±íÊ¾×îÐ¡»¯Ä¿±êº¯Êýf=fun(x)%x0: ÁÐÏòÁ¿£¬±íÊ¾±äÁ¿³õÖµ%xstat£º ÁÐÏòÁ¿£¬xstat(i)=0±íÊ¾x(i)ÎªÁ¬Ðø±äÁ¿£¬1±íÊ¾ÕûÊý£¬2±íÊ¾¹Ì¶¨Öµ%xl£º ÁÐÏòÁ¿£¬±íÊ¾±äÁ¿ÏÂ½ç%xu: ÁÐÏòÁ¿£¬±íÊ¾±äÁ¿ÉÏ½ç%A: ¾ØÕó, ±íÊ¾ÏßÐÔ²»µÈÊ½Ô¼ÊøÏµÊý%B: ÁÐÏòÁ¿, ±íÊ¾ÏßÐÔ²»µÈÊ½Ô¼ÊøÉÏ½ç%Aeq: ¾ØÕó, ±íÊ¾ÏßÐÔµÈÊ½Ô¼ÊøÏµÊý%Beg: ÁÐÏòÁ¿, ±íÊ¾ÏßÐÔ²»µÈÊ½Ô¼ÊøÓÒ¶ËÖµ%nonlcon:MÎÄ¼þÃû£¬±íÊ¾·ÇÏßÐÔÔ¼Êøº¯Êý[C,Ceq]=nonlin(x),ÆäÖÐC(x)Îª²»µÈÊ½Ô¼Êø,% Ceq(x)ÎªµÈÊ½Ô¼Êø%setts: Ëã·¨ÉèÖÃ%errmsq: ·µ»Ø´íÎóÌáÊ¾%Z: ·µ»ØÄ¿±êº¯Êý×îÐ¡Öµ%X: ·µ»Ø×îÓÅ½â%%ÀýÌâ% max x1*x2*x3% -x1+2*x2+2*x3>=0% x1+2*x2+2*x3<=72% 10<=x2<=20% x1-x2=10% ÏÈÐ´ Mº¯Êýdiscfun.m% function f=discfun(x)% f=-x(1)*x(2)*x(3);%Çó½â% 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,opti ons2,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 2.0 (R11) 09-Oct-1998% From this toolbox fmincon.m 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 erro r 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 bnb18.m).% 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.% E.C. Kuipers% e-mail E.C.Kuipers@cpedu.rug.nl% FI-Lab% Applied Physics% Rijksuniversiteit Groningen% To get rid of bugs and to stop fmincon from hanging make the following chances:%% In optim/private/nlconst.m ($Revision: 1.20 $ $Date: 1998/08/24 13:46:15 $):% Get EXITFLAG independent of verbosity.% After the lines: disp(' less than 2*options.TolFun 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/qpsub.m ($Revision: 1.21 $ $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 zconvflag]=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]);。

matlab非线性规划

i1
1
giX
m
m
其中 r 称 lngiX或r
i1
i1
gi1X为障碍 r为项障，碍因子
这样问1题）（就转化为求一值系问列题极：
XmDi0nIX,rk 得X（ k rk）
10
内点法的迭代步骤
( 1 ) 给定允许误差 0 ，取 r 1 0 , 0 1 ；
( 2 ) 求出约束集合 D 的一个内点 X 0 D 0 ，令 k 1 ；
7. [x,fval,exitflag]=quaprog(...);
8. [x,fval,exitflag,output]=quaprog(...);
15
例1 min f(x1,x2)=-2x1-6x2+x12-2x1x2+2x22 s.t. x1+x2≤2 -x1+2x2≤2 x1≥0, x2≥0
17
3. 建立主程序.非线性规划求解的函数是fmincon,命令的基本格式如下：
(1) x=fmincon(‘fun’,X0,A,b) (2) x=fmincon(‘fun’,X0,A,b,Aeq,beq) (3) x=fmincon(‘fun’,X0,A,b, Aeq,beq,VLB,VUB)
[3] fmincon函数可能会给出局部最优解，这与初值X0的选取有关。
19
例2 mfi n x 1 2 x 2 1 2 x 1 2 1 2 x 2 2 2x1+3x2 6
s.t x1+4x2 5 x1,x2 0
1、写成标准形式：
mfi n x12x21 2x1 21 2x2 2
2x13x260 s.t. x14x25 0

MATLAB非线性规划问题

MATLAB⾮线性规划问题⼀．⾮线性规划课题实例1 表⾯积为36平⽅⽶的最⼤长⽅体体积。

建⽴数学模型：设x、y、z分别为长⽅体的三个棱长，f为长⽅体体积。

max f = x y (36-2 x y)/2 (x+y)实例2 投资决策问题某公司准备⽤5000万元⽤于A、B两个项⽬的投资，设x1、x2分别表⽰配给项⽬A、B的投资。

预计项⽬A、B的年收益分别为20%和16%。

同时，投资后总的风险损失将随着总投资和单位投资的增加⽽增加，已知总的风险损失为2x12+x22+(x1+x2)2.问应如何分配资⾦，才能使期望的收益最⼤，同时使风险损失为最⼩。

建⽴数学模型：max f=20x1+16x2-λ[2x12+x22+(x1+x2)2]s.t x1+x2≤5000x 1≥0,x2≥0⽬标函数中的λ≥0是权重系数。

由以上实例去掉实际背景，其⽬标函数与约束条件⾄少有⼀处是⾮线性的，称其为⾮线性问题。

⾮线性规划问题可分为⽆约束问题和有约束问题。

实例1为⽆约束问题，实例2为有约束问题。

⼆．⽆约束⾮线性规划问题：求解⽆约束最优化问题的⽅法主要有两类：直接搜索法(Search method)和梯度法(Gradient method)，单变量⽤fminbnd,fminsearch,fminunc;多变量⽤fminsearch,fminnuc 1．fminunc函数调⽤格式：x=fminunc(fun,x0)x=fminunc(fun,x0,options)x=fminunc(fun,x0,options,P1,P2)[x,fval]=fminunc(…)[x,fval, exitflag]=fminunc(…)[x,fval, exitflag,output]=fminunc(…)[x,fval, exitflag,output,grad]=fminunc(…)[x,fval, exitflag,output,grad,hessian]=fminunc(…)说明：fun为需最⼩化的⽬标函数，x0为给定的搜索的初始点。

[线性整数规划matlab程序]matlab非线性规划

第1页共1页
? for i=1:NL
本文格式为 Word 版，下载可任意编辑，页眉双击删除即可。
={k};
if bs(x(i)-round(x(i)))1e-7?
b=BB{k};kk=FN(end);Fra bibliotekbrek;
FN=[FN,kk+1,kk+2];
end
temp_=zeros(1,NL);
? end
本文格式为 Word 版，下载可任意编辑，页眉双击删除即可。
第1页共1页
本文格式为 Word 版，下载可任意编辑，页眉双击删除即可。
[线性整数规划 mtlb 程序]mtlb 非线性规划
(1)={}; BB(1)={b};
线性整数规划的 mtlb 程序〔〕 %%本程序是用分枝定界法求解整数线性规划问题 %%问题的标准形式： %%? min c*x %% s.t. *x=b %% eq*x=beq %%? x 要求是整数 %%程序的编写者：欢乐逍遥侠 %%程序的最终修改时间是：20XX-11-24 function [y,fvl]=BrnchBound(c,,b,eq,beq) NL=length(c); UB=inf; LB=-inf; FN=[0];
k=0;? flg=0; while flg==0; [x,fvl,exitFlg]=linprog(c,,b,eq,beq); if (exitFlg == -2) | (fvl = UB)? ? FN(1)=[]; ? if isempty(FN)==1? flg=1; ? else k=FN(1); ={k}; b=BB{k}; ? end else
temp_(i)=1;
? if (i==NL) (bs(x(i)-round(x(i)))=1e-7)?

MATLAB非线性规划

MATLAB⾮线性规划MATLAB求解⾮线性规划可以使⽤ fmincon 函数，其数学模型可以写成如下形式：x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options)其中，fun是⽬标函数，x0是初始值，A，b 规定线性不等式约束条件，Aeq，beq 规定线性等式约束条件，lb 规定可⾏解的数值下限，ub规定可⾏解的数值上限。

nonlcon是包含⾮线性约束条件（C(x)，Ceq(x)）的函数。

使⽤options所指定的优化选项执⾏最⼩化。

例如，使⽤MATLAB计算如下⾮线性规划。

x0 = [0.5,0];A = [1,-2];b = 1;Aeq = [2,1];beq = 1;x = fmincon(fun,x0,A,b,Aeq,beq)带有边界约束的，例如：fun = @(x)1+x(1)/(1+x(2)) - 3*x(1)*x(2) + x(2)*(1+x(1));lb = [0,0];ub = [1,2];% 没有线性约束，因此将这些参数设置为 []。

A = [];b = [];Aeq = [];beq = [];% 尝试使⽤⼀个位于区域中部的初始点。

x0 = (lb + ub)/2;x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub)带有⾮线性约束的，例如：%% 主函数options=optimset('largescale','off');x = fmincon(@fun,rand(3,1),[],[],[],[],zeros(3,1),[], @nonlcon, options)%% ⽬标函数function f=fun(x)f=sum(x.^2)+8;end%% ⾮线性约束条件function [c,ceq]=nonlcon(x)c=[-x(1)^2+x(2)-x(3)^2x(1)+x(2)^2+x(3)^3-20]; %⾮线性不等式约束ceq=[-x(1)-x(2)^2+2x(2)+2*x(3)^2-3]; %⾮线性等式约束end特别注意：⽬标函数为最⼩化函数，fun是⼀个函数，fun接受向量或数组 x，并返回实数标量 f，即在 x 处计算的⽬标函数值。

MATLAB求解非线性规划

MATLAB基础及求解非线性规划
经济管理学院：李继红
1 Matlab简介及操作环境
1.1 Matlab简介
由美国 DOCTORCLEVER MOLER 于 1980 年开始研制并于1984年推出正式版本。以后陆续推出了4.0，4.2和5.0版本，1999年初推出了功能更为强大的5.3版本，2001年推出6.1版本，2003年推出了6.5版本。2004年6月正式推出7.0版本。 MATLAB是建立在 C 语言基础上的高级语言，并建立了自已独特的语言环境。
plot3(x,y,z,’s’)
plot3(x1,y1,z1,’s1’,x2,y2,z2,’s2’,·· ·) plot3例命：令t将=绘0制:二pi维/5图0形:1的0函*p数i;plot的特性扩展到三维空间。函数格式除p了lo包t括3(第s三in维(t的),信co息s(（t)比,t如);Z方向）之外，与二维函数
1、for循环
li1_6_1.m
• For循环可以按指定的次数重复执行一系列语句。For循环的常见结构：
• for 变量 =表达式循环语句
end
这里循环语句可以是一条或多条，并且可以是变量的函数
• 如：for k=初值: 增值: 终值
例：
％ mzmfor.m
clear all
N=input('请输入矩阵的维数 N:');
6.4 switch-case语句
一般switch-case语句格式为： switch num case n1 command case n2 command case n3 command . . .otherwise Command
li1_6_4.m
7 数据的可视化

用Matlab求解非线性规划

用Matlab 求解非线性规划1．无约束优化问题)(min x f n Rx ∈，其中向量x 的n 个分量i x 都是决策变量，称)(x f 目标函数。

用Matlab 求解：先建立函数文件mbhs.m ，内容是)(x f 的表达式；再回到Matlab 命令区输入决策变量初值数据x0，再命令[x,fmin]=fminunc(@mbhs,x0) 如：)32(m in 22212x x R x +∈的最优解是.)0,0(T x = 用Matlab 计算，函数文件为 function f=mbhs(x)f=2*x(1)^2+3*x(2)^2;再输入初值 x0=[1;1]; 并执行上述命令，结果输出为 x =？ fmin =？略。

2．约束优化问题.),,...,2,1(,0)(),,...,2,1(,0)(..)(min U x L m i x h p i x g t s x f i i Rx n ≤≤===≤∈其中：向量x 的n 个分量i x 都是决策变量，称)(x f 目标函数、)(x g i 等式约束函数、)(x h i 不等式约束函数、L 下界、U 上界。

用Matlab 求解：先把模型写成适用于Matlab 的标准形式.,0)(,0)(,,..)(min U x L x h x g beq x Aeq b Ax t s x f n Rx ≤≤=≤=≤∈ 约束条件中：把线性的式子提炼出来得前两个式子；后三个式子都是列向量。

（如：⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡===⨯⨯)()()([],[],,,11262x g x g x g beq Aeq b A p ）再建立两个函数文件：目标函数mbhs.m ；约束函数yshs.m再回到Matlab 命令区，输入各项数据及决策变量初值数据x0，执行命令[x,fmin]=fmincon(@mbhs,x0,A,b,Aeq,beq,L,U,@yshs)例：单位球1222≤++z y x 内，曲面xy y x z 1.05.022--+=的上方，平面008.0=-++z y x 之上（不是上面），满足上述三个条件的区域记为D ，求函数)1cos()sin(2-+-+-z e z y x e xy xyz 在D 上的最大值、最大值点。

工程优化 6-解非线性规划的MATLAB指令

练习：利用Matlab求解下列最优化问题
（ 1）
1 2 1 2 min f x1 2 x2 x1 x2 2 2 2x1+3x2 6 s.t x1+4x2 5 x1,x2 0
2 2 （2） min f ( x ) e x1 (4 x1 2 x2 4 x1 x2 2 x2 1)
(4) 若修改主程序gying2.m, 取初值为上面的计算结果:
x0=[ 3.0000 5.0000 0.0707 7.0000 0 0.9293 0 0 3.9293 0 6.0000 10.0707 6.3875 4.3943 5.7511 7.1867]’
得结果为:
x=[3.0000 5.0000 0.3094 7.0000 0.0108 0.6798 0 0 3.6906 0 5.9892 10.3202 5.5369 4.9194 5.8291 7.2852]’ fval =103.4760 exitflag = 1 总的吨千米数比上面结果略优. (5) 若再取刚得出的结果为初值, 却计算不出最优解.
(6) 若取初值为: x0=[3 5 4 7 1 0 0 0 0 0 5 11 5.6348 4.8687 7.2479 7.7499]', 则计算结果为: x=[3.0000 5.0000 4.0000 7.0000 1.0000 0 0 0 0 0 5.0000 11.0000 5.6959 4.9285 7.2500 7.7500]’ fval =89.8835 exitflag = 1 总的吨千米数89.8835比上面结果更好.
小结：用Matlab求解非线性规划问题，基本步骤：
1. 首先建立M文件fun.m,定义目标函数 f（x）: function f=fun(x); f= f（x）; 2.若约束条件中有非线性约束:g(x) 0

matlab中的非线性规划求解fmincon函数[修改版]

这个函数的基本形式为x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon，options)其中fun为你要求最小值的函数，可以单写一个文件设置函数，如以上给的例子中。

1.如果fun中有N个变量，如x y z, 或者是X1, X2,X3, 什么的，自己排个顺序，在fun 中统一都是用x(1),x(2）....x(n) 表示的。

2. x0, 表示初始的猜测值，大小要与变量数目相同3. A b 为线性不等约束，A*x <= b， A应为n*n阶矩阵，学过线性代数应不难写出A和b4 Aeq beq为线性相等约束，Aeq*x = beq。

Aeq beq同上可求5 lb ub为变量的上下边界，正负无穷用 -Inf和Inf表示， lb ub应为N阶数组6 nonlcon 为非线性约束，可分为两部分，非线性不等约束 c，非线性相等约束，ceq可按下面的例子设置function [c,ce] = nonlcon1(x)c = -x(1)+x(2)^2-4;ce = []; % no nonlinear equality constraints7，最后是options，可以用OPTIMSET函数设置，见上例具体可见OPTIMSET函数的帮助文件。

对于优化控制，MATLAB提供了18个参数，这些参数的具体意义为：options(1)-参数显示控制（默认值为0）。

等于1时显示一些结果。

    options(2)-优化点x的精度控制(默认值为1e-4)。

options = optimset('TolX',1e-8)    options(3)-优化函数F的精度控制(默认值为1e-4)。

options = optimset('TolFun',1e-10)    options(4)-违反约束的结束标准(默认值为1e-6)。

1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。

2 . 若约束条件中有非线性约束 : G ( X )0
或 Ceq(X)=0,
则建立 M文件 nonlcon.m定义函数 G(X)与 Ceq(X): function [G,Ceq]=nonlcon(X) G=... Ceq=...
17
3. 建立主程序.非线性规划求解的函数是fmincon,命令的基本格式如下： (1) x=fmincon(‘fun’,X0,A,b) (2) x=fmincon(‘fun’,X0,A,b,Aeq,beq) (3) x=fmincon(‘fun’,X0,A,b, Aeq,beq,VLB,VUB)
k
内点法的迭代步骤
(1) 给定允许误差
0 ，取 r1 0 , 0 1 ；
( 2) 求出约束集合 D 的一个内点
( 3) 以
0
X
0
D ，令k 1 ；
0
X
k 1
D 为初始点，求解 min I X , rk ，其中
考虑问题：
设集合 D
0
min f X s .t . g X 0 i
i 1, 2 ,..., m
0
(1)
X | g i X 0 , i 1, 2 , , m ， D 是可行域中
所有严格内点的集合。
构造障碍函数 I X , r ： I X , r f X r ln g i X
11
近似规划法
近似规划法的基本思想：将问题(3)中的目标函数 f X
和约束条件
giX
0 ( i 1,..., m);
h j X
0 ( j 1, , l )
近似为线性函数，并对变量的取值范围加以限制，从
而得到一个近似线性规划问题，再用单纯形法求解之，
把其符合原始条件的最优解作为(3)的解的近似．每得到一个近似解后，都从这点出发，重复以上步骤．这样，通过求解一系列线性规划问题，产生一个由线性规划最优解组成的序列，经验表明，这样的序列往往收敛于非线性规划问题的解。
j 1, , n
高，故需要对变量的取值范围加以限制，所增加的约束条件是：
xj x
k j

k j
求解该线性规划问题，得到最优
解X
k 1
；
( 4)
检验X
k 1
点对原约束是否可行。若X
k 1
对原约束可行，

X D
定义2
对于问题(1)，设 X * D，若存在 0 ,使得对一切，且 X X * ，都有 f X * f X ，则称X*是f(X)在D上的 X D 局部极小值点（局部最优解）．特别地当 X X *时，若 f X * f X ，则称X*是f(X)在D上的严格局部极小值点（严格局部最优解）．定义3 对于问题(1),设 X * D ,对任意的X D ,都有 f X * f X 则称X*是f(X)在D上的全局极小值点（全局最优解）．特别地当 * 时，若 * ，则称X*是f(X)在D上的严格全局极小值点 X X f X f X （严格全局最优解）．返回5
7
SUTM外点法
对一般的非线性规划： min f X

i 1,2,..., m; j 1, 2 ,..., l . (1)
gi X 0 s .t . h j X 0
可设： T X , M
f X M
min 0 , g i X 2
i 1
m
M
h
j 1
l
X 2 j
(2)
将问题（ 1）转化为无约束问题：
X E
min T X , M
n

(3)
其中T(X,M)称为罚函数，M称为罚因子，带M的项称为罚项，这里的罚函数只对不满足约束条件的点实行惩罚：当 D 时，满 X 足各 i X 0 , hi X 0 ，故罚项=0，不受惩罚．当 D 时， g X
(4) x=fmincon(‘fun’,X0,A,b,Aeq,beq,VLB,VUB,’nonlcon’) (5)x=fmincon(‘fun’,X0,A,b,Aeq,beq,VLB,VUB,’nonlcon’,options)
1 2
1
2 x 2
6 x 2
s.t.
1 1
1 x1 2 x 2
2 2
2、输入命令：
0 x1 0 x2
H=[1 -1; -1 2]; c=[-2 ;-6];A=[1 1; -1 2];b=[2;2]; Aeq=[];beq=[]; VLB=[0;0];VUB=[]; [x,z]=quadprog(H,c,A,b,Aeq,beq,VLB,VUB)
min f(x1,x2)=-2x1-6x2+x12-2x1x2+2x22 s.t. x1+x2≤2 -x1+2x2≤2 x1≥0, x2≥0 T 1 - 1 x1 2 x1 1、写成标准形式：min z ( x , x )
例1
0 X D 0
0
X D 的最优解，设为 X
r
*
k
X rk D ；
k
( 4) 检验是否满足足，停止迭代，返回 (3 ）．
m

i 1
ln g i X
k
或r
m

1 giX
k
i 1

，满
X
X ；否则取 rk 1 rk ，令 k k 1 ，
数学建模与数学实验
非线性规划
后勤工程学院数学教研室
1
实验目的
1、直观了解非线性规划的基本内容。
2、掌握用数学软件求解优化问题。
实验内容
1、非线性规划的基本理论。
2、用数学软件求解非线性规划。 3、钢管订购及运输优化模型 4、实验作业。
2
非线性规划
非线性规划的基本概念
*非线性规划的基本解法
返回3
非现性规划的基本概念定义如果目标函数或约束条件中至少有一个是非线性函数时的最优化问题就叫做非线性规划问题．
i 1 m

或
I (X ,r) f (X ) r r 为障碍因子
m
1 gi X
i 1

其中称 r ln g i X
i 1
m

或 r
m
1 gi X
i 1

为障碍项，
X D
min I X , rk
0
这样问题（ 1）就转化为求一系列极

值问题：
得 X （ rk ）
10
( 2)
在点
X
k
处，将 f X ， giX ,h j X 按泰勒级数展开并
取一阶近似，得到近似线性规划问题：
min f X f X

k
f X X
k T k k T i
X
k

k
giX gi X hjX hj X

min T X , M
n
T (X k ,M k ) ；
i m ，使 g i X
3、若存在 i 1
，则取Mk>M(M
k
m
k
k 1
M , 10
)
令k=k+1返回（2），否则，停止迭代．得最优解 X * 计算时也可将收敛性判别准则
gi X
3、运算结果为： x =0.6667 1.3333
MATLAB（youh1）
z = -8.2222
16
标准型为： min F(X) s.t AX<=b Ceq(X)=0
2、一般非线性规划
Aeq X beq
G ( X 0 )
VLB X
VUB
其中X为n维变元向量，G(X)与Ceq(X)均为非线性函数组成的向量，其它变量的含义与线性规划、二次规划中相同.用 Matlab求解上述问题，基本步骤分三步： 1. 首先建立M文件fun.m,定义目标函数F（X）: function f=fun(X); f=F(X);
14 返回
否则，令

k j
j 1, , n ， k = k + 1 ，返回步骤 ( 2 ) ．
1、二次规划
标准型为： Min Z= 1
2
T T
X HX+c X
s.t. AX<=b
Aeq X beq
V LB ≤ X ≤ V UB
用MATLAB软件求解,其输入格式如下: 1. 2. 3. 4. 5. 6. 7. 8. x=quadprog(H,C,A,b); x=quadprog(H,C,A,b,Aeq,beq); x=quadprog(H,C,A,b,Aeq,beq,VLB,VUB); x=quadprog(H,C,A,b, Aeq,beq ,VLB,VUB,X0); x=quadprog(H,C,A,b, Aeq,beq ,VLB,VUB,X0,options); [x,fval]=quaprog(...); [x,fval,exitflag]=quaprog(...); 15 [x,fval,exitflag,output]=quaprog(...);