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

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

2、2：投資問題 某公司有一批資金用于4個工程項目的投資，其投資各項目時所得的凈收益(投入資金锪百分比)如下表： 工程項目收益表 工程項目 A B C D 收益(%) 15 10 8 12 由于某種原因，決定用于項目A的投資不大于其他各項投資之和而用于項目B和C的投資要大于項目D的投資。試確定全文該公司收益最大的投資分配方案。 建立數(shù)學模型： 設x1、 x2 、x3 、x4分別代表用于項目A、B、C、D的投資百分數(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 實例3：運輸問題 有A、B、C三個食品加工廠，負責供給甲、乙、丙、丁四個市場。三個廠每天生產(chǎn)食品箱數(shù)上限如下表： 工廠 A B C 生產(chǎn)數(shù) 60 40 50 四個市場每天的需求量如下表： 市場 甲 乙 丙 丁 需求量 20 35 33 34 從各廠運到各市場的運輸費(元/每箱)由下表給出： 收 點 發(fā) 點 市 場 甲 乙 丙 丁 工 廠 A 2 1 3 2 B 1 3 2 1 C 3 4 1 1 求在基本滿足供需平衡的約束條件下使總運輸費用最小。 建立數(shù)學模型： 設ai j為由工廠i運到市場j的費用，xi j 是由工廠i運到市場j的箱數(shù)。bi是工廠i的產(chǎn)量，dj是市場j的需求量。 b=

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

5、界約束的n維常數(shù)向量。 二線性規(guī)劃問題求最優(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() 說明：x=linprog(f,A,b)返回值x為最優(yōu)解向量。 x=linprog(f,A,b,Aeq,beq) 作有等式約束的問題。若沒有不等式約束，則令A= 、b= 。 x=linprog(f,A,b,Aeq,beq,lb,ub,x0,options) 中l(wèi)b ,ub為變量x的下界和上界，x0為初值點，options為指定優(yōu)化參數(shù)進行最小化。 Options的參數(shù)描述：Display 顯示水平。 選擇off 不顯示輸出；選擇iter顯示每一 步迭代過程的輸出；選擇final 顯示最終結(jié)果。 MaxFunEvals 函數(shù)評價的最大允許次數(shù) Maxiter 最大允許迭代次數(shù) Tol

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

8、朗日乘子。它有以下屬性： lambda.lower-lambda的下界； lambda.upper-lambda的上界； lambda.ineqlin-lambda的線性不等式； lambda.eqlin-lambda的線性等式。 三 舉例 例1：求解線性規(guī)劃問題： max f=2x1+5x2 s.t 先將目標函數(shù)轉(zhuǎ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ī)劃問題： 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) 運行結(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 通過下面將優(yōu)化選項結(jié)構(gòu)options.GradObj設置為on來得到梯度值。 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() 說明：參數(shù)及返回變量同上一函數(shù)。對求解二次以上的問題，fminsearch函數(shù)比fminunc函數(shù)有效。 3 多元非線性最小二乘問題： 非線線性最小二乘問題的數(shù)學模型為： 其中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() 說明：x返回解向量；resnorm返回x處殘差的平方范數(shù)值：sum(fun(x).2)；residual返回x處的殘差值fun(x)；lambda返回包含x處拉格朗日乘子的結(jié)構(gòu)參數(shù)；jacobian返回解x處的fun函數(shù)的雅可比矩陣。 lsqnonlin默認時選擇大型優(yōu)化算法。Lsq

15、nonlin通過將options.LargeScale設置為off來作中型優(yōu)化算法。其采用一維搜索法。 例4求 minf=4(x2-x1)2+(x2-4)2 ，選擇初始點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：求 ，選擇初始點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ī)劃問題： 數(shù)學模型： min F(x) s.t Gi (x) 0 i=1,m Gj (x) =0 j=m+1,n xlxxu 其中：F(x)為多元實值函數(shù)，G(x)為向量值函數(shù)， 在有約束非線性規(guī)劃問題中，通常要將該問題轉(zhuǎn)換為更簡單的子問題，這些子問題可以求并作為迭代過程的基礎(chǔ)。其基于K-T方程解的方法。它的K-T方程可表達為： 方程第一行描述了目標函數(shù)和約束條件在解處梯度的取消。由于梯度取消，需要用拉格朗日乘子i來平衡目

17、標函數(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() 說明：x=fmincon(f,x0,A,b)返回值x為最優(yōu)解向量。其中：x0為初始點。A,b為不等式約束的系數(shù)矩陣和右端列向量。 x=fmincon(f,x0,A,b,Aeq,beq) 作有等式約束的問題。若沒有不等式約束，則令A= 、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ù)進行最小化。 例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：求解： 首先建立目標函數(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 程序：首先建立目標函數(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 當有等式約束時，要放在矩陣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ù)學模型： 其中H為二次型矩陣，A、Aeq分別為不等式約束與等式約束系數(shù)矩陣，f,b,beq,lb,ub,x為向量。 求解二次規(guī)劃問題函數(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() 說明：輸入?yún)?shù)中，x0為初始點；若無等式約束或無不等式約

24、束，就將相應的矩陣和向量設置為空；options為指定優(yōu)化參數(shù)。輸出參數(shù)中，x是返回最優(yōu)解；fval是返回解所對應的目標函數(shù)值；exitflag是描述搜索是否收斂；output是返回包含優(yōu)化信息的結(jié)構(gòu)。Lambda是返回解x入包含拉格朗日乘子的參數(shù)。 例1：求解：二次規(guī)劃問題 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ī)劃問題 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 多目標規(guī)劃模型 多目標規(guī)劃定義為在一組約束下，多個不同的目標函數(shù)進行優(yōu)化設計。 數(shù)學模型：

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

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

28、益最大，同時使風險損失為最小。 建立數(shù)學模型 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)造目標函數(shù)： max f=0.5f1(x) -0.5f2(x) 化最小值問題： min (-f)=- 0.5f1(x) +0.5f2(x) 首先編輯目標函數(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)用單目標規(guī)劃求最小值問

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) 極大極小法 其基本思想是：對于極小化的多目標規(guī)劃，讓其中最大的目標函數(shù)值盡可能地小為此，對每個 xR，我

30、們先求諸目標函數(shù)值fi(x)的最大值，然后再求這些最大值中的最小值。即構(gòu)造單目標規(guī)劃： (4) 目標達到法 對于多目標規(guī)劃： s.t gj (x) 0 j=1, 2, ,n 先設計與目標函數(shù)相應的一組目標值理想化向量 ， 再設為一松弛因子標量。設 為權(quán)值系數(shù)向量。 于是多目標規(guī)劃問題化為： 在Matlab的優(yōu)化工具箱中，fgoalattain函數(shù)用于解決此類問題。 其數(shù)學模型形式為： 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