最優化方法

課程設計報告

班級：________________

姓名:______________

學號：______________

成績：______________

2023年5月21日

目錄

一、摘要錯誤!未定義書簽。

二、單純形算法...........................錯誤!未定義書簽。

1.1單純形算法的基本思緒..................................錯誤!未定義書簽。

1.2算法流程圖.............................................錯誤!未定義書簽。

1.3用matlab編寫源程序...................................錯誤!未定義書簽。

二、黃金分割法錯誤!未定義書簽。

2.1黃金分割法的基本思緒..................................錯誤!未定義書簽。

2.2算法流程圖..............................................錯誤!未定義書簽。

2.3用matlab編寫源程序....................................錯誤!未定義書簽。

2.4黃金分割法應用舉例.....................................錯誤!未定義書簽。

三、最速下降法..........................錯誤!未定義書簽。

3.1最速下降法的基本思緒..................................錯誤!未定義書簽。

3.2算法流程圖.............................................錯誤!未定義書簽。

3.3用matlab編寫源程序....................................錯誤!未定義書簽。

3.4最速下降法應用舉例.....................................錯誤!未定義書簽。

四、處罰函數法……錯誤!未定義書簽。

4.1處罰函數法的基本思緒.....錯誤!未定義書簽。

4.2算法流程圖..............錯誤!未定義書簽。

4.3用matlab編寫源程序….....錯誤!未定義書簽。

4.4處罰函數法應用舉例…......錯誤!未定義書簽。

五、自我總結......錯誤!未定義書簽。

六、參考文獻錯誤!未定義書簽。

一、摘要

運籌學是一門以人機系統的組織、管理為對象，應用數學和計算機等工具來

研究各類有限資源的合理規劃使用并提供優化決策方案的科學。通過對數據的調

查、收集和記錄分析，以及具體模型的建立。收集和記錄上述擬定之模型所需要

的各種基礎數據，并最終將數據整理形成分析和解決問題的具體模型。

最優化理論和方法日益受到重視，己經滲透到生產、管理、商業、軍事、決

策等各個領域，而最優化模型與方法廣泛應用于工業、農業、交通運送、商業、

國防、建筑、通信、政府機關等各個部門及各個領域。隨著著計算機技術的高速

發展，最優化理論與方法的迅速進步為解決實際最優化問題的軟件也在飛速發

展。其中，MATLAB軟件已經成為最優化領域應用最廣的軟件之一。有了MATLAB

這個強大的計算平臺，既可以運用MATLAB優化工具箱(OptimizationToolbox)

中的函數，又可以通過算法變成實現相應的最優化計算。

關鍵詞：優化、線性規劃、黃金分割法、最速下降法、處罰函數法

二、單純形算法

1.1單純形算法的基本思緒

線性規劃問題的可行域是n維向量空間Rn中的多面凸集，其最優值假如存

在必在該凸集的某頂點處達成。頂點所相應的可行解稱為基本可行解。

單純形法的基本思想是：先找出一個基本可行解，對它進行鑒別，看是否是

最優解；若不是，則按照一定法則轉換到另一改善的基本可行解，再鑒別；若仍

不是，則再轉換，按此反復進行。因基本可行解的個數有限，故經有限次轉換必

能得出問題的最優解。假如問題無最優解也可用此法判別。

單純形法是從某一基可行解出發，連續地尋找相鄰的基可行解，直到達成最

優的迭代過程，其實質是解線性方程組。

概述：根據單純形法的原理，在線性規劃問題中，決策變量（控制變量）xl,

x2,…xn的值稱為一個解，滿足所有的約束條件的解稱為可行解。使目的函數達

成最大值（或最小值）的可行解稱為最優解。這樣，一個最優解能在整個由約束

條件所擬定的可行區域內使目的函數達成最大值（或最小值）。

求解線性規劃問題的目的就是要找出最優解。用單純形法求解線性規劃問題

所需的迭代次數重要取決于約束條件的個數。現在一般的線性規劃問題都是應用

單純形法標準軟件在計算機上求解，對于具有106個決策變量和104個約束條件

的線性規劃問題已能在計算機上解得。

求解時也許出現下列情況之一：①存在著一個最優解；②存在著無窮多個最

優解；③不存在最優解，這只在兩種情況下發生，即沒有可行解或各項約束條件

不阻止目的函數的值無限增大（或向負的方向無限增大）。

要縮小對最優解的搜索范圍，就必須結識最優解的一般性質，最優解假如存

在的話，則它必然處在可行區域的邊界上。

任何一項約束條件的邊界方程是用“=”號來替換該約束條件中的“W”

或“2”號而得到的。每一個邊界方程擬定一個超平面。因此，可行區域的邊界

是由那些滿足一個或同時滿足幾個邊界方程(即處在作為邊界的一個或幾個超平

面上)的可行解所組成，并且最優解必在其中。最優解不僅是在可行區域的邊界

上，并且也在這個區域的一個隅角上。一個可行解，假如不處在由另兩個可行解

連接起來的任何線段上，它就是一個角點可行解。假如連接兩個角點可行解的線

段處在可行區域的邊界上，這兩個角點可行解就稱為相鄰的角點可行解。

角點可行解具有下列三個重要性質：

①假如存在著一個最優解，那么它必然是角點可行解。假如存在有多個最優

解，那么至少有兩個最優解必然是相鄰的角點可行解。

②只存在有限個數的角點可行解。

③假如一個角點可行解按目的函數值來衡量時比其所有的相鄰角點可行解

更好一些，那它就比所有其他角點可行解都更好，也就是最優解。

上述這些性質構成單純形法的原理基礎。

最后一個性質的重要性在于它為一個角點可行解是否是最優解提供了一種簡便

的檢查標準，因而毋需列舉所有的可行解。單純形法正是運用了這個性質，只要

檢查少數的角點可行解，并且一旦這個最優性檢查獲得通過就可立即停止運算。

1.2算法流程圖

(1)、擬定初始基可行解

①從線性規劃標準形的系數矩陣中能直接找出m個線性獨立的單位向量；

②對約束條件全為“<=”連接的LP,化為標準形，左端添加松弛變量后即形

成一個單位子矩陣；

③約束條件中具有“仁”或“=”連接的方程，在插入剩余變量后找不到單

位矩陣，則必須采用“人造基”法，

(2)、單純形法的運算環節可歸結為:

①起始環節——在一個角點可行解上開始。

②迭代環節一一移動至一個更好一些的相鄰角點可行解（根據需要反復進行

這一環節）。

③停止法則——在當前角點可行解比所有相鄰角點可行解都更好些時停止。

當前角點可行解就是一個最優解。單純形法的優點及其成功之處在于它只需要

較少的有限次數的迭代，即可找到最優解。

（3）、單純形法的算法流程如下：

①把線性規劃問題的約束方程組表達成典范型方程組，找出基本可行解作為

初始基本可行解。

②若基本可行解不存在，即約束條件有矛盾，則問題無解。

③若基本可行解存在，從初始基本可行解作為起點，根據最優性條件和可行

性條件，引入非基變量取代某一基變量,找出目的函數值更優的另一基本可行解。

④按環節3進行迭代,直到相應檢查數滿足最優性條件（這時目的函數值不

能再改善），即得到問題的最優解。⑤若迭代過程中發現問題的目的函數值無界,

則終止迭代。

1.3用matlab編寫源程序

Matlab程序源代碼：

-----------simplexTab.m子函數-------------------------------------

functionsimplexTab(mat,numFreeVar)

maxRow=length(mat(:,1));

maxCol=length(mat(1

objEntryExcludingMaxPayOff=mat(maxRow,1:maxCol-2);

[objEntbestColToPivot]=min(objEntryExcludingMaxPayOff)；

while(objEnt<0)

lastColExcludingObjEnty=mat(1:(maxRow-1),maxCol);

ithColExcludingObjEnty=niat(l:(maxRow-1),bestColToPivot);

a=lastColExcludingObjEnty./ithColExcludingObjEnty;

[valbestRowToPivot]=min(a);

sprintf(tthebestPivotis%drowand%dcol',bestRowToPivot,bestColToPivot)

dispC單純形表化為：，)；

[mat,[a;0]]

dispC按任意鍵繼續上

pause;

if(val<0)

[sindices]=sort(a);

if(max(a)>0)

count=l;

while(s(count)<0)

count=count+l;

end

bestRowToPivot=indices(count);

end

end

if(length(a)==O)

length(a)

return

end

mat=pivot(mat,bestRowToPivot,bestColToPivot);

objEntryExcludingMaxPayOff=mat(maxRow,1:maxCol-2);

[objEntbestColToPivot]=min(objEntryExcludingMaxPayOff);

end

sprintf(*thebestPicotis%d行and%d列',bestRowToPivot,bestColToPivot)

dispC單純形表化為：，)；

[mat,[a;0]]

dispC運營結束!);

----------interChange.m子函數-------------------------------------

functionnewMat=interChange(mat,row1,row2)

temp=mat(rowl,:);

mat(rowl,:)=mat(row2,:);

mat(row2,:)=temp;

newMat=mat;

----------multiFromRowIbRow.m子函數-----------------------------

functionnewMat=multiFromRowToRow(mat,fromRow,toRow,multiplier)

rG=mat(fromRow,:)*multiplier;

mat(toRow,:)=rG+mat(toRow,:);

newMat=mat;

-----------multMat.m子函數----------------------------------------

fuctionnewMat=multMat(mat,row,mult)

%multiplyarowofamatrixbyanonzeroconstant

mat(row,:)=mat(row,:)*mult;

newMat=mat;

-----------pivot.m子函數--------------------------------------------

functionnewMat=pivot(mat,row,col)

%normalizethisrow

mat(row,:)=mat(row,:)./mat(row,col);%maketheleadinganumbera1

fori-l:length(mat(:,l))

if(r-=row)

mat=multiFromRowToRow(mat,row,r,-mat(r,col));

end

end

newMat=mat;

1.4單純形算法應用舉例

題目：使用單純形法解下面線性規劃問題:

目的函數為：maxf=一七一巧一七

約束條件是：（7毛+3巧+9^<1

s.t.8才］+5X2+44<1

63+9X2+5X3<1

X],巧，X3>0

解：化為標準形式：

maxg=3+巧+巧+。*4+。*5+。*6

(7/+3*2+9*3+才4+0*5+0才6=1

S.t.8x{+5X2+4*3+0*4+*5+0才6=1

《

6/+9巧+5X3+0才4+0才5+才6=1

I才1，*2,*3,才4,才5,才6N0

運用Matlab程序計算：

命令窗口輸入：

mat=[7391001;

8540101;

6950011;

-1-1-10000];

numFreeVar=3;%自由變量個數

simplexTab(mat,numFreeVar)

結果輸出：

初始結果：

ans=thebestPivotis2rowand1col

單純形表化為：

7.00003.00009.00001.0000001.00000.1429

8.00005.00004.000001.000001.00000.1250

6.00009.00005.0000001.00001.00000.1667

-1.0000-1.0000-1.000000000

第一次轉換結果:

ans=thebestPivotis1rowand3col

單純形表化為:

0-1.37505.50001.0000-0.875000.12500.0227

1.00000.62500.500000.125000.12500.2500

05.25002.00000-0.75001.00000.25000.1250

0-0.3750-0.500000.125000.12500

第二次轉換結果：

ans=thebestPivotis1rowand2col

單純形表化為：

0-0.25001.00000.1818-0.159100.0227-0.0909

1.00000.75000-0.09090.204500.11360.1515

05.75000-0.3636-0.43181.00000.20450.0356

0-0.500000.09090.045500.13640

第三次轉換結果：

ans=thebestPicotis3rowand1col

單純形表化為：

001.00000.1660-0.17790.04350.0316-0.0909

1.000000-0.04350.2609-0.13040.08700.1515

01.00000-0.0632-0.07510.17390.03560.0356

0000.05930.00790.08700.15420

運營結束！（并表達程序運營完畢）

最后單純形表為:

X2%局X6btheta

001.00000.1660-0.17790.04350.0316-0.0909

A1.000000-0.04350.2609-0.13040.08700.1515

01.00000-0.0632-0.07510.17390.03560.0356

判別數0000.05930.00790.08700.15420

x二(0.0316,0.0870,0.0356),=0.1542,

f=—-x2-x2=-0.1542

所以該線性規劃的最優解是g=o.1542。

二、黃金分割法

2.1黃金分割法的基本思緒

一維搜索是解函數極小值的方法之一，其解法思想為沿某一已知方向求目的

函數的極小值點。一維搜索的解法很多，這里重要采用黃金分割法(0.618法)。

該方法用不變的區間縮短率0.618代替斐波那契法每次不同的縮短率，從而可以

當作是斐波那契法的近似，實現起來比較容易，也易于人們所接受。

黃金分割法是用于一元函數f(x)在給定初始區間［a,b］內搜索極小點xmin的

一種方法。它是優化計算中的經典算法，以算法簡樸、收斂速度均勻、效果較好

而著稱，是許多優化算法的基礎，但它只合用于一維區間上的凸函數，即只在單

峰區間內才干進行一維尋優，其收斂效率較低。

其基本原理是：

依照“去劣存優”原則、對稱原則、以及等比收縮原則來逐步縮小搜索區間。

具體環節是：

在區間［a,b］內取點：a},a2把［a,b］分為三段。

①假如fQ［)>〃%),令2=@［,al=a2,a2=a+0.618*(b-a)；

②假如/1(4)<￡(%)，令b=%，a—al=b~0.618*(b-a)；

假如-a)/0和恒-%)/九|都大于收斂精度￡重新開始循環。

由于［a,b］為單峰區間，這樣每次可將搜索區間縮小0.618倍，解決后的區

間都將包含極小點的區間縮小，然后在保存下來的區間上作同樣的解決，如此迭

代下去，將使搜索區［a,b］逐步縮小，直到滿足預先給定的精度時，即獲得一維

優化問題的近似最優解。

插入點原理圖如下：

1—A=A2,+A—1=0A=0.618

所謂的黃金分割法，是指將一線段提成兩段的方法，使整段的長與較長段的

長的長度比值等于較長段與較短段的比值，即：

1:2=/I:(1-2)

2.2算法流程圖

黃金分割法的算法環節：

給定a,b(a>b),控制誤差￡>0。

一般取［a,b］為［0,1］,當然也可以根據自己需要調整。假如目的函數對變

量非常靈敏，可以將補償的可行范圍進行范圍縮小，以提高搜索精度；反之，可

以將步長的可行范圍放大，以提高算法的計算速度。

2.3用matlab編寫源程序

Matlab程序源代碼：

-------------golden,m子函數---------------

functionPa=golden(x,p,LowBound,UpBound)

%goldenisfunctionforgoldenlinearsearch

%x(k+1)=x(k)+Pa*p;

%input:

%x:搜索初始點

%p：搜索方向

%LowBoim&UpBound:搜索區間

%output:

%Pa:搜索步長

%disp('itisthegoldenlinereasch1);

%搜索步長范圍［0,1］可以根據需要修改

gold_a=LowBound;

gold_b=UpBound;

z2=gold_a+0.618*(gold_b-gold_a);%%stepl

f2=objfun(x+z2*p);

z1=gold_a+0.382*(gold_b-gold_a);%%step2

fl=objfun(x+z1*p);

in_itera=0;

whileabs(gold_b-gold_a)>0.

%absa=0,b=lgoldenlinesearch

in_itera=in_itera+1;

%fl,f2

iffl<f2

gold_b=z2;

z2=zl;

f2=fl;

z1=gold_a+0.382*(gold_b-gold_a);

fl=objfun(x+zl*p);

continue

elseiffl==f2

gold_a=zl;

gold_b=z2;

z2=gold_a+0.618*(gold_b-gold_a);%%step2

f2=objfun(x+z2*p);

z1=gold_a+0.382*(gold_b-gold_a);%%step2

fl=objfun(x+z1*p);

continue

else

gold_a=zl;

zl=z2;

fl=f2;

z2=gold_a+0.618*(gold_b-gold_a);

f2=objfun(x+z2*p);

end

end

Pa=(gold_a+gold_b)/2;%Pax(k+1)=x(k)+Pa*p;

dispC最優解是：')

%Pa,in_itera

%obj(x+Pa*p)

-----------------objfun.m子函數---------------------------

functionf=objfun(x)

f=x(l)A3+x(2)A2-l0*x(l)*x(2)+1;

2.4黃金分割法應用舉例

題目：運用黃金分割法求解下面函數的最優解：

f=X；+君一10占42+1

解：

運用Matlab程序計算:

命令窗口輸入：

X=[O,O];P=[1,1];

Pa=golden(X,P,0,l)

結果輸出：

最優解是：Pa=1.0000

三、最速下降法

3.1最速下降法的基本思緒

最速下降法的搜索法向是目的函數的負梯度方向，最速下降法從目的函數的

負梯度方向一直前進，直到到達目的函數的最低點。

己知目的函數在X“）點的梯度為：

當求目的函數的最小點時，由于函數沿負梯度方向下降最快，故在點的

探索方向應取該點的負梯度方向，即

顯然，為單位向量。這樣第Z+1次迭代計算所得的新點為

a伏）

X（Z+D=x（*）伏）S（*）=X（k）_______,,k）

■■—

負梯度僅給出了最優化方向，而沒有給出步長的大小，所以也許有各種各樣的最

速下降的過程，它們依賴于的大小。

步長a")有兩種取法：

一種方法是任意給定一個初始步長，使滿足條件：

f(X(k)+a(k}S(k})<f(X(k))

此外一種方法是沿負梯度方向做一維探索，以求解一維最優化問題的最優步

長a,即對目的函數極小，以得到最優步長：

min/(X*)+aS(k))=f(X<k)+a{k}S(k))

a>0

以此最優步長作為由點出發沿該點的負梯度方向探索的步長a")。

這種方法的迭代計算的收斂性，可用以下三式中的任一式或二式作為準則來

進行判斷：

L叫《與

|小巧-/W

用最速下降法求無約束多維極值問題min的算法環節如下：

⑴、取初始點x⑼，精度￡>°,令k=0

(2)、計算搜索方向*)=-9(心>),其中W電力表達函數/Xx)在點燎處

的梯度；

(3)、若卜⑻|?￡，則停止計算；否則，從V出發，沿進行一維搜索，

即求《，使得/'(”+4儼>)=minf(”+/lW))。此處的一維搜索可以用

MATLAB的fminUAC函數；令/旬=x?+4儼>,左=Z+1,轉環節(2)。

3.2算法流程圖

3.3用matlab編寫源程序

-----------BanaFunWithGrad.m子函數--------------------------------

function[f,g]=BanaFunWithGrad(x)%(不含導數解析式)

f=100*(x(2)-x(1)A2)A2+(1-x(1))A2;

g=[100*(4*x(1)A3-4*x(1)*x(2))+2*x(1)-2;100*(2*x(2)-2*x(1)A2)];

-----------BanaFun.m子函數----------------------------------------

functionf=BanaFun(x)%(含導數解析式)

f=100*(x(2)-x(l)A2)A2+(l-x(l))A2;

3.4最速下降法應用舉例

題目：求解函數f=100(X2-+(1-X廣

解：

運用Matlab程序計算：

命令窗口輸入：

OPTIONS^ptimsetC^argeScaleVoff/HessUpdate'/steepdesc^^radobjVonVMaxFu

nEvals,,250,display*,4tef);

x=[-1.9,2];

[x,fval,exitflag,output]=fminunc(@BanaFunWithGrad,x,OPTIONS)

結果輸出：

迭代次數函數計算次數函數值步長一階導數最優性

Gradient's

IterationFunc-countf(x)Step-sizeinfinity-norm

01267.621.23e+003

12214.4160.519

295.836390.0009849313.4

3155.782920.1.58

4185.731270.038797912.9

5245.680230.1.59

6275.630810.036759912.5

7335.58190.1.6

8365.534440.034948212.1

9425.487430.1.61

10455.441730.033324511.7

11515.396410.1.62

12545.352280.031858711.3

13605.308490.1.63

14635.265790.030527111

15695.223380.1.64

16725.181970.029310610.7

17785.140820.1.64

18815.100590.028193710.4

19875.060580.1.65

20905.021440.027163510.1

21964.982490.1.66

22994.944330.02620959.89

231054.906350.1.67

241084.869120.02532289.65

251144.832040.1.68

261174.795660.02449589.42

271234.75940.1.69

281264.723810.0237229.2

291324.688330.1.7

301354.653470.0229968.99

311414.618710.1.7

321444.584530.0223138.79

331504.550430.1.71

341534.51690.02166898.6

351594.483430.1.72

361624.450490.02106018.42

371684.41760.1.73

381714.385220.02048348.24

391774.352890.1.74

401804.321030.0199368.07

411864.28920.1.75

421894.257840.01941547.91

431954.22650.1.76

441984.195590.01891947.75

452044.16470.1.77

462074.134230.01844627.6

472134.103770.001008181.78

482164.073710.01799397.45

492224.043640.001035051.79

502254.013960.01756097.3

512313.984270.001062841.8

522343.954950.01714587.16

532403.925610.001088861.8

542433.894410.01810857.29

552493.863210.001117641.82

Maximumnumberoffunctionevaluationsexceeded;（達成最大的函數值計算次數）

increaseoptions.MaxFunEvals（建議增長最大的函數值計算次數）

x=-0.96340.9372（最后迭代點）

fval=3.8632（最優點相應的函數值）

exitflag=0

output=（函數基本信息）

iterations:56（迭代次數）

funcCount:250（目的函數最大計算次數）

stepsize:0.0011(步長)

firstorderopt:1.8155

algorithm:'medium-scale:Quasi-Newtonlinesearch1

message:11x78charj

從結果可見:在250次最大的目的函數計算次數內，算法沒有到最優點[1,1]o

將最大目的函數計算次數調高到2500次，則：

命令窗口輸入：

OPTIONS=optimset('LargeScale?off,'HessUpdate','steepdesc?gradobj?on?MaxFu

nEvals1,2500,display*,'iter*);

x=[-1.9,21;

[x,fval,exitflag,output]=fminunc(@BanaFunWithGrad,x,OPTIONS)

結果輸出：

迭代次數函數計算次數函數值步長一階導數最優性

Gradients

IterationFunc-countf(x)Step-sizeinfinity-norm

01267.621.23e+003

12214.4160.519

295.836390.0009849313.4

3155.782920.1.58

4185.731270.038797912.9

5245.680230.1.59

6275.630810.036759912.5

7335.58190.1.6

8365.534440.034948212.1

9425.487430.1.61

10455.441730.033324511.7

11515.396410.1.62

12545.352280.031858711.3

13605.308490.1.63

14635.265790.030527111

15695.223380.1.64

16725.181970.029310610.7

17785.140820.1.64

18815.100590.028193710.4

19875.060580.1.65

40014590.0.001117680.0162

40114620.0.01415250.0574

Maximumnumberofiterationsexceeded;

increaseoptions.MaxIter.

x=0.97700.9542

fval=5.3808e-004

exitflag=0

output=

iterations:401

funcCount:1462

stepsize:0.0142

firstorderopt:0.0574

algorithm:'medium-scale:Quasi-Newtonlinesearch*

message:[1x67char]

從以上兩個算法的結果可以看出，最速下降法對bana函數不是十分有效。

當目的函數不可微或者導數求解復雜時，可以無需提供目的函數的導函

數的解析形式，MATLAB可以使用差分的方法求導(下降方向)。

命令窗口輸入相應的代碼為：

OPTIONS=optimset('LargeScale,,,off,,HessUpdate7steepdesc',,MaxFunEvals\250);

x=[-1.9,2];

fx,fval,exitflag,output]=fminunc(@BanaFun,x,OPTIONS)

結果輸出：

Maximumnumberoffunctionevaluationsexceeded;

increaseoptions.MaxFunEvals

x=-1.25531.5640

fval=5.1002

exitflag=0

output=

iterations:19

funcCount:252

stepsize:0.0282

firstorderopt:10.4149

algorithm:'medium-scale:Quasi-Newtonlinesearch*

message:[1x78char]

從計算結果可以明顯看出，這對于最速下降發，有導數解析式會比沒有導數

解析式時，計算質量會有所提高。

四、處罰函數法

4.1處罰函數法的基本思緒

處罰函數法是應用廣泛,非常有效的間接解法，又稱為序列無約束極小化方

法(SUMT法)。該方法通過將原約束優化問題中的等式和不等式約束函數加權解

決后與原目的函數結合，得到新的目的函數(處罰函數)。原問題轉化為新的無約

束優化問題,求解該新的無約束優化問題，間接得