注册 登录  
 加关注
   显示下一条  |  关闭
温馨提示!由于新浪微博认证机制调整,您的新浪微博帐号绑定已过期,请重新绑定!立即重新绑定新浪微博》  |  关闭

点点滴滴

感悟人生,享受精彩。

 
 
 

日志

 
 

第5章 优化问题  

2009-03-30 15:47:33|  分类: 默认分类 |  标签: |举报 |字号 订阅

  下载LOFTER 我的照片书  |
第5章  优化问题
5.1  线性规划问题
线性规划问题是目标函数和约束条件均为线性函数的问题,MATLAB6.0解决的线性规划问题的标准形式为:
min  
sub.to:
       
       
其中f、x、b、beq、lb、ub为向量,A、Aeq为矩阵。
其它形式的线性规划问题都可经过适当变换化为此标准形式。
在MATLAB6.0版中,线性规划问题(Linear Programming)已用函数linprog取代了MATLAB5.x版中的lp函数。当然,由于版本的向下兼容性,一般说来,低版本中的函数在6.0版中仍可使用。
函数  linprog
格式  x = linprog(f,A,b)   %求min  f ' *x   sub.to   线性规划的最优解。
x = linprog(f,A,b,Aeq,beq)   %等式约束 ,若没有不等式约束 ,则A=[ ],b=[ ]。
x = linprog(f,A,b,Aeq,beq,lb,ub)   %指定x的范围 ,若没有等式约束  ,则Aeq=[ ],beq=[ ]
x = linprog(f,A,b,Aeq,beq,lb,ub,x0)   %设置初值x0
x = linprog(f,A,b,Aeq,beq,lb,ub,x0,options)    % options为指定的优化参数
[x,fval] = linprog(…)   % 返回目标函数最优值,即fval= f ' *x。
[x,lambda,exitflag] = linprog(…)   % lambda为解x的Lagrange乘子。
[x, lambda,fval,exitflag] = linprog(…)   % exitflag为终止迭代的错误条件。
[x,fval, lambda,exitflag,output] = linprog(…)   % output为关于优化的一些信息
说明  若exitflag>0表示函数收敛于解x,exitflag=0表示超过函数估值或迭代的最大数字,exitflag<0表示函数不收敛于解x;若lambda=lower 表示下界lb,lambda=upper表示上界ub,lambda=ineqlin表示不等式约束,lambda=eqlin表示等式约束,lambda中的非0元素表示对应的约束是有效约束;output=iterations表示迭代次数,output=algorithm表示使用的运算规则,output=cgiterations表示PCG迭代次数。
例5-1  求下面的优化问题
min  
sub.to      
        
         
 
解:
>>f = [-5; -4; -6];
>>A =  [1 -1  1;3  2  4;3  2  0];
>>b = [20; 42; 30];
>>lb = zeros(3,1);
>>[x,fval,exitflag,output,lambda] = linprog(f,A,b,[],[],lb)
结果为:
x =      %最优解
    0.0000
   15.0000
    3.0000
fval =     %最优值
  -78.0000
exitflag =     %收敛
     1
output =
      iterations: 6   %迭代次数
    cgiterations: 0
       algorithm: 'lipsol'   %所使用规则
lambda =
    ineqlin: [3x1 double]
      eqlin: [0x1 double]
      upper: [3x1 double]
      lower: [3x1 double]
>> lambda.ineqlin
ans =
    0.0000
    1.5000
    0.5000
>> lambda.lower
ans =
    1.0000
    0.0000
    0.0000
表明:不等约束条件2和3以及第1个下界是有效的
5.2  foptions函数
对于优化控制,MATLAB提供了18个参数,这些参数的具体意义为:
options(1)-参数显示控制(默认值为0)。等于1时显示一些结果。
    options(2)-优化点x的精度控制(默认值为1e-4)。
    options(3)-优化函数F的精度控制(默认值为1e-4)。
    options(4)-违反约束的结束标准(默认值为1e-6)。
    options(5)-算法选择,不常用。
    options(6)-优化程序方法选择,为0则为BFCG算法,为1则采用DFP算法。
    options(7)-线性插值算法选择,为0则为混合插值算法,为1则采用立方插算法。
    options(8)-函数值显示 (目标—达到问题中的Lambda )
    options(9)-若需要检测用户提供的梯度,则设为1。
    options(10)-函数和约束估值的数目。
    options(11)-函数梯度估值的个数。
    options(12)-约束估值的数目。
    options(13)-等约束条件的个数。
    options(14)-函数估值的最大次数(默认值是100×变量个数)
    options(15)-用于目标 — 达到问题中的特殊目标。
    options(16)-优化过程中变量的最小有限差分梯度值。
    options(17)- 优化过程中变量的最大有限差分梯度值。
    options(18)-步长设置 (默认为1或更小)。
Foptions已经被optimset和optimget代替,详情请查函数optimset和optimget。
5.3  非线性规划问题
5.3.1  有约束的一元函数的最小值
单变量函数求最小值的标准形式为      sub.to  
在MATLAB5.x中使用fmin函数求其最小值。
函数  fminbnd
格式  x = fminbnd(fun,x1,x2)   %返回自变量x在区间 上函数fun取最小值时x值,fun为目标函数的表达式字符串或MATLAB自定义函数的函数柄。
x = fminbnd(fun,x1,x2,options)   % options为指定优化参数选项
[x,fval] = fminbnd(…)   % fval为目标函数的最小值
[x,fval,exitflag] = fminbnd(…)   %xitflag为终止迭代的条件
[x,fval,exitflag,output] = fminbnd(…)   % output为优化信息
说明  若参数exitflag>0,表示函数收敛于x,若exitflag=0,表示超过函数估计值或迭代的最大数字,exitflag<0表示函数不收敛于x;若参数output=iterations表示迭代次数,output=funccount表示函数赋值次数,output=algorithm表示所使用的算法。
例5-2  计算下面函数在区间(0,1)内的最小值。
 
解:>> [x,fval,exitflag,output]=fminbnd('(x^3+cos(x)+x*log(x))/exp(x)',0,1)
x =
    0.5223
fval =
    0.3974
exitflag =
     1
output =
    iterations: 9
     funcCount: 9
     algorithm: 'golden section search, parabolic interpolation'
例5-3  在[0,5]上求下面函数的最小值
解:先自定义函数:在MATLAB编辑器中建立M文件为:
function f = myfun(x)
f = (x-3).^2 - 1;
保存为myfun.m,然后在命令窗口键入命令:
>> x=fminbnd(@myfun,0,5)
则结果显示为:
x =
     3
5.3.2  无约束多元函数最小值
多元函数最小值的标准形式为
其中:x为向量,如
在MATLAB5.x中使用fmins求其最小值。
命令  利用函数fminsearch求无约束多元函数最小值
函数  fminsearch
格式  x = fminsearch(fun,x0)   %x0为初始点,fun为目标函数的表达式字符串或MATLAB自定义函数的函数柄。
x = fminsearch(fun,x0,options)   % options查optimset
[x,fval] = fminsearch(…)   %最优点的函数值
[x,fval,exitflag] = fminsearch(…)   % exitflag与单变量情形一致
[x,fval,exitflag,output] = fminsearch(…)  %output与单变量情形一致
注意:fminsearch采用了Nelder-Mead型简单搜寻法。
例5-4  求 的最小值点
解:>>X=fminsearch('2*x(1)^3+4*x(1)*x(2)^3-10*x(1)*x(2)+x(2)^2',  [0,0])
结果为
X =
    1.0016    0.8335
或在MATLAB编辑器中建立函数文件
function  f=myfun(x)
f=2*x(1)^3+4*x(1)*x(2)^3-10*x(1)*x(2)+x(2)^2;
保存为myfun.m,在命令窗口键入
>> X=fminsearch ('myfun',  [0,0]) 或 >> X=fminsearch(@myfun,  [0,0])
结果为:
X =
    1.0016    0.8335
命令  利用函数fminunc求多变量无约束函数最小值
函数  fminunc
格式  x = fminunc(fun,x0)   %返回给定初始点x0的最小函数值点
x = fminunc(fun,x0,options)   % options为指定优化参数
[x,fval] = fminunc(…)   %fval最优点x处的函数值
[x,fval,exitflag] = fminunc(…)   % exitflag为终止迭代的条件,与上同。
[x,fval,exitflag,output] = fminunc(…)   %output为输出优化信息
[x,fval,exitflag,output,grad] = fminunc(…)   % grad为函数在解x处的梯度值
[x,fval,exitflag,output,grad,hessian] = fminunc(…)   %目标函数在解x处的海赛(Hessian)值
注意:当函数的阶数大于2时,使用fminunc比fminsearch更有效,但当所选函数高度不连续时,使用fminsearch效果较好。
例5-5  求 的最小值。
>> fun='3*x(1)^2+2*x(1)*x(2)+x(2)^2';
>> x0=[1 1];
>> [x,fval,exitflag,output,grad,hessian]=fminunc(fun,x0)
结果为:
x =
  1.0e-008 *
   -0.7591    0.2665
fval =
  1.3953e-016
exitflag =
     1
output =
       iterations: 3
        funcCount: 16
         stepsize: 1.2353
    firstorderopt: 1.6772e-007
        algorithm: 'medium-scale: Quasi-Newton line search'
grad =
  1.0e-006 *
   -0.1677
    0.0114
hessian =
    6.0000    2.0000
    2.0000    2.0000
或用下面方法:
>> fun=inline('3*x(1)^2+2*x(1)*x(2)+x(2)^2')
fun =
     Inline function:
     fun(x) = 3*x(1)^2+2*x(1)*x(2)+x(2)^2
>> x0=[1 1];
>> x=fminunc(fun,x0)
x =
  1.0e-008 *
   -0.7591    0.2665
5.3.3  有约束的多元函数最小值
非线性有约束的多元函数的标准形式为:
 
sub.to    
 
 
 
 
其中:x、b、beq、lb、ub是向量,A、Aeq为矩阵,C(x)、Ceq(x)是返回向量的函数,f(x)为目标函数,f(x)、C(x)、Ceq(x)可以是非线性函数。
在MATLAB5.x中,它的求解由函数constr实现。
函数  fmincon
格式  x = fmincon(fun,x0,A,b)
x = fmincon(fun,x0,A,b,Aeq,beq)
x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub)
x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon)
x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options)
[x,fval] = fmincon(…)
[x,fval,exitflag] = fmincon(…)
[x,fval,exitflag,output] = fmincon(…)
[x,fval,exitflag,output,lambda] = fmincon(…)
[x,fval,exitflag,output,lambda,grad] = fmincon(…)
[x,fval,exitflag,output,lambda,grad,hessian] = fmincon(…)
参数说明:fun为目标函数,它可用前面的方法定义;
x0为初始值;
A、b满足线性不等式约束 ,若没有不等式约束,则取A=[ ],b=[ ];
Aeq、beq满足等式约束 ,若没有,则取Aeq=[ ],beq=[ ];
lb、ub满足 ,若没有界,可设lb=[ ],ub=[ ];
nonlcon的作用是通过接受的向量x来计算非线性不等约束 和等式约束 分别在x处的估计C和Ceq,通过指定函数柄来使用,如:>>x = fmincon(@myfun,x0,A,b,Aeq,beq,lb,ub,@mycon),先建立非线性约束函数,并保存为mycon.m:function [C,Ceq] = mycon(x)
C = …    % 计算x处的非线性不等约束 的函数值。
Ceq = …   % 计算x处的非线性等式约束 的函数值。
lambda是Lagrange乘子,它体现哪一个约束有效。
output输出优化信息;
grad表示目标函数在x处的梯度;
hessian表示目标函数在x处的Hessiab值。
例5-6  求下面问题在初始点(0,1)处的最优解
min  
sub.to    
 
解:约束条件的标准形式为
sub.to  
 
先在MATLAB编辑器中建立非线性约束函数文件:
function  [c, ceq]=mycon (x)
c=(x(1)-1)^2-x(2);
ceq=[ ];      %无等式约束
然后,在命令窗口键入如下命令或建立M文件:
>>fun='x(1)^2+x(2)^2-x(1)*x(2)-2*x(1)-5*x(2)';     %目标函数
>>x0=[0 1];
>>A=[-2 3];   %线性不等式约束
>>b=6;
>>Aeq=[ ];    %无线性等式约束
>>beq=[ ];
>>lb=[ ];       %x没有下、上界
>>ub=[ ];
>>[x,fval,exitflag,output,lambda,grad,hessian]
=fmincon(fun,x0,A,b,Aeq,beq,lb,ub,@mycon)
则结果为
x =
     3     4
fval =
   -13
exitflag =      %解收敛
     1
output =
       iterations: 2
        funcCount: 9
         stepsize: 1
        algorithm: 'medium-scale: SQP, Quasi-Newton, line-search'
    firstorderopt: [ ]
     cgiterations: [ ]
lambda =
         lower: [2x1 double]   %x下界有效情况,通过lambda.lower可查看。
         upper: [2x1 double]   %x上界有效情况,为0表示约束无效。
         eqlin: [0x1 double]    %线性等式约束有效情况,不为0表示约束有效。
      eqnonlin: [0x1 double]    %非线性等式约束有效情况。
       ineqlin: 2.5081e-008    %线性不等式约束有效情况。
    ineqnonlin: 6.1938e-008    %非线性不等式约束有效情况。
grad =      %目标函数在最小值点的梯度
  1.0e-006 *
   -0.1776
         0
hessian =    %目标函数在最小值点的Hessian值
    1.0000   -0.0000
   -0.0000    1.0000
例5-7  求下面问题在初始点x=(10, 10, 10)处的最优解。
Min  
Sub.to  
解:约束条件的标准形式为
sub.to     
 
>> fun= '-x(1)*x(2)*x(3)';
>> x0=[10,10,10];
>> A=[-1 -2 -2;1 2 2];
>> b=[0;72];
>> [x,fval]=fmincon(fun,x0,A,b)
结果为:
x =
   24.0000   12.0000   12.0000
fval =
       -3456
5.3.4  二次规划问题
二次规划问题(quadratic programming)的标准形式为:
 
sub.to  
    
    
其中,H、A、Aeq为矩阵,f、b、beq、lb、ub、x为向量
其它形式的二次规划问题都可转化为标准形式。
MATLAB5.x版中的qp函数已被6.0版中的函数quadprog取代。
函数  quadprog
格式  x = quadprog(H,f,A,b)   %其中H,f,A,b为标准形中的参数,x为目标函数的最小值。
x = quadprog(H,f,A,b,Aeq,beq)   %Aeq,beq满足等约束条件 。
x = quadprog(H,f,A,b,Aeq,beq,lb,ub)   % lb,ub分别为解x的下界与上界。
x = quadprog(H,f,A,b,Aeq,beq,lb,ub,x0)   %x0为设置的初值
x = quadprog(H,f,A,b,Aeq,beq,lb,ub,x0,options)   % options为指定的优化参数
[x,fval] = quadprog(…)   %fval为目标函数最优值
[x,fval,exitflag] = quadprog(…)   % exitflag与线性规划中参数意义相同
[x,fval,exitflag,output] = quadprog(…)   % output与线性规划中参数意义相同
[x,fval,exitflag,output,lambda] = quadprog(…)   % lambda与线性规划中参数意义相同
例5-8  求解下面二次规划问题
 
sub.to
 
 
 
 
解:
则 , ,
在MATLAB中实现如下:
>>H = [1 -1; -1 2] ;
>>f = [-2; -6];
>>A = [1 1; -1 2; 2 1];
>>b = [2; 2; 3];
>>lb = zeros(2,1);
>>[x,fval,exitflag,output,lambda] = quadprog(H,f,A,b,[ ],[ ],lb)
结果为:
x =     %最优解
    0.6667
    1.3333
fval =     %最优值
   -8.2222
exitflag =    %收敛
     1
output =
       iterations: 3
        algorithm: 'medium-scale: active-set'
    firstorderopt: [ ]
     cgiterations: [ ]
lambda =
      lower: [2x1 double]
      upper: [2x1 double]
      eqlin: [0x1 double]
    ineqlin: [3x1 double]
>> lambda.ineqlin
ans =
    3.1111
    0.4444
         0
>> lambda.lower
ans =
     0
     0
说明  第1、2个约束条件有效,其余无效。
例5-9  求二次规划的最优解
max   f (x1, x2)=x1x2+3
sub.to    x1+x2-2=0
解:化成标准形式:
 
sub.to    x1+x2=2
在Matlab中实现如下:
>>H=[0,-1;-1,0];
>>f=[0;0];
>>Aeq=[1 1];
>>b=2;
>>[x,fval,exitflag,output,lambda] = quadprog(H,f,[ ],[ ],Aeq,b)
结果为:
x =
    1.0000
    1.0000
fval =
   -1.0000
exitflag =
     1
output =
    firstorderopt: 0
       iterations: 1
     cgiterations: 1
        algorithm: [1x58 char]
lambda =
      eqlin: 1.0000
    ineqlin: [ ]
      lower: [ ]
      upper: [ ]
5.4  “半无限”有约束的多元函数最优解
“半无限”有约束多元函数最优解问题的标准形式为
 
sub.to      
 
 
 
 
 

 
其中:x、b、beq、lb、ub都是向量;A、Aeq是矩阵;C(x)、Ceq(x)、 是返回向量的函数,f(x)为目标函数;f(x)、C(x)、Ceq(x)是非线性函数; 为半无限约束, 通常是长度为2的向量。
在MTALAB5.x中,使用函数seminf解决这类问题。
函数  fseminf
格式  x = fseminf(fun,x0,ntheta,seminfcon)
x = fseminf(fun,x0,ntheta,seminfcon,A,b)
x = fseminf(fun,x0,ntheta,seminfcon,A,b,Aeq,beq)
x = fseminf(fun,x0,ntheta,seminfcon,A,b,Aeq,beq,lb,ub)
x = fseminf(fun,x0,ntheta,seminfcon,A,b,Aeq,beq,lb,ub,options)
[x,fval] = fseminf(…)
[x,fval,exitflag] = fseminf(…)
[x,fval,exitflag,output] = fseminf(…)
[x,fval,exitflag,output,lambda] = fseminf(…)
参数说明:x0为初始估计值;
fun为目标函数,其定义方式与前面相同;
A、b由线性不等式约束 确定,没有,则A=[ ],b=[ ];
Aeq、beq由线性等式约束 确定,没有,则Aeq=[ ],beq=[ ];
Lb、ub由变量x的范围 确定;
options为优化参数;
ntheta为半无限约束的个数;
seminfcon用来确定非线性约束向量C和Ceq以及半无限约束的向量K1,K2,…,Kn,通过指定函数柄来使用,如:
x = fseminf(@myfun,x0,ntheta,@myinfcon)
先建立非线性约束和半无限约束函数文件,并保存为myinfcon.m:
function [C,Ceq,K1,K2,…,Kntheta,S] = myinfcon(x,S)
%S为向量w的采样值
% 初始化样本间距
if  isnan(S(1,1)),
S = …    % S 有ntheta行2列
end
w1 = …     %计算样本集
w2 = …     %计算样本集

wntheta = …   % 计算样本集
K1 = …     % 在x和w处的第1个半无限约束值
K2 = …     %在x和w处的第2个半无限约束值

Kntheta = …  %在x和w处的第ntheta个半无限约束值
C = …        % 在x处计算非线性不等式约束值
Ceq = …      % 在x处计算非线性等式约束值
如果没有约束,则相应的值取为“[ ]”,如Ceq=[]
fval为在x处的目标函数最小值;
exitflag为终止迭代的条件;
output为输出的优化信息;
lambda为解x的Lagrange乘子。
例5-10  求下面一维情形的最优化问题
 
sub.to
 
 
 
 
解:将约束方程化为标准形式:
 
 
先建立非线性约束和半无限约束函数文件,并保存为mycon.m:
function [C,Ceq,K1,K2,S] = mycon(X,S)
% 初始化样本间距:
if  isnan(S(1,1)),
   S = [0.2  0; 0.2  0];
end
% 产生样本集:
w1 = 1:S(1,1):100;
w2 = 1:S(2,1):100;
% 计算半无限约束:
K1 = sin(w1*X(1)).*cos(w1*X(2)) - 1/1000*(w1-50).^2 -sin(w1*X(3))-X(3)-1;
K2 = sin(w2*X(2)).*cos(w2*X(1)) - 1/1000*(w2-50).^2 -sin(w2*X(3))-X(3)-1;
% 无非线性约束:
C = [ ]; Ceq=[ ];
% 绘制半无限约束图形
plot(w1,K1,'-',w2,K2,':'),title('Semi-infinite constraints')
然后在MATLAB命令窗口或编辑器中建立M文件:
fun = 'sum((x-0.5).^2)';
x0 = [0.5; 0.2; 0.3];      % Starting guess
[x,fval] = fseminf(fun,x0,2,@mycon)
结果为:
x =
    0.6673
    0.3013
    0.4023
fval =
    0.0770
>>[C,Ceq,K1,K2] = mycon (x,NaN);    % 利用初始样本间距
>>max(K1)
ans =
      -0.0017
>>max(K2)
ans =
      -0.0845
图5-1
例5-11  求下面二维情形的最优化问题
 
sub.to
 
 
 
 
初始点为x0=[0.25, 0.25, 0.25]。
解:先建立非线性和半无限约束函数文件,并保存为mycon.m:
function [C,Ceq,K1,S] = mycon(X,S)
% 初始化样本间距:
if  isnan(s(1,1)),
   s = [2 2];
end
% 设置样本集
w1x = 1:s(1,1):100;
w1y = 1:s(1,2):100;
[wx, wy] = meshgrid(w1x,w1y);
% 计算半无限约束函数值
K1 = sin(wx*X(1)).*cos(wx*X(2))-1/1000*(wx-50).^2 -sin(wx*X(3))-X(3)+…
sin(wy*X(2)).*cos(wx*X(1))-1/1000*(wy-50).^2-sin(wy*X(3))-X(3)-1.5;
% 无非线性约束
C = [ ]; Ceq=[ ];
%作约束曲面图形
m = surf(wx,wy,K1,'edgecolor','none','facecolor','interp');
camlight headlight
title('Semi-infinite constraint')
drawnow
然后在MATLAB命令窗口下键入命令:
>>fun = 'sum((x-0.2).^2)';
>>x0 = [0.25, 0.25, 0.25];
>>[x,fval] = fseminf(fun,x0,1,@mycon)
结果为(如图)
x =
    0.2926    0.1874    0.2202
fval =
    0.0091
>>[c,ceq,K1] = mycon(x,[0.5,0.5]);  % 样本间距为0.5
>>max(max(K1))
ans =
    -0.0027
5.5  极小化极大(Minmax)问题
极小化极大问题的标准形式为
 
sub.to    
 
 
 
 
其中:x、b、beq、lb、ub是向量,A、Aeq为矩阵,C(x)、Ceq(x)和F(x)是返回向量的函数,F(x)、C(x)、Ceq(x)可以是非线性函数。
在MATLAB5.x中,它的求解由函数minmax实现。
函数  fminimax
格式  x = fminimax(fun,x0)
x = fminimax(fun,x0,A,b)
x = fminimax(fun,x0,A,b,Aeq,beq)
x = fminimax(fun,x0,A,b,Aeq,beq,lb,ub)
x = fminimax(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon)
x = fminimax(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options)
[x,fval,maxfval] = fminimax(…)
[x,fval,maxfval,exitflag] = fminimax(…)
[x,fval,maxfval,exitflag,output] = fminimax(…)
[x,fval,maxfval,exitflag,output,lambda] = fminimax(…)
参数说明:fun为目标函数;
x0为初始值;
A、b满足线性不等约束 ,若没有不等约束,则取A=[ ],b=[ ];
Aeq、beq满足等式约束 ,若没有,则取Aeq=[ ],beq=[ ];
lb、ub满足 ,若没有界,可设lb=[ ],ub=[ ];
nonlcon的作用是通过接受的向量x来计算非线性不等约束 和等式约束 分别在x处的值C和Ceq,通过指定函数柄来使用,如:>>x = fminimax(@myfun,x0,A,b,Aeq,beq,lb,ub,@mycon),先建立非线性约束函数,并保存为mycon.m:function [C,Ceq] = mycon(x)
C = …    % 计算x处的非线性不等约束 的函数值。
Ceq = …   % 计算x处的非线性等式约束 的函数值。
options为指定的优化参数;
fval为最优点处的目标函数值;
maxfval为目标函数在x处的最大值;
exitflag为终止迭代的条件;
lambda是Lagrange乘子,它体现哪一个约束有效。
output输出优化信息。
例5-12  求下列函数最大值的最小化问题
 
其中:
 
 
 
 
解:先建立目标函数文件,并保存为myfun.m:function f = myfun(x)
f(1)= 2*x(1)^2+x(2)^2-48*x(1)-40*x(2)+304;
f(2)= -x(1)^2 - 3*x(2)^2;
f(3)= x(1) + 3*x(2) -18;
f(4)= -x(1)- x(2);
f(5)= x(1) + x(2) - 8;
然后,在命令窗口键入命令:
x0 = [0.1; 0.1];       % 初始值
[x,fval] = fminimax(@myfun,x0)
结果为:
x =
    4.0000
    4.0000
fval =
    0.0000  -64.0000   -2.0000   -8.0000   -0.0000
例5-13  求上述问题的绝对值的最大值最小化问题。
目标函数为:
解:先建立目标函数文件(与上例相同)
然后,在命令窗口或编辑器中建立M文件:
>>x0 = [0.1; 0.1];        % 初始点
>>options = optimset('MinAbsMax',5);   % 指定绝对值的最小化
>>[x,fval] = fminimax(@myfun,x0,[ ],[ ],[ ],[ ],[ ],[ ],[ ],options)
则结果为
x =
    4.9256
    2.0796
fval =
   37.2356  -37.2356   -6.8357   -7.0052   -0.9948
5.6  多目标规划问题
多目标规划是指在一组约束下,对多个不同目标函数进行优化。它的一般形式为
 
sub.to    
其中: 。
在同一约束下,当目标函数处于冲突状态时,不存在最优解x使所有目标函数同时达到最优。此时,我们使用有效解,即如果不存在 ,使得 ,i=1,2,…m, 则称x*为有效解。
在MATLAB中,多目标问题的标准形式为
 
sub.to     
 
 
 
 
 
其中:x、b、beq、lb、ub是向量;A、Aeq为矩阵;C(x)、Ceq(x)和F(x)是返回向量的函数;F(x)、C(x)、Ceq(x)可以是非线性函数;weight为权值系数向量,用于控制对应的目标函数与用户定义的目标函数值的接近程度;goal为用户设计的与目标函数相应的目标函数值向量; 为一个松弛因子标量;F(x)为多目标规划中的目标函数向量。
在MATLAB5.x中,它的最优解由attgoal函数实现。
函数  fgoalattain
格式  x = fgoalattain(fun,x0,goal,weight)
x = fgoalattain(fun,x0,goal,weight,A,b)
x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq)
x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,lb,ub)
x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,lb,ub,nonlcon)
x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,lb,ub,nonlcon,options)
[x,fval] = fgoalattain(…)
[x,fval,attainfactor] = fgoalattain(…)
[x,fval,attainfactor,exitflag] = fgoalattain(…)
[x,fval,attainfactor,exitflag,output] = fgoalattain(…)
[x,fval,attainfactor,exitflag,output,lambda] = fgoalattain(…)
参数说明:
x0为初始解向量;
fun为多目标函数的文件名字符串,其定义方式与前面fun的定义方式相同;
goal为用户设计的目标函数值向量;
weight为权值系数向量,用于控制目标函数与用户自定义目标值的接近程度;
A、b满足线性不等式约束 ,没有时取A=[ ],b=[ ];
Aeq、beq满足线性等式约束 ,没有时取Aeq=[ ],beq=[ ];
lb、ub为变量的下界和上界: ;
nonlcon的作用是通过接受的向量x来计算非线性不等约束 和等式约束 分别在x处的值C和Ceq,通过指定函数柄来使用。
如:>>x = fgoalattain(@myfun,x0,goal,weight,A,b,Aeq,beq,lb,ub,@mycon),先建立非线性约束函数,并保存为mycon.m:function  [C,Ceq] = mycon(x)
C = …    % 计算x处的非线性不等式约束 的函数值。
Ceq = …   % 计算x处的非线性等式约束 的函数值。
options为指定的优化参数;
fval为多目标函数在x处的值;
attainfactor为解x处的目标规划因子;
exitflag为终止迭代的条件;
output为输出的优化信息;
lambda为解x处的Lagrange乘子
例5-14    控制系统输出反馈器设计。
设如下线性系统
 
 
其中:            
要求设计输出反馈控制器K,使闭环系统
 
 
在复平面实轴上点[-5,-3,-1]的左侧有极点,并要求 
解:上述问题就是要求解矩阵K,使矩阵(A+BKC)的极点为[-5,-3,-1],这是一个多目标规划问题。
先建立目标函数文件,保存为eigfun.m:
function F = eigfun(K,A,B,C)
F = sort(eig(A+B*K*C));    % 估计目标函数值
然后,输入参数并调用优化程序:
A = [-0.5 0 0; 0 -2 10; 0 1 -2];
B = [1 0; -2 2; 0 1];
C = [1 0 0; 0 0 1];
K0 = [-1 -1; -1 -1];       % 初始化控制器矩阵
goal = [-5 -3 -1];         % 为闭合环路的特征值(极点)设置目标值向量
weight = abs(goal)        % 设置权值向量
lb = -4*ones(size(K0));    % 设置控制器的下界
ub = 4*ones(size(K0));     % 设置控制器的上界
options = optimset('Display','iter');   % 设置显示参数:显示每次迭代的输出
[K,fval,attainfactor] = fgoalattain(@eigfun,K0,goal,weight,[],[],[],[],lb,ub,[],options,A,B,C)
结果为:
weight =
     5     3     1
                    Attainment                  Directional
 Iter   F-count       factor         Step-size     derivative      Procedure
    1      6          1.885            1            1.03    
    2     13          1.061            1          -0.679    
    3     20         0.4211            1          -0.523    Hessian modified 
    4     27       -0.06352            1          -0.053    Hessian modified twice 
    5     34        -0.1571            1          -0.133    
    6     41        -0.3489            1        -0.00768    Hessian modified 
    7     48        -0.3643            1      -4.25e-005    Hessian modified 
    8     55        -0.3645            1        -0.00303    Hessian modified twice 
    9     62        -0.3674            1         -0.0213    Hessian modified 
   10     69        -0.3806            1        0.00266    
   11     76        -0.3862            1      -2.73e-005    Hessian modified twice 
   12     83        -0.3863            1      -1.22e-013    Hessian modified twice 
Optimization terminated successfully:
    Search direction less than 2*options. TolX and maximum constraint violation is less than options.TolCon
Active Constraints:
     1
     2
     4
     9
    10
K =
   -4.0000   -0.2564
   -4.0000   -4.0000
fval =
   -6.9313
   -4.1588
   -1.4099
attainfactor =
   -0.3863
5.7  最小二乘最优问题
5.7.1  约束线性最小二乘
有约束线性最小二乘的标准形式为
 
sub.to   
 
 
其中:C、A、Aeq为矩阵;d、b、beq、lb、ub、x是向量。
在MATLAB5.x中,约束线性最小二乘用函数conls求解。
函数  lsqlin  
格式  x = lsqlin(C,d,A,b)   %求在约束条件 下,方程Cx = d的最小二乘解x。
x = lsqlin(C,d,A,b,Aeq,beq)   %Aeq、beq满足等式约束 ,若没有不等式约束,则设A=[ ],b=[ ]。
x = lsqlin(C,d,A,b,Aeq,beq,lb,ub)   %lb、ub满足 ,若没有等式约束,则Aeq=[ ],beq=[ ]。
x = lsqlin(C,d,A,b,Aeq,beq,lb,ub,x0)   % x0为初始解向量,若x没有界,则lb=[ ],ub=[ ]。
x = lsqlin(C,d,A,b,Aeq,beq,lb,ub,x0,options)   % options为指定优化参数
[x,resnorm] = lsqlin(…)   % resnorm=norm(C*x-d)^2,即2-范数。
[x,resnorm,residual] = lsqlin(…)   %residual=C*x-d,即残差。
[x,resnorm,residual,exitflag] = lsqlin(…)   %exitflag为终止迭代的条件
[x,resnorm,residual,exitflag,output] = lsqlin(…)   % output表示输出优化信息
[x,resnorm,residual,exitflag,output,lambda] = lsqlin(…)   % lambda为解x的Lagrange乘子
例5-15  求解下面系统的最小二乘解
系统:Cx=d
约束: ;
先输入系统系数和x的上下界:
C = [0.9501    0.7620    0.6153    0.4057;…
    0.2311    0.4564    0.7919    0.9354;…
    0.6068    0.0185    0.9218    0.9169;…
    0.4859    0.8214    0.7382    0.4102;…
    0.8912    0.4447    0.1762    0.8936];
d = [ 0.0578; 0.3528; 0.8131; 0.0098; 0.1388];
A =[ 0.2027    0.2721    0.7467    0.4659;…
    0.1987    0.1988    0.4450    0.4186;…
    0.6037    0.0152    0.9318    0.8462];
b =[ 0.5251; 0.2026; 0.6721];
lb = -0.1*ones(4,1);
ub = 2*ones(4,1);
然后调用最小二乘命令:
[x,resnorm,residual,exitflag,output,lambda] = lsqlin(C,d,A,b,[ ],[ ],lb,ub);
结果为:
x =
   -0.1000
   -0.1000
    0.2152
    0.3502
resnorm =
    0.1672
residual =
    0.0455
    0.0764
   -0.3562
    0.1620
    0.0784
exitflag =
     1      %说明解x是收敛的
output =
       iterations: 4
        algorithm: 'medium-scale: active-set'
    firstorderopt: []
     cgiterations: []
lambda =
      lower: [4x1 double]
      upper: [4x1 double]
      eqlin: [0x1 double]
    ineqlin: [3x1 double]
通过lambda.ineqlin可查看非线性不等式约束是否有效。
5.7.2  非线性数据(曲线)拟合
非线性曲线拟合是已知输入向量xdata和输出向量ydata,并且知道输入与输出的函数关系为ydata=F(x, xdata),但不知道系数向量x。今进行曲线拟合,求x使得下式成立:
 
在MATLAB5.x中,使用函数curvefit解决这类问题。
函数  lsqcurvefit
格式  x = lsqcurvefit(fun,x0,xdata,ydata)
x = lsqcurvefit(fun,x0,xdata,ydata,lb,ub)
x = lsqcurvefit(fun,x0,xdata,ydata,lb,ub,options)
[x,resnorm] = lsqcurvefit(…)
[x,resnorm,residual] = lsqcurvefit(…)
[x,resnorm,residual,exitflag] = lsqcurvefit(…)
[x,resnorm,residual,exitflag,output] = lsqcurvefit(…)
[x,resnorm,residual,exitflag,output,lambda] = lsqcurvefit(…)
[x,resnorm,residual,exitflag,output,lambda,jacobian] =lsqcurvefit(…)
参数说明:
x0为初始解向量;xdata,ydata为满足关系ydata=F(x, xdata)的数据;
lb、ub为解向量的下界和上界 ,若没有指定界,则lb=[ ],ub=[ ];
options为指定的优化参数;
fun为拟合函数,其定义方式为:x = lsqcurvefit(@myfun,x0,xdata,ydata),
其中myfun已定义为     function F = myfun(x,xdata)
F = …      % 计算x处拟合函数值fun的用法与前面相同;
resnorm=sum ((fun(x,xdata)-ydata).^2),即在x处残差的平方和;
residual=fun(x,xdata)-ydata,即在x处的残差;
exitflag为终止迭代的条件;
output为输出的优化信息;
lambda为解x处的Lagrange乘子;
jacobian为解x处拟合函数fun的jacobian矩阵。
例5-16  求解如下最小二乘非线性拟合问题
已知输入向量xdata和输出向量ydata,且长度都是n,拟合函数为
 
即目标函数为
其中:
初始解向量为x0=[0.3, 0.4, 0.1]。
解:先建立拟合函数文件,并保存为myfun.m
function F = myfun(x,xdata)
F = x(1)*xdata.^2 + x(2)*sin(xdata) + x(3)*xdata.^3;
然后给出数据xdata和ydata
>>xdata = [3.6 7.7 9.3 4.1 8.6 2.8 1.3 7.9 10.0 5.4];
>>ydata = [16.5 150.6 263.1 24.7 208.5 9.9 2.7 163.9 325.0 54.3];
>>x0 = [10, 10, 10];    %初始估计值
>>[x,resnorm] = lsqcurvefit(@myfun,x0,xdata,ydata)
结果为:
Optimization terminated successfully:
Relative function value changing by less than OPTIONS.TolFun
x =
0.2269    0.3385    0.3021
resnorm =
     6.2950
5.7.3  非线性最小二乘
非线性最小二乘(非线性数据拟合)的标准形式为
 
其中:L为常数
在MATLAB5.x中,用函数leastsq解决这类问题,在6.0版中使用函数lsqnonlin。

则目标函数可表达为
其中:x为向量,F(x)为函数向量。
函数  lsqnonlin
格式  x = lsqnonlin(fun,x0)   %x0为初始解向量;fun为 ,i=1,2,…,m,fun返回向量值F,而不是平方和值,平方和隐含在算法中,fun的定义与前面相同。
x = lsqnonlin(fun,x0,lb,ub)    %lb、ub定义x的下界和上界: 。
x = lsqnonlin(fun,x0,lb,ub,options)   %options为指定优化参数,若x没有界,则lb=[ ],ub=[ ]。
[x,resnorm] = lsqnonlin(…)    % resnorm=sum(fun(x).^2),即解x处目标函数值。
[x,resnorm,residual] = lsqnonlin(…)   % residual=fun(x),即解x处fun的值。
[x,resnorm,residual,exitflag] = lsqnonlin(…)    %exitflag为终止迭代条件。
[x,resnorm,residual,exitflag,output] = lsqnonlin(…)   %output输出优化信息。
[x,resnorm,residual,exitflag,output,lambda] = lsqnonlin(…)   %lambda为Lagrage乘子。
[x,resnorm,residual,exitflag,output,lambda,jacobian] =lsqnonlin(…)   %fun在解x处的Jacobian矩阵。
例5-17  求下面非线性最小二乘问题 初始解向量为x0=[0.3, 0.4]。
解:先建立函数文件,并保存为myfun.m,由于lsqnonlin中的fun为向量形式而不是平方和形式,因此,myfun函数应由 建立:
     k=1,2,…,10
function  F = myfun(x)
k = 1:10;
F = 2 + 2*k-exp(k*x(1))-exp(k*x(2));
然后调用优化程序:
x0 = [0.3 0.4];
[x,resnorm] = lsqnonlin(@myfun,x0)
结果为:
Optimization terminated successfully:
Norm of the current step is less than OPTIONS.TolX
x =
    0.2578    0.2578
resnorm =    %求目标函数值
  124.3622
5.7.4  非负线性最小二乘
非负线性最小二乘的标准形式为:
 
sub.to    
其中:矩阵C和向量d为目标函数的系数,向量x为非负独立变量。
在MATLAB5.x中,用函数nnls求解这类问题,在6.0版中则用函数lsqnonneg。
函数  lsqnonneg
格式  x = lsqnonneg(C,d)   %C为实矩阵,d为实向量
x = lsqnonneg(C,d,x0)   % x0为初始值且大于0
x = lsqnonneg(C,d,x0,options)   % options为指定优化参数
[x,resnorm] = lsqnonneg(…)   % resnorm=norm (C*x-d)^2
[x,resnorm,residual] = lsqnonneg(…)   %residual=C*x-d
[x,resnorm,residual,exitflag] = lsqnonneg(…)
[x,resnorm,residual,exitflag,output] = lsqnonneg(…)
[x,resnorm,residual,exitflag,output,lambda] = lsqnonneg(…)
例5-18  一个最小二乘问题的无约束与非负约束解法的比较。
先输入数据:
>>C = [ 0.0372  0.2869; 0.6861  0.7071; 0.6233  0.6245; 0.6344  0.6170];
>>d = [0.8587; 0.1781; 0.0747; 0.8405];
>> [Cd, lsqnonneg(C,d)]
ans =
   -2.5627         0
    3.1108    0.6929
注意:1。当问题为无约束线性最小二乘问题时,使用MATLAB下的“”运算即可以解决。2.对于非负最小二乘问题,调用lsqnonneg(C,d)求解。
5.8  非线性方程(组)求解
5.8.1  非线性方程的解
非线性方程的标准形式为f(x)=0
函数  fzero
格式  x = fzero (fun,x0)   %用fun定义表达式f(x),x0为初始解。
x = fzero (fun,x0,options)
[x,fval] = fzero(…)     %fval=f(x)
[x,fval,exitflag] = fzero(…)
[x,fval,exitflag,output] = fzero(…)
说明  该函数采用数值解求方程f(x)=0的根。
例5-19  求 的根
解:>> fun='x^3-2*x-5';
>> z=fzero(fun,2)    %初始估计值为2
结果为
z =
       2.0946
5.8.2  非线性方程组的解
非线性方程组的标准形式为:F(x) = 0
其中:x为向量,F(x)为函数向量。
函数  fsolve
格式  x = fsolve(fun,x0)   %用fun定义向量函数,其定义方式为:先定义方程函数function F = myfun (x)。
F =[表达式1;表达式2;…表达式m]   %保存为myfun.m,并用下面方式调用:x = fsolve(@myfun,x0),x0为初始估计值。
x = fsolve(fun,x0,options)
[x,fval] = fsolve(…)     %fval=F(x),即函数值向量
[x,fval,exitflag] = fsolve(…)
[x,fval,exitflag,output] = fsolve(…)
[x,fval,exitflag,output,jacobian] = fsolve(…)   % jacobian为解x处的Jacobian阵。
其余参数与前面参数相似。
例5-20   求下列系统的根
 
解:化为标准形式
 
设初值点为x0=[-5  -5]。
先建立方程函数文件,并保存为myfun.m:
function F = myfun(x)
F = [2*x(1) - x(2) - exp(-x(1));
      -x(1) + 2*x(2) - exp(-x(2))];
然后调用优化程序
x0 = [-5; -5];           % 初始点
options=optimset('Display','iter');   % 显示输出信息
[x,fval] = fsolve(@myfun,x0,options)
结果为
                                      Norm of    First-order
Iteration   Func-count       f(x)         step        optimality   CG-iterations
    1          4        47071.2             1    2.29e+004         0
    2          7        6527.47        1.45207    3.09e+003      1
    3         10        918.372        1.49186          418     1
    4         13         127.74        1.55326         57.3       1
    5         16        14.9153        1.57591         8.26      1
    6         19       0.779051        1.27662         1.14       1
    7         22     0.00372453       0.484658       0.0683       1
    8         25    9.21617e-008      0.0385552     0.000336       1
    9         28    5.66133e-017    0.000193707    8.34e-009      1
Optimization terminated successfully:
 Relative function value changing by less than OPTIONS.TolFun
x =
    0.5671
    0.5671
fval =
  1.0e-008 *
   -0.5320
   -0.5320
例5-21  求矩阵x使其满足方程 ,并设初始解向量为x=[1, 1; 1, 1]。
解:先编写M文件:
function F = myfun(x)
F = x*x*x-[1,2;3,4];
然后调用优化程序求解:
>>x0 = ones(2,2);     %初始解向量
>>options = optimset('Display','off');    %不显示优化信息
>>[x,Fval,exitflag] = fsolve(@myfun,x0,options)
则结果为
x =
   -0.1291    0.8602
    1.2903    1.1612
Fval =
  1.0e-003 *
    0.1541   -0.1163
    0.0109   -0.0243
exitflag =
     1
  评论这张
 
阅读(230)| 评论(0)
推荐 转载

历史上的今天

评论

<#--最新日志,群博日志--> <#--推荐日志--> <#--引用记录--> <#--博主推荐--> <#--随机阅读--> <#--首页推荐--> <#--历史上的今天--> <#--被推荐日志--> <#--上一篇,下一篇--> <#-- 热度 --> <#-- 网易新闻广告 --> <#--右边模块结构--> <#--评论模块结构--> <#--引用模块结构--> <#--博主发起的投票-->
 
 
 
 
 
 
 
 
 
 
 
 
 
 

页脚

网易公司版权所有 ©1997-2017