- 相關(guān)推薦
matlab線性規(guī)劃
MATLAB 優(yōu)化問(wèn)題
1 線性規(guī)劃問(wèn)題
線性規(guī)劃問(wèn)題是目標(biāo)函數(shù)和約束條件均為線性函數(shù)的問(wèn)題,MATLAB6.0解決的線性規(guī)劃問(wèn)題的標(biāo)準(zhǔn)形式為:
min f(x)x?R
sub.to:A?x?b
Aeq?x?beq
lb?x?ub n
其中f、x、b、beq、lb、ub為向量,A、Aeq為矩陣。
其它形式的線性規(guī)劃問(wèn)題都可經(jīng)過(guò)適當(dāng)變換化為此標(biāo)準(zhǔn)形式。
在MATLAB6.0版中,線性規(guī)劃問(wèn)題(Linear Programming)已用函數(shù)linprog取代了函數(shù) linprog
格式 x = linprog(f,A,b) %求min f ' *x sub.to A?x?b線性規(guī)劃的最優(yōu)解。
x = linprog(f,A,b,Aeq,beq) %等式約束Aeq?x?beq,若沒有不等式約束
A?x?b,則A=[ ],b=[ ]。
x = linprog(f,A,b,Aeq,beq,lb,ub) %指定x的范圍lb?x?ub,若沒有等式約束
Aeq?x?beq ,則Aeq=[ ],beq=[ ]
x = linprog(f,A,b,Aeq,beq,lb,ub,x0) %設(shè)置初值x0
x = linprog(f,A,b,Aeq,beq,lb,ub,x0,options) % options為指定的優(yōu)化參數(shù)
[x,fval] = linprog(?) % 返回目標(biāo)函數(shù)最優(yōu)值,即fval= f ' *x。
[x,lambda,exitflag] = linprog(?) % lambda為解x的Lagrange乘子。
[x, lambda,fval,exitflag] = linprog(?) % exitflag為終止迭代的錯(cuò)誤條件。
[x,fval, lambda,exitflag,output] = linprog(?) % output為關(guān)于優(yōu)化的一些信息
說(shuō)明 若exitflag>0表示函數(shù)收斂于解x,exitflag=0表示超過(guò)函數(shù)估值或迭代的最大數(shù)字,exitflag
例5-1 求下面的優(yōu)化問(wèn)題
min ?5x1?4x2?6x3
sub.to x1?x2?x3?20
3x1?2x2?4x3?42
3x1?2x2?30
.
0?x1,0?x2,0?x3
解:
>>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)
結(jié)果為:
x = %最優(yōu)解
0.0000
15.0000
3.0000
fval = %最優(yōu)值
-78.0000
exitflag = %收斂
1
output =
iterations: 6 %迭代次數(shù)
cgiterations: 0
algorithm: 'lipsol' %所使用規(guī)則
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個(gè)下界是有效的
2 非線性規(guī)劃問(wèn)題
2.1 有約束的一元函數(shù)的最小值
單變量函數(shù)求最小值的標(biāo)準(zhǔn)形式為minf(x) sub.to x1?x?x2 x
在MATLAB中使用fmin函數(shù)求其最小值。
函數(shù) fminbnd
格式 x = fminbnd(fun,x1,x2) %返回自變量x在區(qū)間x1?x?x2上函數(shù)fun取最小值
時(shí)x值,fun為目標(biāo)函數(shù)的表達(dá)式字符串或MATLAB
自定義函數(shù)的函數(shù)柄。
x = fminbnd(fun,x1,x2,options) % options為指定優(yōu)化參數(shù)選項(xiàng)
[x,fval] = fminbnd(?) % fval為目標(biāo)函數(shù)的最小值
[x,fval,exitflag] = fminbnd(?) %xitflag為終止迭代的條件
[x,fval,exitflag,output] = fminbnd(?) % output為優(yōu)化信息
說(shuō)明 若參數(shù)exitflag>0,表示函數(shù)收斂于x,若exitflag=0,表示超過(guò)函數(shù)估計(jì)值或迭代的最大數(shù)字,exitflag
例5-2 計(jì)算下面函數(shù)在區(qū)間(0,1)內(nèi)的最小值。
x3?cosx?xlogx f(x)?e解:>> [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]上求下面函數(shù)的最小值f(x)?(x?3)3?1
解:先自定義函數(shù):在MATLAB編輯器中建立M文件為:
function f = myfun(x)
f = (x-3).^2 - 1;
保存為myfun.m,然后在命令窗口鍵入命令:
>> x=fminbnd(@myfun,0,5)
則結(jié)果顯示為:
x =
3
2.2 無(wú)約束多元函數(shù)最小值
多元函數(shù)最小值的標(biāo)準(zhǔn)形式為minf(x) x
其中:x為向量,如x?[x1,x2,?,xn]
在MATLAB中使用fmins求其最小值。
命令 利用函數(shù)fminsearch求無(wú)約束多元函數(shù)最小值
函數(shù) fminsearch
格式 x = fminsearch(fun,x0) %x0為初始點(diǎn),fun為目標(biāo)函數(shù)的表達(dá)式字符串或
MATLAB自定義函數(shù)的函數(shù)柄。
x = fminsearch(fun,x0,options) % options查optimset
.
[x,fval] = fminsearch(?) %最優(yōu)點(diǎn)的函數(shù)值
[x,fval,exitflag] = fminsearch(?) % exitflag與單變量情形一致
[x,fval,exitflag,output] = fminsearch(?) %output與單變量情形一致
注意:fminsearch采用了Nelder-Mead型簡(jiǎn)單搜尋法。
32例5-4 求y?2x1?4x1x32?10x1x2?x2的最小值點(diǎn)
解:>>X=fminsearch('2*x(1)^3+4*x(1)*x(2)^3-10*x(1)*x(2)+x(2)^2', [0,0])
結(jié)果為
X =
1.0016 0.8335
或在MATLAB編輯器中建立函數(shù)文件
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])
結(jié)果為:
X =
1.0016 0.8335
命令 利用函數(shù)fminunc求多變量無(wú)約束函數(shù)最小值
函數(shù) fminunc
格式 x = fminunc(fun,x0) %返回給定初始點(diǎn)x0的最小函數(shù)值點(diǎn)
x = fminunc(fun,x0,options) % options為指定優(yōu)化參數(shù)
[x,fval] = fminunc(?) %fval最優(yōu)點(diǎn)x處的函數(shù)值
[x,fval,exitflag] = fminunc(?) % exitflag為終止迭代的條件,與上同。
[x,fval,exitflag,output] = fminunc(?) %output為輸出優(yōu)化信息
[x,fval,exitflag,output,grad] = fminunc(?) % grad為函數(shù)在解x處的梯度值
[x,fval,exitflag,output,grad,hessian] = fminunc(?) %目標(biāo)函數(shù)在解x處的海賽
(Hessian)值
注意:當(dāng)函數(shù)的階數(shù)大于2時(shí),使用fminunc比f(wàn)minsearch更有效,但當(dāng)所選函數(shù)高度不連續(xù)時(shí),使用fminsearch效果較好。
2?2x1x2?x2例5-5 求f(x)?3x12的最小值。
>> 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)
結(jié)果為:
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
2.3 有約束的多元函數(shù)最小值
非線性有約束的多元函數(shù)的標(biāo)準(zhǔn)形式為:
minf(x) x
sub.to C(x)?0
Ceq(x)?0
A?x?b
Aeq?x?beq
lb?x?ub
其中:x、b、beq、lb、ub是向量,A、Aeq為矩陣,C(x)、Ceq(x)是返回向量的函數(shù),f(x)為目標(biāo)函數(shù),f(x)、C(x)、Ceq(x)可以是非線性函數(shù)。
在MATLAB5.x中,它的求解由函數(shù)constr實(shí)現(xiàn)。
函數(shù) 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(?)
參數(shù)說(shuō)明:fun為目標(biāo)函數(shù),它可用前面的方法定義;
x0為初始值;
A、b滿足線性不等式約束A?x?b,若沒有不等式約束,則取A=[ ],b=[ ];
Aeq、beq滿足等式約束Aeq?x?beq,若沒有,則取Aeq=[ ],beq=[ ];
lb、ub滿足lb?x?ub,若沒有界,可設(shè)lb=[ ],ub=[ ];
nonlcon的作用是通過(guò)接受的向量x來(lái)計(jì)算非線性不等約束C(x)?0和等式
約束Ceq(x)?0分別在x處的估計(jì)C和Ceq,通過(guò)指定函數(shù)柄來(lái)使用,
如:>>x = fmincon(@myfun,x0,A,b,Aeq,beq,lb,ub,@mycon),先建立非
線性約束函數(shù),并保存為mycon.m:function [C,Ceq] = mycon(x)
C = ? % 計(jì)算x處的非線性不等約束C(x)?0的函數(shù)值。
Ceq = ? % 計(jì)算x處的非線性等式約束Ceq(x)?0的函數(shù)值。
lambda是Lagrange乘子,它體現(xiàn)哪一個(gè)約束有效。
output輸出優(yōu)化信息;
grad表示目標(biāo)函數(shù)在x處的梯度;
hessian表示目標(biāo)函數(shù)在x處的Hessiab值。
例5-6 求下面問(wèn)題在初始點(diǎn)(0,1)處的最優(yōu)解
2min x1?x22?x1x2?2x1?5x2
sub.to ?(x1?1)2?x2?0
2x1?3x2?6?0
解:約束條件的標(biāo)準(zhǔn)形式為
sub.to (x1?1)2?x2?0
?2x1?3x2?6
先在MATLAB編輯器中建立非線性約束函數(shù)文件:
function [c, ceq]=mycon (x)
c=(x(1)-1)^2-x(2);
ceq=[ ]; %無(wú)等式約束
然后,在命令窗口鍵入如下命令或建立M文件:
>>fun='x(1)^2+x(2)^2-x(1)*x(2)-2*x(1)-5*x(2)'; %目標(biāo)函數(shù)
>>x0=[0 1];
>>A=[-2 3]; %線性不等式約束
>>b=6;
>>Aeq=[ ]; %無(wú)線性等式約束
>>beq=[ ];
>>lb=[ ]; %x沒有下、上界
>>ub=[ ];
>>[x,fval,exitflag,output,lambda,grad,hessian]
=fmincon(fun,x0,A,b,Aeq,beq,lb,ub,@mycon)
則結(jié)果為
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下界有效情況,通過(guò)lambda.lower可查看。
upper: [2x1 double] %x上界有效情況,為0表示約束無(wú)效。
eqlin: [0x1 double] %線性等式約束有效情況,不為0表示約束有效。
eqnonlin: [0x1 double] %非線性等式約束有效情況。
ineqlin: 2.5081e-008 %線性不等式約束有效情況。
ineqnonlin: 6.1938e-008 %非線性不等式約束有效情況。
grad = %目標(biāo)函數(shù)在最小值點(diǎn)的梯度
1.0e-006 *
-0.1776
hessian = %目標(biāo)函數(shù)在最小值點(diǎn)的Hessian值
1.0000 -0.0000
-0.0000 1.0000
例5-7 求下面問(wèn)題在初始點(diǎn)x=(10, 10, 10)處的最優(yōu)解。
Min f(x)??x1x2x3
Sub.to 0?x1?2x2?2x3?72
解:約束條件的標(biāo)準(zhǔn)形式為
sub.to ?x1?2x2?2x3?0 x1?2x2?2x3?72
>> 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)
結(jié)果為:
x =
24.0000 12.0000 12.0000
fval =
-3456
2.4 二次規(guī)劃問(wèn)題
二次規(guī)劃問(wèn)題(quadratic programming)的標(biāo)準(zhǔn)形式為:
minx?Hx?f?x sub.to A?x?b
?x?beq Aeq
lb?x?ub
其中,H、A、Aeq為矩陣,f、b、beq、lb、ub、x為向量
.
其它形式的二次規(guī)劃問(wèn)題都可轉(zhuǎn)化為標(biāo)準(zhǔn)形式。
MATLAB中二次規(guī)劃應(yīng)用函數(shù)quadprog。
函數(shù) quadprog
格式 x = quadprog(H,f,A,b) %其中H,f,A,b為標(biāo)準(zhǔn)形中的參數(shù),x為目標(biāo)函數(shù)的最小
值。
x = quadprog(H,f,A,b,Aeq,beq) %Aeq,beq滿足等約束條件Aeq?x?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為設(shè)置的初值
x = quadprog(H,f,A,b,Aeq,beq,lb,ub,x0,options) % options為指定的優(yōu)化參數(shù)
[x,fval] = quadprog(?) %fval為目標(biāo)函數(shù)最優(yōu)值
[x,fval,exitflag] = quadprog(?) % exitflag與線性規(guī)劃中參數(shù)意義相同
[x,fval,exitflag,output] = quadprog(?) % output與線性規(guī)劃中參數(shù)意義相同
[x,fval,exitflag,output,lambda] = quadprog(?) % lambda與線性規(guī)劃中參數(shù)意義
相同 例5-8 求解下面二次規(guī)劃問(wèn)題
min
sub.to 2f(x)?x1?x22?x1x2?2x1?6x2 x1?x2?2
?x1?2x2?2
2x1?x2?3
0?x1,0?x2
解:f(x)?x?Hx?f?x ?1?1??x1???2?f?x?則H??,,???6??x? ?12???2???
在MATLAB中實(shí)現(xiàn)如下:
>>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)
結(jié)果為:
x = %最優(yōu)解
0.6667
1.3333
fval = %最優(yōu)值
-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
>> lambda.lower
ans =
說(shuō)明 第1、2個(gè)約束條件有效,其余無(wú)效。
例5-9 求二次規(guī)劃的最優(yōu)解
max f (x1, x2)=x1x2+3
sub.to x1+x2-2=0
解:化成標(biāo)準(zhǔn)形式:
?0?1??x1??x1?minf(x1x2)??x1x2?3?(x1x2)?????(0,0)???3 ??10??x2??x2?
sub.to x1+x2=2
在Matlab中實(shí)現(xiàn)如下:
>>H=[0,-1;-1,0];
>>f=[0;0];
>>Aeq=[1 1];
>>b=2;
>>[x,fval,exitflag,output,lambda] = quadprog(H,f,[ ],[ ],Aeq,b)
結(jié)果為:
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: [ ]
.
3 “半無(wú)限”有約束的多元函數(shù)最優(yōu)解
“半無(wú)限”有約束多元函數(shù)最優(yōu)解問(wèn)題的標(biāo)準(zhǔn)形式為
minxf(x)
sub.to C(x)?0
Ceq(x)?0
A?x?b
Aeq?x?beq
K1(x,w1)?0
K2(x,w2)?0
?
Kn(x,wn)?0
其中:x、b、beq、lb、ub都是向量;A、Aeq是矩陣;C(x)、Ceq(x)、Ki(x,wi)是返回向量的函數(shù),f(x)為目標(biāo)函數(shù);f(x)、C(x)、Ceq(x)是非線性函數(shù);Ki(x,wi)為半無(wú)限約束,w1,w2,?,wn通常是長(zhǎng)度為2的向量。
在MTALAB5.x中,使用函數(shù)seminf解決這類問(wèn)題。
函數(shù) 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(?)
參數(shù)說(shuō)明:x0為初始估計(jì)值;
fun為目標(biāo)函數(shù),其定義方式與前面相同;
A、b由線性不等式約束A?x?b確定,沒有,則A=[ ],b=[ ];
Aeq、beq由線性等式約束Aeq?x?beq確定,沒有,則Aeq=[ ],beq=[ ]; Lb、ub由變量x的范圍lb?x?ub確定;
options為優(yōu)化參數(shù);
ntheta為半無(wú)限約束的個(gè)數(shù);
seminfcon用來(lái)確定非線性約束向量C和Ceq以及半無(wú)限約束的向量K1,
K2,?,Kn,通過(guò)指定函數(shù)柄來(lái)使用,如:
x = fseminf(@myfun,x0,ntheta,@myinfcon)
先建立非線性約束和半無(wú)限約束函數(shù)文件,并保存為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 = ? %計(jì)算樣本集
w2 = ? %計(jì)算樣本集
?
wntheta = ? % 計(jì)算樣本集
K1 = ? % 在x和w處的第1個(gè)半無(wú)限約束值
K2 = ? %在x和w處的第2個(gè)半無(wú)限約束值
?
Kntheta = ? %在x和w處的第ntheta個(gè)半無(wú)限約束值
C = ? % 在x處計(jì)算非線性不等式約束值
Ceq = ? % 在x處計(jì)算非線性等式約束值
如果沒有約束,則相應(yīng)的值取為“[ ]”,如Ceq=[]
fval為在x處的目標(biāo)函數(shù)最小值;
exitflag為終止迭代的條件;
output為輸出的優(yōu)化信息;
lambda為解x的Lagrange乘子。
例5-10 求下面一維情形的最優(yōu)化問(wèn)題
minxf(x)?(x1?0.5)2?(x2?0.5)2?(x3?0.5)2
sub.to
K1(x,w1)?sin(w1x1)cos(w1x2)?(w1?50)2?sin(w1x3)?x3?1 K2(x,w2)?sin(w2x2)cos(w2x1)?(w2?50)2?sin(w2x3)?x3?1 1?w1?100
1?w2?100
解:將約束方程化為標(biāo)準(zhǔn)形式:
K1(x,w1)?sin(w1x1)cos(w1x2)?(w1?50)2?sin(w1x3)?x3?1?0 K2(x,w2)?sin(w2x2)cos(w2x1)?(w2?50)2?sin(w2x3)?x3?1?0 先建立非線性約束和半無(wú)限約束函數(shù)文件,并保存為mycon.m:
function [C,Ceq,K1,K2,S] = mycon(X,S)
% 初始化樣本間距:
.
if isnan(S(1,1)),
S = [0.2 0; 0.2 0];
end
% 產(chǎn)生樣本集:
w1 = 1:S(1,1):100;
w2 = 1:S(2,1):100;
% 計(jì)算半無(wú)限約束:
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;
% 無(wú)非線性約束:
C = [ ]; Ceq=[ ];
% 繪制半無(wú)限約束圖形
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)
結(jié)果為:
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 求下面二維情形的最優(yōu)化問(wèn)題
minxf(x)?(x1?0.2)2?(x2?0.2)2?(x3?0.2)2
sub.to
K1(x,w)?sin(w1x1)cos(w2x2)?(w1?50)2?sin(w1x3)?x3??
sin(w2x2)cos(w1x1)?(w2?50)2?sin(w2x3)?x3?1.5 1?w1?100
1?w2?100
初始點(diǎn)為x0=[0.25, 0.25, 0.25]。
解:先建立非線性和半無(wú)限約束函數(shù)文件,并保存為mycon.m:
function [C,Ceq,K1,S] = mycon(X,S)
% 初始化樣本間距:
if isnan(s(1,1)),
s = [2 2];
end
% 設(shè)置樣本集
w1x = 1:s(1,1):100;
w1y = 1:s(1,2):100;
[wx, wy] = meshgrid(w1x,w1y);
% 計(jì)算半無(wú)限約束函數(shù)值
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;
% 無(wú)非線性約束
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)
結(jié)果為(如圖)
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-2
4 極小化極大(Minmax)問(wèn)題
極小化極大問(wèn)題的標(biāo)準(zhǔn)形式為
minmaxx{Fi}{Fi(x)}
sub.to C(x)?0
Ceq(x)?0
A?x?b
.
Aeq?x?beq
lb?x?ub
其中:x、b、beq、lb、ub是向量,A、Aeq為矩陣,C(x)、Ceq(x)和F(x)是返回向量的函數(shù),F(xiàn)(x)、C(x)、Ceq(x)可以是非線性函數(shù)。
在MATLAB5.x中,它的求解由函數(shù)minmax實(shí)現(xiàn)。
函數(shù) 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(?)
參數(shù)說(shuō)明:fun為目標(biāo)函數(shù);
x0為初始值;
A、b滿足線性不等約束A?x?b,若沒有不等約束,則取A=[ ],b=[ ];
Aeq、beq滿足等式約束Aeq?x?beq,若沒有,則取Aeq=[ ],beq=[ ];
lb、ub滿足lb?x?ub,若沒有界,可設(shè)lb=[ ],ub=[ ];
nonlcon的作用是通過(guò)接受的向量x來(lái)計(jì)算非線性不等約束C(x)?0和等式約束
Ceq(x)?0分別在x處的值C和Ceq,通過(guò)指定函數(shù)柄來(lái)使用,如:>>x =
fminimax(@myfun,x0,A,b,Aeq,beq,lb,ub,@mycon),先建立非線性約束函數(shù),
并保存為mycon.m:function [C,Ceq] = mycon(x)
C = ? % 計(jì)算x處的非線性不等約束C(x)?0的函數(shù)值。
Ceq = ? % 計(jì)算x處的非線性等式約束Ceq(x)?0的函數(shù)值。
options為指定的優(yōu)化參數(shù);
fval為最優(yōu)點(diǎn)處的目標(biāo)函數(shù)值;
maxfval為目標(biāo)函數(shù)在x處的最大值;
exitflag為終止迭代的條件;
lambda是Lagrange乘子,它體現(xiàn)哪一個(gè)約束有效。
output輸出優(yōu)化信息。
例5-12 求下列函數(shù)最大值的最小化問(wèn)題
[f1(x), f2(x), f3(x), f4(x), f5(x)]
2?x2其中:f1(x)?2x12?48x1?40x2?304
2f2(x)??x22?3x2
f3(x)?x1?3x2?18
f4(x)??x1?x2
f5(x)?x1?x2?8
解:先建立目標(biāo)函數(shù)文件,并保存為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)
結(jié)果為:
x =
4.0000
4.0000
fval =
0.0000 -64.0000 -2.0000 -8.0000 -0.0000
例5-13 求上述問(wèn)題的絕對(duì)值的最大值最小化問(wèn)題。
目標(biāo)函數(shù)為:[|f1(x)|, |f2(x)|, |f3(x)|, |f4(x)|, |f5(x)|]
解:先建立目標(biāo)函數(shù)文件(與上例相同)
然后,在命令窗口或編輯器中建立M文件:
>>x0 = [0.1; 0.1]; % 初始點(diǎn)
>>options = optimset('MinAbsMax',5); % 指定絕對(duì)值的最小化
>>[x,fval] = fminimax(@myfun,x0,[ ],[ ],[ ],[ ],[ ],[ ],[ ],options)
則結(jié)果為
x =
4.9256
2.0796
fval =
37.2356 -37.2356 -6.8357 -7.0052 -0.9948
5 多目標(biāo)規(guī)劃問(wèn)題
多目標(biāo)規(guī)劃是指在一組約束下,對(duì)多個(gè)不同目標(biāo)函數(shù)進(jìn)行優(yōu)化。它的一般形式為
min[f1(x),f2(x),?,fm(x)]
j?1,2,?,p sub.to gj(x)?0
其中:x?(x1,x2,?,xn)。
在同一約束下,當(dāng)目標(biāo)函數(shù)處于沖突狀態(tài)時(shí),不存在最優(yōu)解x使所有目標(biāo)函數(shù)同時(shí)達(dá)到最優(yōu)。此時(shí),我們使用有效解,即如果不存在x?S,使得fi(x)?fi(x*),i=1,2,?m, 則稱x*為有效解。
在MATLAB中,多目標(biāo)問(wèn)題的標(biāo)準(zhǔn)形式為
.
minimize? x,?
sub.to F(x)?weight???goal
C(x)?0
Ceq(x)?0
A?x?b
Aeq?x?beq
lb?x?ub
其中:x、b、beq、lb、ub是向量;A、Aeq為矩陣;C(x)、Ceq(x)和F(x)是返回向量的函數(shù);F(x)、C(x)、Ceq(x)可以是非線性函數(shù);weight為權(quán)值系數(shù)向量,用于控制對(duì)應(yīng)的目標(biāo)函數(shù)與用戶定義的目標(biāo)函數(shù)值的接近程度;goal為用戶設(shè)計(jì)的與目標(biāo)函數(shù)相應(yīng)的目標(biāo)函數(shù)值向量;?為一個(gè)松弛因子標(biāo)量;F(x)為多目標(biāo)規(guī)劃中的目標(biāo)函數(shù)向量。
在MATLAB5.x中,它的最優(yōu)解由attgoal函數(shù)實(shí)現(xiàn)。
函數(shù) 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(?)
參數(shù)說(shuō)明:
x0為初始解向量;
fun為多目標(biāo)函數(shù)的文件名字符串,其定義方式與前面fun的定義方式相同;
goal為用戶設(shè)計(jì)的目標(biāo)函數(shù)值向量;
weight為權(quán)值系數(shù)向量,用于控制目標(biāo)函數(shù)與用戶自定義目標(biāo)值的接近程度;
A、b滿足線性不等式約束A?x?b,沒有時(shí)取A=[ ],b=[ ];
Aeq、beq滿足線性等式約束Aeq?x?beq,沒有時(shí)取Aeq=[ ],beq=[ ];
lb、ub為變量的下界和上界:lb?x?ub;
nonlcon的作用是通過(guò)接受的向量x來(lái)計(jì)算非線性不等約束C(x)?0和等式約束
Ceq(x)?0分別在x處的值C和Ceq,通過(guò)指定函數(shù)柄來(lái)使用。
如:>>x = fgoalattain(@myfun,x0,goal,wei
http://www.ishadingyu.com/news/559EE48D172564E3.html ght,A,b,Aeq,beq,lb,ub,@mycon),先建立非線性約束函數(shù),并保存為mycon.m:function [C,Ceq] = mycon(x)
C = ? % 計(jì)算x處的非線性不等式約束C(x)?0的函數(shù)值。
Ceq = ? % 計(jì)算x處的非線性等式約束Ceq(x)?0的函數(shù)值。
options為指定的優(yōu)化參數(shù);
fval為多目標(biāo)函數(shù)在x處的值;
attainfactor為解x處的目標(biāo)規(guī)劃因子;
exitflag為終止迭代的條件;
output為輸出的優(yōu)化信息;
lambda為解x處的Lagrange乘子
例5-14 控制系統(tǒng)輸出反饋器設(shè)計(jì)。
設(shè)如下線性系統(tǒng)
??Ax?Bu x
y?Cx
0???o.50?10??100??210? B???22? C??其中:A??0? ????001????1?2??0??01??
要求設(shè)計(jì)輸出反饋控制器K,使閉環(huán)系統(tǒng)
??(A?BKC)x?Bu x
y?Cx
在復(fù)平面實(shí)軸上點(diǎn)[-5,-3,-1]的左側(cè)有極點(diǎn),并要求 ?4?Kij?4(i,j?1,2)
解:上述問(wèn)題就是要求解矩陣K,使矩陣(A+BKC)的極點(diǎn)為[-5,-3,-1],這是一個(gè)多目標(biāo)規(guī)劃問(wèn)題。
先建立目標(biāo)函數(shù)文件,保存為eigfun.m:
function F = eigfun(K,A,B,C)
F = sort(eig(A+B*K*C)); % 估計(jì)目標(biāo)函數(shù)值
然后,輸入?yún)?shù)并調(diào)用優(yōu)化程序:
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]; % 為閉合環(huán)路的特征值(極點(diǎn))設(shè)置目標(biāo)值向量
weight = abs(goal) % 設(shè)置權(quán)值向量
lb = -4*ones(size(K0)); % 設(shè)置控制器的下界
ub = 4*ones(size(K0)); % 設(shè)置控制器的上界
options = optimset('Display','iter'); % 設(shè)置顯示參數(shù):顯示每次迭代的輸出
[K,fval,attainfactor] = fgoalattain(@eigfun,K0,goal,weight,[],[],[],[],lb,ub,[],options,A,B,C)
結(jié)果為:
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
6 最小二乘最優(yōu)問(wèn)題
6.1 約束線性最小二乘
有約束線性最小二乘的標(biāo)準(zhǔn)形式為
minCx?dxsub.to A?x?b
Aeq?x?beq
lb?x?ub 22
其中:C、A、Aeq為矩陣;d、b、beq、lb、ub、x是向量。
在MATLAB5.x中,約束線性最小二乘用函數(shù)conls求解。
函數(shù) lsqlin
格式 x = lsqlin(C,d,A,b) %求在約束條件A?x?b下,方程Cx = d的最小二乘解x。
x = lsqlin(C,d,A,b,Aeq,beq) %Aeq、beq滿足等式約束Aeq?x?beq,若沒有不
等式約束,則設(shè)A=[ ],b=[ ]。
x = lsqlin(C,d,A,b,Aeq,beq,lb,ub) %lb、ub滿足lb?x?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為指定優(yōu)化參數(shù)
[x,resnorm] = lsqlin(?) % resnorm=norm(C*x-d)^2,即2-范數(shù)。
[x,resnorm,residual] = lsqlin(?) %residual=C*x-d,即殘差。
[x,resnorm,residual,exitflag] = lsqlin(?) %exitflag為終止迭代的條件
[x,resnorm,residual,exitflag,output] = lsqlin(?) % output表示輸出優(yōu)化信息
[x,resnorm,residual,exitflag,output,lambda] = lsqlin(?) % lambda為解x的
Lagrange乘子
例5-15 求解下面系統(tǒng)的最小二乘解
系統(tǒng):Cx=d
約束:A?x?b;lb?x?ub
先輸入系統(tǒng)系數(shù)和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);
然后調(diào)用最小二乘命令:
[x,resnorm,residual,exitflag,output,lambda] = lsqlin(C,d,A,b,[ ],[ ],lb,ub);
結(jié)果為:
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 %說(shuō)明解x是收斂的
output =
iterations: 4
algorithm: 'medium-scale: active-set'
firstorderopt: []
cgiterations: []
lambda =
lower: [4x1 double]
upper: [4x1 double]
.
eqlin: [0x1 double]
ineqlin: [3x1 double]
通過(guò)lambda.ineqlin可查看非線性不等式約束是否有效。
6.2 非線性數(shù)據(jù)(曲線)擬合
非線性曲線擬合是已知輸入向量xdata和輸出向量ydata,并且知道輸入與輸出的函數(shù)關(guān)系為ydata=F(x, xdata),但不知道系數(shù)向量x。今進(jìn)行曲線擬合,求x使得下式成立:
minF(x,xdata)?ydatax2
2??(F(x,xdatai)?ydatai)2 i
在MATLAB5.x中,使用函數(shù)curvefit解決這類問(wèn)題。
函數(shù) 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(?)
參數(shù)說(shuō)明:
x0為初始解向量;xdata,ydata為滿足關(guān)系ydata=F(x, xdata)的數(shù)據(jù);
lb、ub為解向量的下界和上界lb?x?ub,若沒有指定界,則lb=[ ],ub=[ ]; options為指定的優(yōu)化參數(shù);
fun為擬合函數(shù),其定義方式為:x = lsqcurvefit(@myfun,x0,xdata,ydata),
其中myfun已定義為 function F = myfun(x,xdata)
F = ? % 計(jì)算x處擬合函數(shù)值fun的用法與前面相同;
resnorm=sum ((fun(x,xdata)-ydata).^2),即在x處殘差的平方和;
residual=fun(x,xdata)-ydata,即在x處的殘差;
exitflag為終止迭代的條件;
output為輸出的優(yōu)化信息;
lambda為解x處的Lagrange乘子;
jacobian為解x處擬合函數(shù)fun的jacobian矩陣。
例5-16 求解如下最小二乘非線性擬合問(wèn)題
已知輸入向量xdata和輸出向量ydata,且長(zhǎng)度都是n,擬合函數(shù)為
ydata(i)?x(1)?xdata(i)2?x(2)?sin(xdata(i))?x(3)?xdata(i)3
即目標(biāo)函數(shù)為min?(F(x,xdatai)?ydatai)2 xi?1
其中:F(x,xdata)?x(1)?xdata2?x(2)?sin(xdata)?x(3)?xdata3
n
初始解向量為x0=[0.3, 0.4, 0.1]。
解:先建立擬合函數(shù)文件,并保存為myfun.m
function F = myfun(x,xdata)
F = x(1)*xdata.^2 + x(2)*sin(xdata) + x(3)*xdata.^3;
然后給出數(shù)據(jù)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]; %初始估計(jì)值
>>[x,resnorm] = lsqcurvefit(@myfun,x0,xdata,ydata)
結(jié)果為:
Optimization terminated successfully:
Relative function value changing by less than OPTIONS.TolFun
x =
0.2269 0.3385 0.3021
resnorm =
6.2950
6.3 非線性最小二乘
非線性最小二乘(非線性數(shù)據(jù)擬合)的標(biāo)準(zhǔn)形式為
minf(x)?f1(x)2?f2(x)2???fm(x)2?L x
其中:L為常數(shù)
在MATLAB5.x中,用函數(shù)leastsq解決這類問(wèn)題,在6.0版中使用函數(shù)lsqnonlin。
?f1(x)??f(x)?2? 設(shè)F(x)???????f(x)m??
則目標(biāo)函數(shù)可表達(dá)為minF(x)x2
2??fi(x)2 i
其中:x為向量,F(xiàn)(x)為函數(shù)向量。
函數(shù) lsqnonlin
格式 x = lsqnonlin(fun,x0) %x0為初始解向量;fun為fi(x),i=1,2,?,m,fun返回向
量值F,而不是平方和值,平方和隱含在算法中,fun的定義與前面相同。
x = lsqnonlin(fun,x0,lb,ub) %lb、ub定義x的下界和上界:lb?x?ub。
x = lsqnonlin(fun,x0,lb,ub,options) %options為指定優(yōu)化參數(shù),若x沒有界,則
lb=[ ],ub=[ ]。
[x,resnorm] = lsqnonlin(?) % resnorm=sum(fun(x).^2),即解x處目標(biāo)函數(shù)值。
[x,resnorm,residual] = lsqnonlin(?) % residual=fun(x),即解x處fun的值。
[x,resnorm,residual,exitflag] = lsqnonlin(?) %exitflag為終止迭代條件。
[x,resnorm,residual,exitflag,output] = lsqnonlin(?) %output輸出優(yōu)化信息。
[x,resnorm,residual,exitflag,output,lambda] = lsqnonlin(?) %lambda為L(zhǎng)agrage
.
乘子。
[x,resnorm,residual,exitflag,output,lambda,jacobian] =lsqnonlin(?) %fun在解x
處的Jacobian矩陣。
例5-17 求下面非線性最小二乘問(wèn)題?(2?2k?ekx1?ekx2)2初始解向量為x0=[0.3,
k?110
0.4]。
解:先建立函數(shù)文件,并保存為myfun.m,由于lsqnonlin中的fun為向量形式而不是平方和形式,因此,myfun函數(shù)應(yīng)由fi(x)建立:
fk(x)?2?2k?ekx1?ekx2 k=1,2,…,10
function F = myfun(x)
k = 1:10;
F = 2 + 2*k-exp(k*x(1))-exp(k*x(2));
然后調(diào)用優(yōu)化程序:
x0 = [0.3 0.4];
[x,resnorm] = lsqnonlin(@myfun,x0)
結(jié)果為:
Optimization terminated successfully:
Norm of the current step is less than OPTIONS.TolX
x =
0.2578 0.2578
resnorm = %求目標(biāo)函數(shù)值
124.3622
6.4 非負(fù)線性最小二乘
非負(fù)線性最小二乘的標(biāo)準(zhǔn)形式為:
minCx?dxsub.to x?0 22
其中:矩陣C和向量d為目標(biāo)函數(shù)的系數(shù),向量x為非負(fù)獨(dú)立變量。
在MATLAB5.x中,用函數(shù)nnls求解這類問(wèn)題,在6.0版中則用函數(shù)lsqnonneg。 函數(shù) lsqnonneg
格式 x = lsqnonneg(C,d) %C為實(shí)矩陣,d為實(shí)向量
x = lsqnonneg(C,d,x0) % x0為初始值且大于0
x = lsqnonneg(C,d,x0,options) % options為指定優(yōu)化參數(shù)
[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 一個(gè)最小二乘問(wèn)題的無(wú)約束與非負(fù)約束解法的比較。
先輸入數(shù)據(jù):
>>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];
>> [C\d, lsqnonneg(C,d)]
ans =
-2.5627 0
3.1108 0.6929
注意:1。當(dāng)問(wèn)題為無(wú)約束線性最小二乘問(wèn)題時(shí),使用MATLAB下的“\”運(yùn)算即可以解決。2.對(duì)于非負(fù)最小二乘問(wèn)題,調(diào)用lsqnonneg(C,d)求解。
7 非線性方程(組)求解
7.1 非線性方程的解
非線性方程的標(biāo)準(zhǔn)形式為f(x)=0
函數(shù) fzero
格式 x = fzero (fun,x0) %用fun定義表達(dá)式f(x),x0為初始解。
x = fzero (fun,x0,options)
[x,fval] = fzero(?) %fval=f(x)
[x,fval,exitflag] = fzero(?)
[x,fval,exitflag,output] = fzero(?)
說(shuō)明 該函數(shù)采用數(shù)值解求方程f(x)=0的根。
例5-19 求x3?2x?5?0的根
解:>> fun='x^3-2*x-5';
>> z=fzero(fun,2) %初始估計(jì)值為2
結(jié)果為
z =
2.0946
7.2 非線性方程組的解
非線性方程組的標(biāo)準(zhǔn)形式為:F(x) = 0
其中:x為向量,F(xiàn)(x)為函數(shù)向量。
函數(shù) fsolve
格式 x = fsolve(fun,x0) %用fun定義向量函數(shù),其定義方式為:先定義方程函數(shù)
function F = myfun (x)。
F =[表達(dá)式1;表達(dá)式2;?表達(dá)式m] %保存為myfun.m,并用下面方式調(diào)用:
x = fsolve(@myfun,x0),x0為初始估計(jì)值。
x = fsolve(fun,x0,options)
[x,fval] = fsolve(?) %fval=F(x),即函數(shù)值向量
[x,fval,exitflag] = fsolve(?)
[x,fval,exitflag,output] = fsolve(?)
.
[x,fval,exitflag,output,jacobian] = fsolve(?) % jacobian為解x處的Jacobian陣。 其余參數(shù)與前面參數(shù)相似。
例5-20 求下列系統(tǒng)的根
2x1?x2?e?x1
?x1?2x2?e?x2
解:化為標(biāo)準(zhǔn)形式
2x1?x2?e?x1?0
?x1?2x2?e?x2?0
設(shè)初值點(diǎn)為x0=[-5 -5]。
先建立方程函數(shù)文件,并保存為myfun.m:
function F = myfun(x)
F = [2*x(1) - x(2) - exp(-x(1));
-x(1) + 2*x(2) - exp(-x(2))];
然后調(diào)用優(yōu)化程序
x0 = [-5; -5]; % 初始點(diǎn)
options=optimset('Display','iter'); % 顯示輸出信息
[x,fval] = fsolve(@myfun,x0,options)
結(jié)果為
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
?12?例5-21 求矩陣x使其滿足方程x?x?x???,并設(shè)初始解向量為x=[1, 1; 1, 1]。 34??
解:先編寫M文件:
function F = myfun(x)
F = x*x*x-[1,2;3,4];
然后調(diào)用優(yōu)化程序求解:
>>x0 = ones(2,2); %初始解向量
>>options = optimset('Display','off'); %不顯示優(yōu)化信息 >>[x,Fval,exitflag] = fsolve(@myfun,x0,options) 則結(jié)果為
x =
-0.1291 0.8602
1.2903 1.1612
Fval =
1.0e-003 *
0.1541 -0.1163
0.0109
exitflag =
1
. -0.0243
【matlab線性規(guī)劃】相關(guān)文章:
Matlab用于發(fā)動(dòng)機(jī)故障預(yù)報(bào)及Delphi對(duì)Matlab的調(diào)用05-01
利用線性規(guī)劃思想解題04-29
基于Matlab的斜拋運(yùn)動(dòng)研究04-29
MATLAB環(huán)境下地面模型構(gòu)建04-27
基于MATLAB的ABS控制仿真研究04-27
基于matlab的物料大小分級(jí)算法的實(shí)現(xiàn)05-02
楊氏雙縫干涉實(shí)驗(yàn)的MatLab模擬04-29