數(shù)學(xué)建模課件之非線性規(guī)劃問(wèn)題的算法

數(shù)學(xué)建模課件之非線性規(guī)劃問(wèn)題的算法第一頁(yè)，共五十八頁(yè)，2022年，8月28日一、非線性規(guī)劃問(wèn)題的幾種求解方法

1.罰函數(shù)法（外點(diǎn)法）

基本思想：利用目標(biāo)函數(shù)和約束函數(shù)構(gòu)造輔助函數(shù)：第二頁(yè)，共五十八頁(yè)，2022年，8月28日要求構(gòu)造的函數(shù)具有這樣的性質(zhì)：當(dāng)點(diǎn)x位于可行域以外時(shí)，取值很大，而離可行域越遠(yuǎn)則越大；當(dāng)點(diǎn)在可行域內(nèi)時(shí)，函數(shù)因此可以將前面的有約束規(guī)劃問(wèn)題轉(zhuǎn)換為下列無(wú)約束規(guī)劃模型：其中稱(chēng)為罰項(xiàng)，稱(chēng)為罰因子，稱(chēng)為罰函數(shù)。第三頁(yè)，共五十八頁(yè)，2022年，8月28日的定義一般如下：函數(shù)一般定義如下：第四頁(yè)，共五十八頁(yè)，2022年，8月28日算法步驟如何將此算法模塊化:第五頁(yè)，共五十八頁(yè)，2022年，8月28日求解非線性規(guī)劃模型例子罰項(xiàng)函數(shù)：無(wú)約束規(guī)劃目標(biāo)函數(shù)：第六頁(yè)，共五十八頁(yè)，2022年，8月28日globallamada%主程序main2.m,罰函數(shù)方法x0=[11];lamada=2;c=10;e=1e-5;k=1;whilelamada*fun2p(x0)>=e x0=fminsearch('fun2min',x0);lamada=c*lamada;k=k+1;enddisp(‘最優(yōu)解’),disp(x0)disp('k='),disp(k)

程序1：主程序main2.m第七頁(yè)，共五十八頁(yè)，2022年，8月28日程序2：計(jì)算的函數(shù)fun2p.mfunctionr=fun2p(x)%罰項(xiàng)函數(shù)r=((x(1)-1)^3-x(2)*x(2))^2;

第八頁(yè)，共五十八頁(yè)，2022年，8月28日程序3：輔助函數(shù)程序fun2min.mfunctionr=fun2min(x)%輔助函數(shù)globallamadar=x(1)^2+x(2)^2+lamada*fun2p(x);第九頁(yè)，共五十八頁(yè)，2022年，8月28日運(yùn)行輸出：最優(yōu)解

k=33第十頁(yè)，共五十八頁(yè)，2022年，8月28日練習(xí)題：

1、用外點(diǎn)法求解下列模型2、將例子程序改寫(xiě)為一個(gè)較為通用的罰函數(shù)法程序。（考慮要提供哪些參數(shù)）第十一頁(yè)，共五十八頁(yè)，2022年，8月28日2.內(nèi)點(diǎn)法（障礙函數(shù)法）僅適合于不等式約束的最優(yōu)化問(wèn)題其中都是連續(xù)函數(shù)，將模型的定義域記為第十二頁(yè)，共五十八頁(yè)，2022年，8月28日構(gòu)造輔助函數(shù)為了保持迭代點(diǎn)含于可行域內(nèi)部，我們定義障礙函數(shù)第十三頁(yè)，共五十八頁(yè)，2022年，8月28日3.問(wèn)題轉(zhuǎn)化為一個(gè)無(wú)約束規(guī)劃由于很小，則函數(shù)取值接近于f(x)，所以原問(wèn)題可以歸結(jié)為如下規(guī)劃問(wèn)題的近似解：第十四頁(yè)，共五十八頁(yè)，2022年，8月28日第十五頁(yè)，共五十八頁(yè)，2022年，8月28日練習(xí)題：

請(qǐng)用內(nèi)點(diǎn)法算法求解下列問(wèn)題：

第十六頁(yè)，共五十八頁(yè)，2022年，8月28日小結(jié)講解了兩個(gè)求解有約束非線性最小化規(guī)劃特點(diǎn)：易于實(shí)現(xiàn)，方法簡(jiǎn)單；沒(méi)有用到目標(biāo)函數(shù)的導(dǎo)數(shù)問(wèn)題的轉(zhuǎn)化技巧（近似為一個(gè)無(wú)約束規(guī)劃）第十七頁(yè)，共五十八頁(yè)，2022年，8月28日4、其它求解算法

（1）間接法

（2）直接法

直接搜索法以梯度法為基礎(chǔ)的間接法無(wú)約束規(guī)劃的Matlab求解函數(shù)數(shù)學(xué)建模案例分析（截?cái)嗲懈?，飛機(jī)排隊(duì)）第十八頁(yè)，共五十八頁(yè)，2022年，8月28日（1）間接法 在非線性最優(yōu)化問(wèn)題當(dāng)中，如果目標(biāo)函數(shù)能以解析函數(shù)表示，可行域由不等式約束確定，則可以利用目標(biāo)函數(shù)和可行域的已知性質(zhì)，在理論上推導(dǎo)出目標(biāo)函數(shù)為最優(yōu)值的必要條件，這種方法就稱(chēng)為間接法（也稱(chēng)為解析法）。一般要用到目標(biāo)函數(shù)的導(dǎo)數(shù)。第十九頁(yè)，共五十八頁(yè)，2022年，8月28日（2）直接法

直接法是一種數(shù)值方法這種方法的基本思想是迭代，通過(guò)迭代產(chǎn)生一個(gè)點(diǎn)序列{X(k)}，使之逐步接近最優(yōu)點(diǎn)。只用到目標(biāo)函數(shù)。如黃金分割法、Fibonacci、隨機(jī)搜索法。第二十頁(yè)，共五十八頁(yè)，2022年，8月28日（3）迭代法一般步驟注意：數(shù)值求解最優(yōu)化問(wèn)題的計(jì)算效率取決于確定搜索方向P

(k)和步長(zhǎng)的效率。第二十一頁(yè)，共五十八頁(yè)，2022年，8月28日最速下降法（steepestdescentmethod）由法國(guó)數(shù)學(xué)家Cauchy于1847年首先提出。在每次迭代中，沿最速下降方向（負(fù)梯度方向）進(jìn)行搜索，每步沿負(fù)梯度方向取最優(yōu)步長(zhǎng)，因此這種方法稱(chēng)為最優(yōu)梯度法。特點(diǎn)：方法簡(jiǎn)單，只以一階梯度的信息確定下一步的搜索方向，收斂速度慢；越是接近極值點(diǎn)，收斂越慢；它是其它許多無(wú)約束、有約束最優(yōu)化方法的基礎(chǔ)。該法一般用于最優(yōu)化開(kāi)始的幾步搜索。第二十二頁(yè)，共五十八頁(yè)，2022年，8月28日以梯度法為基礎(chǔ)的最優(yōu)化方法求f(x)在En中的極小點(diǎn)思想：方向?qū)?shù)是反映函數(shù)值沿某一方向的變化率問(wèn)題方向?qū)?shù)沿梯度方向取得最大值基礎(chǔ)：方向?qū)?shù)、梯度第二十三頁(yè)，共五十八頁(yè)，2022年，8月28日通過(guò)一系列一維搜索來(lái)實(shí)現(xiàn)。本方法的核心問(wèn)題是選擇搜索方向。搜索方向的不同則形成不同的最優(yōu)化方法。第二十四頁(yè)，共五十八頁(yè)，2022年，8月28日最速下降法算法：

第二十五頁(yè)，共五十八頁(yè)，2022年，8月28日算法說(shuō)明可通過(guò)一維無(wú)約束搜索方法求解第二十六頁(yè)，共五十八頁(yè)，2022年，8月28日例子：用最速下降法解下列問(wèn)題分析：1、編寫(xiě)一個(gè)梯度函數(shù)程序fun1gra.m2、求(可以調(diào)用函數(shù)fminsearch)函數(shù)fungetlamada.m3、最速下降法主程序main1.m初始條件第二十七頁(yè)，共五十八頁(yè)，2022年，8月28日第一步：計(jì)算梯度程序fun1gra.mfunctionr=fun1gra(x)%最速下降法求解示例%函數(shù)f(x)=2*x1^2+x2^2的梯度的計(jì)算%r(1)=4*x(1);r(2)=2*x(2);第二十八頁(yè)，共五十八頁(yè)，2022年，8月28日第二步：求最優(yōu)的目標(biāo)函數(shù)functionr=fungetlamada(lamada)%關(guān)于lamada的一元函數(shù)，求最優(yōu)步長(zhǎng)globalx0d=fun1gra(x0);r=2*(x0(1)-lamada*d(1))^2+(x0(2)-lamada*d(2))^2;%注意負(fù)號(hào)表示是負(fù)梯度第二十九頁(yè)，共五十八頁(yè)，2022年，8月28日第三步：主程序main1.m%最速下降方法實(shí)現(xiàn)一個(gè)非線性最優(yōu)化問(wèn)題%minf(x)=2*x1^2+x2^2globalx0x0=[11];yefi=0.0001;k=1;d=-fun1gra(x0);lamada=1;第三十頁(yè)，共五十八頁(yè)，2022年，8月28日主程序main1.m（續(xù)）whilesqrt(sum(d.^2))>=yefilamada=fminsearch(‘fungetlamada’,lamada);%求最優(yōu)步長(zhǎng)lamadax0=x0-lamada*fun1gra(x0);%計(jì)算x0d=fun1gra(x0);%計(jì)算梯度

k=k+1;%迭代次數(shù)enddisp('x='),disp(x0),disp('k='),disp(k),disp('funobj='),disp(2*x0(1)^2+x0(2)^2)第三十一頁(yè)，共五十八頁(yè)，2022年，8月28日三、Matlab求解有約束非線性規(guī)劃第三十二頁(yè)，共五十八頁(yè)，2022年，8月28日1.用fmincon函數(shù)求解形如下面的有約束非線性規(guī)劃模型一般形式：第三十三頁(yè)，共五十八頁(yè)，2022年，8月28日用Matlab求解有約束非線性最小化問(wèn)題求解非線性規(guī)劃問(wèn)題的Matlab函數(shù)為：fmincon

1.約束中可以有等式約束 2.可以含線性、非線性約束均可第三十四頁(yè)，共五十八頁(yè)，2022年，8月28日輸入?yún)?shù)語(yǔ)法：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=fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options,P1,P2,...)第三十五頁(yè)，共五十八頁(yè)，2022年，8月28日輸入?yún)?shù)的幾點(diǎn)說(shuō)明模型中如果沒(méi)有A,b,Aeq,beq,lb,ub的限制，則以空矩陣[]作為參數(shù)傳入；nonlcon：如果包含非線性等式或不等式約束，則將這些函數(shù)

編寫(xiě)為一個(gè)Matlab函數(shù)，nonlcon就是定義這些函數(shù)的程序文件名；不等式約束c(x)<=0等式約束ceq(x)=0.如果nonlcon=‘mycon’;則myfun.m定義如下function[c,ceq]=mycon(x)c=...

%計(jì)算非線性不等式約束在點(diǎn)x處的函數(shù)值ceq=...

%計(jì)算機(jī)非線性等式約束在點(diǎn)x處的函數(shù)值

第三十六頁(yè)，共五十八頁(yè)，2022年，8月28日對(duì)參數(shù)nonlcon的進(jìn)一步示例2個(gè)不等式約束，2個(gè)等式約束3個(gè)決策變量x1，x2，x3如果nonlcon以‘mycon1’作為參數(shù)值，則程序mycon1.m如下第三十七頁(yè)，共五十八頁(yè)，2022年，8月28日對(duì)照約束條件編寫(xiě)myfun1.mfunction[c,ceq]=mycon1(x)c(1)=x(1)*x(1)+x(2)*x(2)+x(3)*x(3)-100c(2)=60-x(1)*x(1)+10*x(3)*x(3)ceq(1)=x(1)+x(2)*x(2)+x(3)-80ceq(2)=x(1)^3+x(2)*x(2)+x(3)-80第三十八頁(yè)，共五十八頁(yè)，2022年，8月28日nonlcon的高級(jí)用法允許提供非線性約束條件中函數(shù)的梯度設(shè)置方法：options=optimset('GradConstr','on')

如果提供非線性約束條件中函數(shù)梯度，nonlcon的函數(shù)必須如下格式：第三十九頁(yè)，共五十八頁(yè)，2022年，8月28日參數(shù)nonlcon的函數(shù)一般格式如下function[c,ceq,GC,GCeq]=mycon(x)c=...

%計(jì)算非線性不等式約束在點(diǎn)x處的函數(shù)值

ceq=...

%計(jì)算機(jī)非線性等式約束在點(diǎn)x處的函數(shù)值

ifnargout>2

%nonlcon如果四個(gè)輸出參數(shù)

GC=...

%不等式約束的梯度

GCeq=...

%等式約束的梯度end第四十頁(yè)，共五十八頁(yè)，2022年，8月28日輸出參數(shù)語(yǔ)法：[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(...)運(yùn)用步驟：將自己的模型轉(zhuǎn)化為上面的形式寫(xiě)出對(duì)應(yīng)的參數(shù)調(diào)用函數(shù)第四十一頁(yè)，共五十八頁(yè)，2022年，8月28日f(shuō)mincon應(yīng)用求解示例：請(qǐng)問(wèn)：1、結(jié)合fmincon函數(shù)，需要提供哪些參數(shù)第四十二頁(yè)，共五十八頁(yè)，2022年，8月28日第一步：編寫(xiě)一個(gè)M文件返回目標(biāo)函數(shù)f在點(diǎn)x處的值函數(shù)程序functionf=myfun(x)f=-x(1)*x(2)*x(3);函數(shù)myfun.m第四十三頁(yè)，共五十八頁(yè)，2022年，8月28日第二步：為了調(diào)用MATLAB函數(shù)，必須將模型中的約束轉(zhuǎn)化為如下形式(<=)。

這里：A=[-1-2-2;122];b=[072]’;這是2個(gè)線性約束，形如第四十四頁(yè)，共五十八頁(yè)，2022年，8月28日第三步：提供一個(gè)搜索起點(diǎn)，然后調(diào)用相應(yīng)函數(shù)，程序如下：%給一個(gè)初始搜索點(diǎn) x0=[10;10;10];[x,fval]=fmincon('myfun',x0,A,b)第四十五頁(yè)，共五十八頁(yè)，2022年，8月28日主程序（整體）：A=[-1-2-2;122];b=[072]’;%給一個(gè)初始搜索點(diǎn) x0=[10;10;10];[x,fval]=fmincon('myfun',x0,A,b)第四十六頁(yè)，共五十八頁(yè)，2022年，8月28日最后得到如下結(jié)果：

x=24.000012.000012.0000

fval=-3.4560e+03第四十七頁(yè)，共五十八頁(yè)，2022年，8月28日2.非負(fù)條件下線性最小二乘lsqnonneg

適合如下模型：

注意：約束只有非負(fù)約束第四十八頁(yè)，共五十八頁(yè)，2022年，8月28日語(yǔ)法：x=lsqnonneg(c,d)x=lsqnonneg(c,d,x0)x=lsqnonneg(c,d,x0,options)

第四十九頁(yè)，共五十八頁(yè)，2022年，8月28日3.有約束線性最小二乘lsqlin

適合如下模型：注意：約束有線性等式、不等式約束第五十頁(yè)，共五十八頁(yè)，2022年，8月28日語(yǔ)法：x=lsqlin(C,d,A,b)x=lsqlin(C,d,A,b,Aeq,beq)x=lsqlin(C,d,A,b,Aeq,beq,lb,ub)x=lsqlin(C,d,A,b,Aeq,beq,lb,ub,x0)x=lsqlin(C,d,A,b,Aeq,beq,lb,ub,x0,options)[x,resnorm]=lsqlin(...)[x,resnorm,residual]=lsqlin(...)[x,resnorm,residual,exitflag]=lsqlin(...)[x,resnorm,residual,exitflag,output]=lsqlin(...)[x,resnorm,residual,exitflag,output,lambda]=lsqlin(...)第五十一頁(yè)，共五十八頁(yè)，2022年，8月28日4.非線性最小二乘lsqnonlin

適合模型：第五十二頁(yè)，共五十八頁(yè)，2022年，8月28日語(yǔ)法：x=lsqnonlin(fun,x0)x=lsqnonlin(fun,x0,lb,ub)x=lsqnonlin(fun,x0,lb,ub,options)x=lsqnonlin(fun,x0,options,P1,P2,...)[x,resnorm]=lsqnonlin(...)[x,resnorm,residual]=lsqnonlin(...)[x,resnorm,resid