MATLAB優(yōu)化應(yīng)用非線性規(guī)劃(共22頁(yè))

1、MATLAB優(yōu)化應(yīng)用 1 線性規(guī)劃模型 一、線性規(guī)劃課題： 實(shí)例1：生產(chǎn)計(jì)劃問(wèn)題 假設(shè)某廠計(jì)劃生產(chǎn)甲、乙兩種產(chǎn)品，現(xiàn)庫(kù)存主要材料有A類3600公斤，B類2000公斤，C類3000公斤。每件甲產(chǎn)品需用材料A類9公斤，B類4公斤，C類3公斤。每件乙產(chǎn)品，需用材料A類4公斤，B類5公斤，C類10公斤。甲單位產(chǎn)品的利潤(rùn)70元，乙單位產(chǎn)品的利潤(rùn)120元。問(wèn)如何安排生產(chǎn)，才能使該廠所獲的利潤(rùn)最大。 建立數(shù)學(xué)模型： 設(shè)x1、x2分別為生產(chǎn)甲、乙產(chǎn)品的件數(shù)。f為該廠所獲總潤(rùn)。 max f=70x1+120x2 s.t 9x1+4x23600 4x1+5x22000 3x1+10x23000 x1,x20 實(shí)例

2、2：投資問(wèn)題 某公司有一批資金用于4個(gè)工程項(xiàng)目的投資，其投資各項(xiàng)目時(shí)所得的凈收益(投入資金锪百分比)如下表： 工程項(xiàng)目收益表 工程項(xiàng)目 A B C D 收益(%) 15 10 8 12 由于某種原因，決定用于項(xiàng)目A的投資不大于其他各項(xiàng)投資之和而用于項(xiàng)目B和C的投資要大于項(xiàng)目D的投資。試確定全文該公司收益最大的投資分配方案。 建立數(shù)學(xué)模型： 設(shè)x1、 x2 、x3 、x4分別代表用于項(xiàng)目A、B、C、D的投資百分?jǐn)?shù)。 max f=0.15x1+0.1x2+0.08 x3+0.12 x4 s.t x1-x2- x3- x40 x2+ x3- x40 x1+x2+x3+ x4=1 xj0 j=1,2,

3、3,4 實(shí)例3：運(yùn)輸問(wèn)題 有A、B、C三個(gè)食品加工廠，負(fù)責(zé)供給甲、乙、丙、丁四個(gè)市場(chǎng)。三個(gè)廠每天生產(chǎn)食品箱數(shù)上限如下表： 工廠 A B C 生產(chǎn)數(shù) 60 40 50 四個(gè)市場(chǎng)每天的需求量如下表： 市場(chǎng) 甲 乙 丙 丁 需求量 20 35 33 34 從各廠運(yùn)到各市場(chǎng)的運(yùn)輸費(fèi)(元/每箱)由下表給出： 收 點(diǎn) 發(fā) 點(diǎn) 市 場(chǎng) 甲 乙 丙 丁 工 廠 A 2 1 3 2 B 1 3 2 1 C 3 4 1 1 求在基本滿足供需平衡的約束條件下使總運(yùn)輸費(fèi)用最小。 建立數(shù)學(xué)模型： 設(shè)ai j為由工廠i運(yùn)到市場(chǎng)j的費(fèi)用，xi j 是由工廠i運(yùn)到市場(chǎng)j的箱數(shù)。bi是工廠i的產(chǎn)量，dj是市場(chǎng)j的需求量。 b=

4、 ( 60 40 50 ) d= ( 20 35 33 34 ) s.t x i j0 當(dāng)我們用MATLAB軟件作優(yōu)化問(wèn)題時(shí)，所有求maxf 的問(wèn)題化為求min(-f )來(lái)作。約束g i (x)0，化為 -g i0來(lái)作。 上述實(shí)例去掉實(shí)際背景，歸結(jié)出規(guī)劃問(wèn)題：目標(biāo)函數(shù)和約束條件都是變量x的線性函數(shù)。 形如： (1) min f T X s.t A Xb Aeq X =beq lbXub 其中X為n維未知向量，f T=f1,f2,fn為目標(biāo)函數(shù)系數(shù)向量，小于等于約束系數(shù)矩陣A為mn矩陣，b為其右端m維列向量，Aeq為等式約束系數(shù)矩陣，beq為等式約束右端常數(shù)列向量。lb,ub為自變量取值上界與下

5、界約束的n維常數(shù)向量。 二線性規(guī)劃問(wèn)題求最優(yōu)解函數(shù)： 調(diào)用格式： x=linprog(f,A,b) x=linprog(f,A,b,Aeq,beq) x=linprog(f,A,b,Aeq,beq,lb,ub) x=linprog(f,A,b,Aeq,beq,lb,ub,x0) x=linprog(f,A,b,Aeq,beq,lb,ub,x0,options) x,fval=linprog() x, fval, exitflag=linprog() x, fval, exitflag, output=linprog() x, fval, exitflag, output, lambda=lin

6、prog() 說(shuō)明：x=linprog(f,A,b)返回值x為最優(yōu)解向量。 x=linprog(f,A,b,Aeq,beq) 作有等式約束的問(wèn)題。若沒有不等式約束，則令A(yù)= 、b= 。 x=linprog(f,A,b,Aeq,beq,lb,ub,x0,options) 中l(wèi)b ,ub為變量x的下界和上界，x0為初值點(diǎn)，options為指定優(yōu)化參數(shù)進(jìn)行最小化。 Options的參數(shù)描述：Display 顯示水平。 選擇off 不顯示輸出；選擇iter顯示每一 步迭代過(guò)程的輸出；選擇final 顯示最終結(jié)果。 MaxFunEvals 函數(shù)評(píng)價(jià)的最大允許次數(shù) Maxiter 最大允許迭代次數(shù) Tol

7、X x處的終止容限 x,fval=linprog() 左端 fval 返回解x處的目標(biāo)函數(shù)值。 x,fval,exitflag,output,lambda=linprog(f,A,b, Aeq,beq,lb,ub,x0) 的輸出部分： exitflag 描述函數(shù)計(jì)算的退出條件：若為正值，表示目標(biāo)函數(shù)收斂于解x處；若為負(fù)值，表示目標(biāo)函數(shù)不收斂；若為零值，表示已經(jīng)達(dá)到函數(shù)評(píng)價(jià)或迭代的最大次數(shù)。 output 返回優(yōu)化信息：output.iterations表示迭代次數(shù)；output.algorithm表示所采用的算法；outprt.funcCount表示函數(shù)評(píng)價(jià)次數(shù)。 lambda 返回x處的拉格

8、朗日乘子。它有以下屬性： lambda.lower-lambda的下界； lambda.upper-lambda的上界； lambda.ineqlin-lambda的線性不等式； lambda.eqlin-lambda的線性等式。 三 舉例 例1：求解線性規(guī)劃問(wèn)題： max f=2x1+5x2 s.t 先將目標(biāo)函數(shù)轉(zhuǎn)化成最小值問(wèn)題：min(-f)=- 2x1-5x2 程序： f=-2 -5; A=1 0;0 1;1 2; b=4;3;8; x,fval=linprog(f,A,b) f=fval*(-1) 結(jié)果： x = 2 3 fval = -19.0000 maxf = 19 例2：min

9、f=5x1-x2+2x3+3x4-8x5 s.t -2x1+x2-x3+x4-3x56 2x1+x2-x3+4x4+x57 0xj15 j=1,2,3,4,5 程序： f=5 -1 2 3 -8; A=-2 1 -1 1 -3;2 1 -1 4 1; b=6;7; lb=0 0 0 0 0; ub=15 15 15 15 15; x,fval=linprog(f,A,b,lb,ub) 結(jié)果：x = 0.0000 0.0000 8.0000 0.0000 15.0000 minf = -104 例3：求解線性規(guī)劃問(wèn)題： minf=5x1+x2+2x3+3x4+x5 s.t -2x1+x2-x3+

10、x4-3x51 2x1+3x2-x3+2x4+x5-2 0xj1 j=1,2,3,4,5 程序： f=5 1 2 3 1; A=-2 1 -1 1 -3;2 3 -1 2 1; b=1;-2; lb=0 0 0 0 0; ub=1 1 1 1 1; x,fval,exitflag,output,lambda=linprog(f,A,b,lb,ub) 運(yùn)行結(jié)果： Exiting: One or more of the residuals, duality gap, or total relative error has grown 100000 times greater than its mi

11、nimum value so far: the primal appears to be infeasible (and the dual unbounded). (The dual residual 1 g(1)=8*x(1)+5*x(2); g(2)=5*x(1)+4*x(2); end 通過(guò)下面將優(yōu)化選項(xiàng)結(jié)構(gòu)options.GradObj設(shè)置為on來(lái)得到梯度值。 options=optimset(Gradobj,on); x0=1,1; x,fval,exitflag=fminunc(ff3,x0,options) 結(jié)果： x = 1.0e-015 * -0.2220 -0.2220 f

12、val = 5.4234e-031 exitflag = 1 2 minsearch函數(shù) 調(diào)用格式： x=fminsearch(fun,x0) x=fminsearch(fun,x0,options) x=fminsearch(fun,x0,options,P1,P2) x,fval=fminsearch() x,fval, exitflag=fminsearch() x,fval, exitflag,output=fminsearch() x,fval, exitflag,output,grad=fminsearch() x,fval, exitflag,output,grad,hessia

13、n=fminsearch() 說(shuō)明：參數(shù)及返回變量同上一函數(shù)。對(duì)求解二次以上的問(wèn)題，fminsearch函數(shù)比f(wàn)minunc函數(shù)有效。 3 多元非線性最小二乘問(wèn)題： 非線線性最小二乘問(wèn)題的數(shù)學(xué)模型為： 其中L為常數(shù)。 調(diào)用格式： x=lsqnonlin(fun,x0) x=lsqnonlin(fun,x0,lb,ub) x=lsqnonlin(fun,x0,options) x=lsqnonlin(fun,x0,options,P1,P2) x,resnorm=lsqnonlin() x,resnorm, residual,exitflag=lsqnonlin() x,resnorm, res

14、idual , exitflag,output=lsqnonlin() x,resnorm, residual,exitflag, output,lambda=lsqnonlin() x,resnorm, r esidual,exitflag, output,lambda,jacobian=lsqnonlin() 說(shuō)明：x返回解向量；resnorm返回x處殘差的平方范數(shù)值：sum(fun(x).2)；residual返回x處的殘差值fun(x)；lambda返回包含x處拉格朗日乘子的結(jié)構(gòu)參數(shù)；jacobian返回解x處的fun函數(shù)的雅可比矩陣。 lsqnonlin默認(rèn)時(shí)選擇大型優(yōu)化算法。Lsq

15、nonlin通過(guò)將options.LargeScale設(shè)置為off來(lái)作中型優(yōu)化算法。其采用一維搜索法。 例4求 minf=4(x2-x1)2+(x2-4)2 ，選擇初始點(diǎn)x0(1,1) 程序： f =4*(x(2)-x(1)2+(x(2)-4)2 x,reshorm=lsqnonlin(f,1,1) 結(jié)果： x = 3.9896 3.9912 reshorm = 5.0037e-009 例5：求 ，選擇初始點(diǎn)x0(0.2,0.3) 求解：先編輯ff5.m文件： function f=ff5(x) k=1:10; f=2+2*k-exp(k*x(1)-exp(k*x(2); 然后作程序：x0=0

16、.2,0.3; x,resnorm=lsqnonlin(ff5,x0) 結(jié)果 ： x = 0.2578 0.2578 resnorm = 124.3622 二 有約束非線性規(guī)劃問(wèn)題： 數(shù)學(xué)模型： min F(x) s.t Gi (x) 0 i=1,m Gj (x) =0 j=m+1,n xlxxu 其中：F(x)為多元實(shí)值函數(shù)，G(x)為向量值函數(shù)， 在有約束非線性規(guī)劃問(wèn)題中，通常要將該問(wèn)題轉(zhuǎn)換為更簡(jiǎn)單的子問(wèn)題，這些子問(wèn)題可以求并作為迭代過(guò)程的基礎(chǔ)。其基于K-T方程解的方法。它的K-T方程可表達(dá)為： 方程第一行描述了目標(biāo)函數(shù)和約束條件在解處梯度的取消。由于梯度取消，需要用拉格朗日乘子i來(lái)平衡目

17、標(biāo)函數(shù)與約束梯度間大小的差異。 調(diào)用格式: x=fmincon(f,x0,A,b) x=fmincon(f,x0,A,b,Aeq,beq) x=fmincon(f,x0,A,b,Aeq,beq,lb,ub) x=fmincon(f,x0,A,b,Aeq,beq,lb,ub,nonlcon) x=fmincon(f,x0,A,b,Aeq,beq,lb,ub,nonlcon,options) x,fval=fmincon() x, fval, exitflag=fmincon() x, fval, exitflag, output=fmincon() x, fval, exitflag, outp

18、ut, lambda=fmincon() 說(shuō)明：x=fmincon(f,x0,A,b)返回值x為最優(yōu)解向量。其中：x0為初始點(diǎn)。A,b為不等式約束的系數(shù)矩陣和右端列向量。 x=fmincon(f,x0,A,b,Aeq,beq) 作有等式約束的問(wèn)題。若沒有不等式約束，則令A(yù)= 、b= 。 x=fmincon(f, x0,A,b,Aeq,beq,lb,ub, nonlcon ,options) 中l(wèi)b ,ub為變量x的下界和上界；nonlcon=fun,由M文件fun.m給定非線性不等式約束c (x) 0和等式約束g(x)=0；options為指定優(yōu)化參數(shù)進(jìn)行最小化。 例6：求解：min 100(

19、x2-x12 )2+(1-x1)2 s.t x12; x22 程序：首先建立ff6.m文件： function f=ff6(x) f=100*(x(2)-x(2)2)2+(1-x(1)2; 然后在工作空間鍵入程序： x0=1.1,1.1; A=1 0;0 1; b=2;2; x,fval=fmincon(ff6,x0,A,b) 結(jié)果： x = 1.0000 1.0000 fval = 3.1936e-011 例7：求解： 首先建立目標(biāo)函數(shù)文件ff7.m文件： function f=ff7(x) f=-x(1)*x(2)*x(3) 然后將約束條件改寫成如下不等式： -x1-2x2-2x30 x1

20、+2x2+2x372 在工作空間鍵入程序： A=-1 -2 -2;1 2 2; b=0;72; x0=10;10;10; x,fval=fmincon(ff71,x0,A,b) 結(jié)果： x = 24.0000 12.0000 12.0000 fval = -3456 例8求解：minf=ex1(6x12+3x22+2x1x2+4x2+1) s.t x1x2-x1-x2+10 -2x1x2-50 程序：首先建立目標(biāo)函數(shù)文件ff8.m文件： function f=ff8(x) f=exp(x(1)*(6*x(1)2+3*x(2)2+2*x(1)*x(2)+4*x(2)+1); 再建立非線性的約束條

21、件文件：ff8g.m function c,g=ff8g(x) c(1)=x(1)*x(2)-x(1)-x(2)+1； c(2)=-2*x(1)*x(2)-5； g=; 然后在工作空間鍵入程序： x0=1,1; nonlcon=ff8g x, fval =fmincon(ff8,x0, nonlcon) 結(jié)果： x = -2.5000 1.0000 fval = 3.3244 exitflag = 1 當(dāng)有等式約束時(shí)，要放在矩陣g的位置，如上例中加等式約束： x(1)+2*x(1)=0 程序：首先建立 fun1.m文件： functionc,g=ff8g1(x) c(1)=x(1)*x(2)-

22、x(1)-x(2)+1； c(2)=-2*x(1)*x(2)-5； g(1)=x(1)+2*x(2); 然后在工作空間鍵入程序： x0=-1,1; nonlcon=ff8g1; x, fval,exitflag =fmincon(ff8,x0, nonlcon) 結(jié)果： x = -2.2361 1.1180 fval = 3.6576 exitflag = 1 3 二次規(guī)劃模型 數(shù)學(xué)模型： 其中H為二次型矩陣，A、Aeq分別為不等式約束與等式約束系數(shù)矩陣，f,b,beq,lb,ub,x為向量。 求解二次規(guī)劃問(wèn)題函數(shù)為quadprog( ) 調(diào)用格式： X= quadprog(H,f,A,b)

23、X= quadprog(H,f,A,b,Aeq,beq) X= quadprog(H,f,A,b,Aeq,beq,lb,ub) X= quadprog(H,f,A,b,Aeq,beq,lb,ub,x0) X= quadprog(H,f,A,b,Aeq,beq,lb,ub,x0,options) x,fval= quadprog() x,fval,exitflag= quadprog() x,fval,exitflag,output= quadprog() x,fval,exitflag,output,lambda= quadprog() 說(shuō)明：輸入?yún)?shù)中，x0為初始點(diǎn)；若無(wú)等式約束或無(wú)不等式約

24、束，就將相應(yīng)的矩陣和向量設(shè)置為空；options為指定優(yōu)化參數(shù)。輸出參數(shù)中，x是返回最優(yōu)解；fval是返回解所對(duì)應(yīng)的目標(biāo)函數(shù)值；exitflag是描述搜索是否收斂；output是返回包含優(yōu)化信息的結(jié)構(gòu)。Lambda是返回解x入包含拉格朗日乘子的參數(shù)。 例1：求解：二次規(guī)劃問(wèn)題 min f(x)= x1-3x2+3x12+4x22-2x1x2 s.t 2x1+x22 -x1+4x23 程序： f=1;-3 H=6 -2;-2 8 A=2 1;-1 4 b=2;3 X,fval,exitflag=quadprog(H,f,A,b) 結(jié)果： X = -0.0455 0.3636 fval = -0.

25、5682 exitflag = 1 例2：求解：二次規(guī)劃問(wèn)題 min +x12+2x22-2x1x2-4x1-12x2 s.t x1+x22 -x1+2x22 2x1+x23 0x1, 0x2 程序： H=2 -2;-2 4; f=-4;-12; A=1 1;-1 2;2 1; b=2;2;3; lb=zeros(2,1); x,fval,exitflag=quadprog(H,f,A,b,lb) 結(jié)果： x = 0.6667 1.3333 fval = -16.4444 exitflag = 1 4 多目標(biāo)規(guī)劃模型 多目標(biāo)規(guī)劃定義為在一組約束下，多個(gè)不同的目標(biāo)函數(shù)進(jìn)行優(yōu)化設(shè)計(jì)。 數(shù)學(xué)模型：

26、s.t gj (x) 0 j=1, 2, ,k 其中x=(x1 ,x2 , ,xn)為一個(gè)n維向量；fi(x)為目標(biāo)函數(shù)，i=1, 2, ,m; gj (x)為系統(tǒng)約束, j=1, 2, ,k。 當(dāng)目標(biāo)函數(shù)處于沖突狀態(tài)時(shí)，不存在最優(yōu)解使所有目標(biāo)函數(shù)同時(shí)達(dá)到最優(yōu)。于是我們尋求有效解(又稱非劣解或非支配解或帕累托解) 定義：若 ( )的鄰域內(nèi)不存在x，使得( +x)，且 則稱 為有效解。 多目標(biāo)規(guī)劃問(wèn)題的幾種常用解法： (1) 主要目標(biāo)法 其基本思想是：在多目標(biāo)問(wèn)題中，根據(jù)問(wèn)題的實(shí)際情況，確定一個(gè)目標(biāo)為主要目標(biāo)，而把其余目標(biāo)作為次要目標(biāo)，并且根據(jù)經(jīng)驗(yàn)，選取一定的界限值。這樣就可以把次要目標(biāo)作為約束

27、來(lái)處理，于是就將原來(lái)的多目標(biāo)問(wèn)題轉(zhuǎn)化為一個(gè)在新的約束下的單目標(biāo)最優(yōu)化問(wèn)題。 (2) 線性加權(quán)和法 其基本思想是：按照多目標(biāo)fi(x) (i=1, 2, ,m)的重要程度，分別乘以一組權(quán)系數(shù)j(j=1, 2, ,m)然后相加作為目標(biāo)函數(shù)而構(gòu)成單目標(biāo)規(guī)劃問(wèn)題。即 ，其中 例1：某鋼鐵廠準(zhǔn)備用5000萬(wàn)用于A、B兩個(gè)項(xiàng)目的技術(shù)改造投資。設(shè)x1、x2分別表示分配給項(xiàng)目A、B的投資。據(jù)專家預(yù)估計(jì)，投資項(xiàng)目A、B的年收益分別為70%和66%。同時(shí)，投資后總的風(fēng)險(xiǎn)損失將隨著總投資和單項(xiàng)投資的增加而增加，已知總的風(fēng)險(xiǎn)損失為0.02x12+0.01x22+0.04(x1+x2)2，問(wèn)應(yīng)如何分配資金才能使期望的收

28、益最大，同時(shí)使風(fēng)險(xiǎn)損失為最小。 建立數(shù)學(xué)模型 max f1(x)=70x1+66x2 min f2(x)= 0.02x12+0.01x22+0.04(x1+x2)2 s.t x1+x25000 0x1, 0x2 線性加權(quán)構(gòu)造目標(biāo)函數(shù)： max f=0.5f1(x) -0.5f2(x) 化最小值問(wèn)題： min (-f)=- 0.5f1(x) +0.5f2(x) 首先編輯目標(biāo)函數(shù)M文件ff11.m function f=ff11(x) f=-0.5*(70*x(1)+66*x(2)+0.5*(0.02*x(1)2+0.01*x(2)2+0.04*(x(1)+x(2)2); 調(diào)用單目標(biāo)規(guī)劃求最小值問(wèn)

29、題的函數(shù) x0=1000,1000 A=1 1; b=5000; lb=zeros(2,1); x,fval, exitflag=fmincon(ff11,x0, A,b,lb,) f1=70*x(1)+66*x(2) f2=0.02*x(1)2+0.01*x(2)2+0.04*(x(1)+x(2)2 結(jié)果：x = 307.1428 414.2857 fval = -1.2211e+004 exitflag = 1 f1 = 4.8843e+004 f2 = 2.4421e+004 (3) 極大極小法 其基本思想是：對(duì)于極小化的多目標(biāo)規(guī)劃，讓其中最大的目標(biāo)函數(shù)值盡可能地小為此，對(duì)每個(gè) xR，我

30、們先求諸目標(biāo)函數(shù)值fi(x)的最大值，然后再求這些最大值中的最小值。即構(gòu)造單目標(biāo)規(guī)劃： (4) 目標(biāo)達(dá)到法 對(duì)于多目標(biāo)規(guī)劃： s.t gj (x) 0 j=1, 2, ,n 先設(shè)計(jì)與目標(biāo)函數(shù)相應(yīng)的一組目標(biāo)值理想化向量 ， 再設(shè)為一松弛因子標(biāo)量。設(shè) 為權(quán)值系數(shù)向量。 于是多目標(biāo)規(guī)劃問(wèn)題化為： 在Matlab的優(yōu)化工具箱中，fgoalattain函數(shù)用于解決此類問(wèn)題。 其數(shù)學(xué)模型形式為： min F(x)-weight goal c(x) 0 ceq(x)=0 A xb Aeq x=beq lbxub 其中，x,weight,goal,b,beq,lb和ub為向量，A和Aeq為矩陣，c(x),ceq(x)和F(x)為函數(shù)， 調(diào)用格式： x=fgoalattain(F,x0,goal,weight) x=fgoalattain(F,x0,goal,weight,A,b) x=fgoalattain(F,x0,goal,weight,A,b,Aeq,beq) x=fgoalattain(F,x0,goal,weight,A,b,Aeq,beq,lb,ub) x=fgoalattain(F,x0,goal,wei