數(shù)學(xué)建模最優(yōu)化方法建模及實現(xiàn)

1、 實驗實驗07 07 最優(yōu)化方法建模及實現(xiàn)最優(yōu)化方法建模及實現(xiàn) 實驗?zāi)康膶嶒災(zāi)康膶嶒瀮?nèi)容實驗內(nèi)容3、基于最優(yōu)化方法建模及實現(xiàn)、論文寫作。、基于最優(yōu)化方法建模及實現(xiàn)、論文寫作。1、了解最優(yōu)化問題的基本內(nèi)容。、了解最優(yōu)化問題的基本內(nèi)容。2、用數(shù)學(xué)軟件包、用數(shù)學(xué)軟件包matlab求解求解(非非)線性規(guī)劃問題。線性規(guī)劃問題。4、實驗題目：鋼管的訂購與運輸。、實驗題目：鋼管的訂購與運輸。1、基礎(chǔ)知識、例子。、基礎(chǔ)知識、例子。3、建模案例：投資的收益與風(fēng)險、建模案例：投資的收益與風(fēng)險2、掌握線性規(guī)劃及非線性規(guī)劃建模及其、掌握線性規(guī)劃及非線性規(guī)劃建模及其MATLAB實現(xiàn)。實現(xiàn)。 最優(yōu)化問題最優(yōu)化問題v優(yōu)化問

2、題，一般是指用“最好”的方式，使用或分配有限的資源，即勞動力、原材料、機器、資金等，使得費用最小或利潤最大.v建立優(yōu)化問題的數(shù)學(xué)模型 1) 確定問題的決策變量 2) 構(gòu)造模型的目標(biāo)函數(shù)和允許取值的范圍，常用一組不等式來表示. T1min(max)( ),(,)s.t.( )0,1,2,nizf xxxxgxim或（1）（2）由（1）、（2）組成的模型屬于約束優(yōu)化，若只有（1）式就是無約束優(yōu)化，f(x)稱為目標(biāo)函數(shù)，gi(x)稱為約束條件若目標(biāo)函數(shù)f(x)和約束條件g(x)都是線性函數(shù)，則稱該模型是線性規(guī)劃. 線性規(guī)劃模型線性規(guī)劃模型例1 、生產(chǎn)炊事用具需要兩種資源勞動力和原材料，某公司制定生產(chǎn)

3、計劃，生產(chǎn)三種不同的產(chǎn)品，生產(chǎn)管理部門提供的數(shù)據(jù)如下 A B C勞動力（小時/件）736原材料（千克/件）445利潤（元/件）423 每天供應(yīng)原材料200kg，每天可使用的勞動力為150h. 建立線性規(guī)劃模型，使總收益最大，并求各種產(chǎn)品的日產(chǎn)量. 解解 第一步,確定決策變量. 用 分別表示A, B, C三種產(chǎn)品的日產(chǎn)量 第二步, 約束條件 原材料： 勞動力： 第三步，確定目標(biāo)函數(shù) AxBxCx445200ABCxxx736150ABCxxx423ABCZxxx 例2 一家廣告公司想在電視、廣播上做廣告，其目的是盡可能多的招來顧客，下面是調(diào)查結(jié)果： 電視無線電 廣播雜志白天最佳時間一次廣告費用（

4、千元）40753015受每次廣告影響的顧客數(shù)（千人）400900500200受每次廣告影響的女顧客數(shù)（千人）300400200100 這家公司希望廣告費用不超過800（千元）還要求：1）至少要有200萬婦女收看廣告；2）電視廣告費用不超過500（千元） 3）電視廣告白天至少播出3次，最佳時間至少播出2次；4）通過廣播、雜志做的廣告要重復(fù)5到10次.123412341234121234m ax400900500200407530158003004002001002000. .40755003,2,510,510Zxxxxxxxxxxxxs txxxxxx令 分別白天，最佳電視、廣播、雜志廣告次數(shù)

5、 1234,xxxx， 例例3： 任務(wù)分配問題：某車間有甲、乙兩臺機床，可用于加工三種工件。假定這兩臺車床的可用臺時數(shù)分別為800和900，三種工件的數(shù)量分別為400、600和500，且已知用三種不同車床加工單位數(shù)量不同工件所需的臺時數(shù)和加工費用如下表。問怎樣分配車床的加工任務(wù)，才能既滿足加工工件的要求，又使加工費用最低？ 解解 設(shè)在甲車床上加工工件1、2、3的數(shù)量分別為x1、x2、x3，在乙車床上加工工件1、2、3的數(shù)量分別為x4、x5、x6?？山⒁韵戮€性規(guī)劃模型： 解答 例例4： 某廠每日8小時的產(chǎn)量不低于1800件。為了進(jìn)行質(zhì)量控制，計劃聘請兩種不同水平的檢驗員。一級檢驗員的標(biāo)準(zhǔn)為：速

6、度25件/小時，正確率98%，計時工資4元/小時；二級檢驗員的標(biāo)準(zhǔn)為：速度15小時/件，正確率95%，計時工資3元/小時。檢驗員每錯檢一次，工廠要損失2元。為使總檢驗費用最省，該工廠應(yīng)聘一級、二級檢驗員各幾名？解解 設(shè)需要一級和二級檢驗員的人數(shù)分別為x1、x2人,則應(yīng)付檢驗員的工資為：212124323848xxxx因檢驗員錯檢而造成的損失為：21211282)%5158%2258(xxxx 故目標(biāo)函數(shù)為：故目標(biāo)函數(shù)為：2121213640)128()2432(minxxxxxxz約束條件為：0, 0180015818002581800158258212121xxxxxx 線性規(guī)劃模型：線性規(guī)

7、劃模型：213640minxxz0, 01594535 . .212121xxxxxxts 解答返 回 線性規(guī)劃模型的一般形式線性規(guī)劃模型的一般形式.,.,.,.,.minnixnibxat sxcuinkikikniii2102111 目標(biāo)函數(shù)和所有的約束條件都是決策變量目標(biāo)函數(shù)和所有的約束條件都是決策變量的線性函數(shù)。的線性函數(shù)。min. . ucxAxbstvlbxvub矩矩陣陣形形式式： 實際問題中實際問題中的優(yōu)化模型的優(yōu)化模型mixgtsxxxxfzMaxMiniTn, 2 , 1, 0)(. .),(),()(1或x決策變量決策變量f(x)目標(biāo)函數(shù)目標(biāo)函數(shù)gi(x) 0約束條件約束條

8、件數(shù)學(xué)規(guī)劃數(shù)學(xué)規(guī)劃線性規(guī)劃線性規(guī)劃(LP)二次規(guī)劃二次規(guī)劃(QP)非線性規(guī)劃非線性規(guī)劃(NLP)純整數(shù)規(guī)劃純整數(shù)規(guī)劃(PIP)混合整數(shù)規(guī)劃混合整數(shù)規(guī)劃(MIP)整數(shù)規(guī)劃整數(shù)規(guī)劃(IP)0-1整數(shù)規(guī)劃整數(shù)規(guī)劃一般整數(shù)規(guī)劃一般整數(shù)規(guī)劃連續(xù)規(guī)劃連續(xù)規(guī)劃 優(yōu)化模型的分類優(yōu)化模型的分類 線性規(guī)劃問題的求解在理論上有單純形法，在實際建模中常用以下解法： 1. 圖解法 2. LINGO 軟件包； 3. Excel中的規(guī)劃求解； 4. MATLAB軟件包. min z=cX bAXts. .1、模型：命令：x=linprog（c，A，b） 2、模型：min z=cX bAXts. .beqXAeq命令：x=l

9、inprog（c，A，b，Aeq,beq） 或或 x=linprog（c，A，b，Aeq,beq,x0）或或 x,fval=linprog（c，A，b，Aeq,beq）注意注意：若沒有不等式： 存在，則令A(yù)= ，b= .bAX 用用MATLAB優(yōu)化工具箱解線性優(yōu)化工具箱解線性linear規(guī)劃規(guī)劃 3、模型：min z=cX bAXts. .beqXAeqVLBXVUB 命令：1 x=linprog（c，A，b，Aeq,beq, VLB，VUB） 2 x=linprog（c，A，b，Aeq,beq, VLB，VUB, X0） 注意：1 若沒有等式約束: , 則令A(yù)eq= , beq= . 2其中

10、X0表示初始點 beqXAeq4、命令：x,fval=linprog() 返回最優(yōu)解返回最優(yōu)解及及處的目標(biāo)函數(shù)值處的目標(biāo)函數(shù)值fval. 321436minxxxz 1231231232380120. .3005020 xxxxxxstxxx 解解: 編寫編寫M文件文件xxgh1.m如下：如下： c=6 3 4; A=1,2,-3;0 1 0; b=80;50; Aeq=1 1 1; beq=120; vlb=30,0,20; vub=; x,fval=linprog(c,A,b,Aeq,beq,vlb,vub)To Matlab (xxgh1)123123123123min(6 3 4)12

11、38001050. . 1 1 112030020 xzxxxxxxstxxxxx 例例5 解解 編寫編寫M文件文件xxgh2.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; Aeq=; beq=; vlb=0;0;0;0;0;0; vub=;x,fval=linprog(c,A,b,Aeq,beq,vlb,vub) To Matlab (

12、xxgh2)max 6543216 . 064. 072. 032. 028. 04 . 0 xxxxxxz 85003. 003. 003. 001. 001. 001. 0. .654321xxxxxxts 70005. 002. 041xx 10005. 002. 052xx 90008. 003. 063xx 6, 2 , 10jxj 例例6 S.t.Xz8121110913min 9008003 . 12 . 15 . 000000011 . 14 . 0X500600400100100010010001001X ，0654321xxxxxxX改寫為： 問題問題例例3的解答 編寫編寫

13、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;beq=400 600 500;vlb = zeros(6,1);vub=;x,fval = linprog(f,A,b,Aeq,beq,vlb,vub)To Matlab (xxgh3) x = 0.0000 600.0000 0.0000 400.0000 0.0000 500.0000fval =1.3800e+004計算結(jié)果：計算結(jié)果

14、： 即在甲機床上加工600個工件2,在乙機床上加工400個工件1、500個工件3，可在滿足條件的情況下使總加工費最小為13800。 問題問題 213640minxxz s.t. )45(3521xx改寫為：例例4的解答 編寫編寫M文件文件xxgh4.m如下：如下：c = 40;36;A=-5 -3;b=-45;Aeq=;beq=;vlb = zeros(2,1);vub=9;15; %調(diào)用linprog函數(shù)：x,fval = linprog(c,A,b,Aeq,beq,vlb,vub)To Matlab (xxgh4) 結(jié)果為：結(jié)果為： x = 9.0000 0.0000 fval =360即

15、只需聘用9個一級檢驗員。 注：注：本問題應(yīng)還有一個約束條件：x1、x2取整數(shù)。故它是一個整數(shù)線性規(guī)劃整數(shù)線性規(guī)劃問題。這里把它當(dāng)成一個線性規(guī)劃來解，求得其最優(yōu)解剛好是整數(shù)：x1=9，x2=0，故它就是該整數(shù)規(guī)劃的最優(yōu)解。若用線性規(guī)劃解法求得的最優(yōu)解不是整數(shù)，將其取整后不一定是相應(yīng)整數(shù)規(guī)劃的最優(yōu)解，這樣的整數(shù)規(guī)劃應(yīng)用專門的方法求解。返 回 1） 首先建立M文件fun.m,定義目標(biāo)函數(shù)F（X）:function f=fun(X);f=F(X); 其中X為n維變元向量，G(X)與Ceq(X)均為非線性函數(shù)組成的向量，其它變量的含義與線性規(guī)劃、二次規(guī)劃中相同.用Matlab求解上述問題，基本步驟分三步

16、：二、非線性規(guī)劃問題及其二、非線性規(guī)劃問題及其MatlabMatlab (). .()0()0MinF XAXbAeq XbeqstG XCeq XVLBXVUB 3） 建立主程序.非線性規(guī)劃求解的函數(shù)是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) (4) x=fmincon(fun,X0,A,b,Aeq,beq,VLB,VUB,nonlcon) (5) x=fmincon(fun,X0,A,b,

17、Aeq,beq,VLB,VUB,nonlcon,options) (6) x,fval= fmincon(.) (7) x,fval,exitflag= fmincon(.) (8) x,fval,exitflag,output= fmincon(.)輸出極值點M文件迭代的初值參數(shù)說明變量上下限 注意：注意：1 fmincon函數(shù)提供了大型優(yōu)化算法和中型優(yōu)化算法。默認(rèn)時，若在fun函數(shù)中提供了梯度（options參數(shù)的GradObj設(shè)置為on），并且只有上下界存在或只有等式約束，fmincon函數(shù)將選擇大型算法。當(dāng)既有等式約束又有梯度約束時，使用中型算法。2 fmincon函數(shù)的中型算法使用的

18、是序列二次規(guī)劃法。在每一步迭代中求解二次規(guī)劃子問題，并用BFGS法更新拉格朗日Hessian矩陣。3 fmincon函數(shù)可能會給出局部最優(yōu)解，這與初值X0的選取有關(guān)。 1.寫成標(biāo)準(zhǔn)形式寫成標(biāo)準(zhǔn)形式： s.t. 00546322121xxxx2100 xx22212121212minxxxxf22212121212minxxxxf 2x1+3x2 6 s.t x1+4x2 5 x1,x2 0例例1 2.先建立先建立M-文件文件 fun1.m: function f=fun1(x); f=-x(1)-2*x(2)+(1/2)*x(1)2+(1/2)*x(2)2MATLAB(youh1)3.再建立主

19、程序再建立主程序youh1.m： 4.運算結(jié)果為：運算結(jié)果為： x = 0.7647 1.0588 fval = -2.0294 1先建立先建立M文件文件 fun2.m,定義目標(biāo)函數(shù)定義目標(biāo)函數(shù): function f=fun2(x) f=exp(x(1)*(4*x(1)2+2*x(2)2+4*x(1)*x(2)+2*x(2)+1);) 12424()(22122211xxxxxexfx x1+x2=0 s.t. 1.5+x1x2 - x1 - x2 0 -x1x2 10 02再建立再建立M文件文件mycon2.m定義非線性約束：定義非線性約束： function g,ceq=mycon2(x

20、) g=1.5+x(1)*x(2)-x(1)-x(2);-x(1)*x(2)-10; ceq= ;例例2 3主程序主程序youh2.m為為:x0=-1;1;A=;b=;Aeq=1 1;beq=0;vlb=;vub=;x,fval=fmincon(fun2,x0,A,b,Aeq,beq,vlb,vub,mycon2)MATLAB(youh2)4. 運算結(jié)果為運算結(jié)果為： x = -1.2247 1.2247 fval = 1.8951 12221122221212min2250s.t.7005, 010fXxxgXxxgXxxxx 1先建立先建立M-文件文件fun3.m定義目標(biāo)函數(shù)定義目標(biāo)函數(shù):

21、 function f=fun3(x) f=-2*x(1)-x(2);2再建立再建立M文件文件mycon3.m定義非線性約束：定義非線性約束： function g,ceq=mycon3(x) g=x(1)2+x(2)2-25;x(1)2-x(2)2-7; ceq=;例例3 3. 主程序主程序youh3.m為為: x0=3;2.5; VLB=0 0;VUB=5 10; x,fval,exitflag,output =fmincon(fun3,x0,VLB,VUB,mycon3)MATLAB(youh3(fun3) 4. 運算結(jié)果為運算結(jié)果為: x = 4.0000 3.0000fval =-1

22、1.0000exitflag = 1output = iterations: 4 funcCount: 17 stepsize: 1 algorithm: 1x44 char firstorderopt: cgiterations: 返回返回 建模案例：投資的收益和風(fēng)險建模案例：投資的收益和風(fēng)險(1998A)(1998A) 二、基本假設(shè)和符號規(guī)定二、基本假設(shè)和符號規(guī)定 1.1.總體總體風(fēng)險用所投資的風(fēng)險用所投資的S Si i中最大的一個風(fēng)險來衡量，即：中最大的一個風(fēng)險來衡量，即： iiqxinmax|1,2, 三、模型的建立與分析三、模型的建立與分析 3.3.建立模型建立模型00101 2ma

23、x ()minmax (). ., ,niiiiiiniiiirp xq xp xMs txin 雙目標(biāo)模型雙目標(biāo)模型為：為： 4.4.模型簡化模型簡化00101 2max (). .(), ,niiiiiiniiiirp xq xaMs tp xMxin 即模型為即模型為： 四、模型四、模型1 1的求解的求解siri （%）qi （%）pi （%）ui （元）（元）S0（銀行）5000S1282.51103S2211.52198S3235.54.552S4252.66.540將將n=4,M=1,及平均收益率及平均收益率ri, 風(fēng)險損失率風(fēng)險損失率qi,費率費率 pi代入模型代入模型1得：得：

24、00101 2max (). .(), ,niiiiiiniiiirp xq xaMs tp xMxin 由于由于a是任意給定的風(fēng)險度，到底怎樣給定沒有一個準(zhǔn)是任意給定的風(fēng)險度，到底怎樣給定沒有一個準(zhǔn)則，不同的投資者有不同的風(fēng)險度。我們從則，不同的投資者有不同的風(fēng)險度。我們從a=0開始，以開始，以步長步長a=0.001進(jìn)行循環(huán)搜索，編制程序如下：進(jìn)行循環(huán)搜索，編制程序如下： a=0;while(1.1-a)1 c=-0.05 -0.27 -0.19 -0.185 -0.185; Aeq=1 1.01 1.02 1.045 1.065; beq=1; A=0 0.025 0 0 0;0 0 0.

25、015 0 0;0 0 0 0.055 0;0 0 0 0 0.026; b=a;a;a;a; vlb=0,0,0,0,0;vub=; x,val=linprog(c,A,b,Aeq,beq,vlb,vub); ax=x Q=-val plot(a,Q,.);axis(0 0.1 0 0.5);hold on a=a+0.001;end xlabel(a),ylabel(Q)To Matlab（xxgh5）模型模型1 1的的MATLABMATLAB程序：程序： a=0.006計算結(jié)果：計算結(jié)果： 4. 在a=0.006=0.006附近有一個轉(zhuǎn)折點，在這一點左邊，風(fēng)險增加很少時，利潤增長很快。在

26、這一點右邊，風(fēng)險增加很大時，利潤增長很緩慢，所以對于風(fēng)險和收益沒有特殊偏好的投資者來說，應(yīng)該選擇曲線的拐點作為最優(yōu)投資組合， 大約是a*=0.6%，Q*=20% ，所對應(yīng)投資方案為: 風(fēng)險度風(fēng)險度 收益收益 x0 x1 x2 x3 x4 0.0060 0.2019 0 0.2400 0.4000 0.1091 0.2212 3. 3. 曲線上的任一點都表示該風(fēng)險水平的最大可能收益和該收益要求的最小風(fēng)險。對于不同風(fēng)險的承受能力，選擇該風(fēng)險水平下的最優(yōu)投資組合。當(dāng)投資越分散時，投資者承擔(dān)的風(fēng)險越小，這與題意一致。 即: 冒險的投資者會出現(xiàn)集中投資的情況，保守的投資者則盡量分散投資。1. 1. 風(fēng)險

27、大，收益也大。模型模型1 1的結(jié)果分析的結(jié)果分析 此模型又可改寫為min y001s. t . 0,00,1,niiiiniiiiiirpxkpxMxqyxyinL模型模型2 2的求解：的求解： min y. . . 01234012341230 050 270 190 1850 185 101 10210451065=10 0250 0150 055xxxxxkxxxxxxyxyxy . (),40 0260 0，1， ，40ixyxiy由于k是任意給定的盈利，到底怎樣給定沒有一個準(zhǔn)則，不同的投資者有不同的盈利.我們從k=0.05開始，以步長k=0.01進(jìn)行循環(huán)搜索，編制程序如下：模型模型2

28、 2的求解：的求解： k=0.05while k0.26/1.01;C= 0 0 0 0 0 1;A= 0 0.025 0 0 0 -1;0 0 0.015 0 0 -1;0 0 0 0.055 0 -1;0 0 0 0 0.026 -1;B= 0;0;0;0;Aeq= 0.05 0.27 0.19 0.185 0.185,0;1 1.01 1.02 1.045 1.065,0;Beq= k;1;Vlb= 0;0;0;0;0;0;% or Vlb= zeros(6,1);Vub= ;x,fval= linprog(C,A,B,Aeq,Beq,Vlb,Vub);模型模型2 2的的MATLABMA

29、TLAB求解：求解： kQ=fval x=xplot(k, Q, m.)axis(0 0.5 0 0.05)xlabel(收益k) ylabel(最小風(fēng)險度Q)title(最小風(fēng)險度Q隨收益R的變化趨勢圖)hold onk=k+0.01;grid onend模型模型2 2的的MATLABMATLAB求解：求解： 00.050.10.150.20.250.30.350.40.450.500.0050.010.0150.020.0250.030.0350.040.0450.05收 益 k最小風(fēng)險度Q最 小 風(fēng) 險 度 Q隨 收 益 R的 變 化 趨 勢 圖模型模型2 2的結(jié)果分析：的結(jié)果分析： 此

30、模型又可改寫為此模型又可改寫為0min (sy-(1-s)niiiirpx01s. t . 0,0niiiiiipxMxqyxy模型模型3 3的求解：的求解： 012340.050min (.270.190.1850.185sy-(1-s)()xxxxx. . . . (),012341234 101 10210451065=10 0250 0150 0550 0260 0，1， ，40ixxxxxxyxyxyxyxiy模型模型3 3的求解：的求解： s=0while s1;C= -0.05*(1-s),-0.27*(1-s),-0.19*(1-s),-0.185*(1-s),-0.185*(1-s),s;A= 0 0.025 0 0 0 -1;0 0 0.015 0 0 -1;0 0 0 0.055 0 -1;0 0 0 0 0.026 -1;B= 0;0;0;0;Aeq= 1 1.01 1.02 1.045 1.065,0;Beq= 1;Vlb= 0;0;0;0;0;0;% or Vlb= zeros(6,1);Vub= ;x,fval= linprog(C,A,B,Aeq,Beq,Vlb,Vub);模型模型3 3的的MATLABMATLAB程序：程序： sQ=x(6) x=xplot(s,Q,r.)axis(0 1 0 0.025)xlabel(權(quán)重s)ylabel(