最近在上最優理論這門課，剛開始是線性規划部分，主要的方法就是單純形方法，學完之后做了一下大M算法和分段法的仿真，拿出來與大家分享一下。單純形方法是求解線性規划問題的一種基本方法。
線性規划就是在一系列不等式約束下求目標函數最大值或最小值的問題，要把數學中的線性規划問題用計算機來解決，首先要確定一個標准形式。
將所給的線性規划問題化為標准形式：
s.t.是英文subject to 的簡寫，意思是受約束，也就是說第一個方程受到后面兩個方程的約束。對於求最大值問題可以將目標函數加負號轉換為最小值問題。
對於求最大值問題可以將目標函數加負號轉換為最小值問題。
其他的問題就是將實際問題中的不等式約束改為等式約束，主要方法是引進松弛變量和剩余變量，以及將自有變量轉換為非負變量。
①對於不等式
，引入松弛變量將其變為等式形式如下：
②對於不等式
，引入剩余變量將其變為等式形式如下：
③若變量為自有變量（可取正、負或零，符號無限制），則引入兩個非負變量將其表示如下：
關於線性規划問題的解：
確定了標准形式，我們就針對這個標准形式討論一下線性規划問題的解。線性規划問題的解能滿足標准形式中約束條件的向量X的值，但只有最優解才能使目標函數值最小。
對於上文中的標准形式，約束矩陣A是一個m*n維矩陣，且m<n，所以一定可以從A中找到一個滿秩m*m矩陣。這個矩陣就稱作矩陣A的一個基陣，矩陣A就可以寫作 [B N] , 相應的解 x 也可以寫成 x=(xB,xN)’，那么 Ax=b 就變為
，左式兩端同乘B矩陣的逆，得到
。由此引出下列名詞：
基陣：非奇異矩陣（滿秩矩陣、可逆矩陣）B
基向量：基陣B由m個線性無關的向量組成，稱之為基向量
基變量：向量xB各分量，與基向量對應的xB中的m個分量成為基變量
非基變量：向量xN各分量
基本解：令xN各分量為0，由
得到的解
稱為基陣B對應的基本解
基本可行解：當
成立時，稱基本解為基本可行解，因為只有滿足所有分量不小於0，才符合標准形式中的約束條件（最后一條
）。
單純形表
如上圖所示，在做單純形表時，我們要找到這么一個滿秩矩陣B，而且要通過行變換把它化為單位陣，同時把這個單位陣上方對應的向量c中元素變為0。也就是說，在標准的單純形表中，在表的第一行中，基變量對應的元素全為0，非基變量對應的元素稱之為檢驗數。這時便找到了此問題的一個基本可行解，此時單純形表最右邊一列的各數從上到下為此基本解對應的目標函數值f和基本解的基變量的值b’（非基變量為0）。
舉個例子：
| S |
X1(非基變量) |
X2(非基變量) |
X3 |
X4 |
X5 |
|
| f |
5（檢驗數） |
10（檢驗數） |
0 |
0 |
0 |
0（目標函數值） |
| X3 |
1/14 |
1/7 |
1 |
0 |
0 |
1（基變量X3值） |
| X4 |
1/7 |
1/12 |
0 |
1 |
0 |
1（基變量X4值） |
| X5 |
1 |
1 |
0 |
0 |
1 |
8（基變量X5值） |
這是一個已經變換好的單純形表，紅色部分是b’，也就是此時基本可行解中基變量分別為X3,X4,X5，他們的值分別為1,1,8，對應的基本可行解就是（0,0,1,1,8）。可以看出，此時目標函數值為0。
單純形方法基本步驟：
初始的單純形表已經給出了線性規划問題的一個基本可行解，接下來要做的就是判斷這個解是否是最優解，如果是的話不用繼續找啦，如果不是的話就找一個比這個解更好的解，再進行判斷，直到找到最優解。但有的問題是沒有最優解的，所以還要判定問題是否無解。
1） 將所給的線性規划問題化為標准形式。
2）找出一個初始的基陣B，作出單純形表，作為程序的輸入（A,C,f,b’）。
3）測試所有的檢驗數，記錄檢驗數中的正數，若全部小於等於0，則已經找到最優解，計算終止。否則轉至4）。
4）測試所有為正的檢驗數，若在單純性表中，其所在的列中其他元素全部小於等於0，則此問題無最優解，計算終止，否則轉至5）。
5）找出檢驗數中的最大值，此最大值元素所在列為A(i,:)，令約束條件中約束向量b與A(i,:)的比值為向量r，向量r中為正的最小值為h，計算過程如下圖。
| S單純形表 |
X1 |
X2 |
X3 |
X4 |
X5 |
|
| f |
5 |
10（最大值） |
0 |
0 |
0 |
0 |
| X3 |
1/14 |
1/7 |
1 |
0 |
0 |
1 |
| X4 |
1/7 |
1/12 |
0 |
1 |
0 |
1 |
| X5 |
1 |
1 |
0 |
0 |
1 |
8 |
表格中黃色部分組成的向量點除（對應元素相除）紅色部分，得到向量（7,12,8），那么7就是我們要找的那個元素，此時記錄元素大小h和坐標(i,j)，注意是在S表中行號和列號，此處是2和2（如果有多個相等的最小值則任取一個即可）。
這個元素1/7就是所謂的轉軸元（或稱基本元），找到他之后要圍繞他進行一系列的行變換，稱之為換基。步驟如下：
①使轉軸元變為1，方法很簡單，就是讓本行所有元素同時除以轉軸元1/7。
②把轉軸元所在列的其他元素都變為0，做法是通過一個循環，遍歷每一個行（自身所在行除外），每行中與轉軸元同列的元素為a，令每行減去轉軸元所在行的第a倍即可。轉至3）。
理論部分到此為止，如有疏漏歡迎留言（參考書目為黃平的《最優化理論》）
matlab仿真程序編寫:
Simplex.m
%Simplex Method function [x,y]=Simplex(f,A,b) %輸入f是檢驗數的數組，1*n維 %輸入A是約束矩陣， m*n維 %輸入b是約束向量， 1*m維 %輸出x是解向量 %輸出y是最優解 %判斷輸入維數是否相符 %做初始單純形表 format rat %將結果以分數表示 S=[f 0;A b']; [n,m]=size(S); %判斷檢驗數 r<=0 r=find(f>0); len=length(r); %有大於0的檢驗數則進入循環 while(len) %檢查非負檢驗數所對列向量元素是否都小於等於0 for k=1:length(r) d=find(S(:,r(k))>0); if(length(d)+1==2) error('無最優解！！！') %break; end end %找到檢驗數中最大值 [Rk,j]=max(S(1,1:m-1)); %b與最大值所在列的比值 rb=S(2:n,m)./S(2:n,j); %把比值中的負數都變無窮 for p=1:length(rb) if(rb(p)<0) rb(p)=Inf; end end [h,i]=min(rb);%列向量比值最小值 % i+1為轉軸元行號(在S中)，j為轉軸元列號 i=i+1; %進行換基,轉軸元置1 S(i,:)=S(i,:)./S(i,j); %轉軸元所在列其他元素都置0 for k=1:n if(k~=i) S(k,:)=S(k,:)-S(i,:)*S(k,j); end end %判斷檢驗數 r<=0 r=find(S(1,1:m-1)>0); len=length(r); end %檢驗數全部非正，找到最優解 %非基變量置0 x=zeros(1,m-1); for i=1:m-1 %找到基變量 j=find(S(:,i)==1);%每列中1的個數 k=find(S(:,i)==0);%每列中0的個數 if((length(j)+1==2)&&(length(k)+1==n)) %i為基本元列號，j是行號 x(i)=S(j,m); end end y=S(1,m);%最優解 S end
下面附帶一個測試程序(test.m)：
clear all clc f=[4 3 0 0 0]; A=[1 1 1 0 0; 1 2 0 1 0; 3 2 0 0 1]; b=[50 80 140]; [x,y]=Simplex(f,A,b) f=[1 1 0 0]; A=[-2 1 1 0; 1 -1 0 1]; b=[4 2]; [x,y]=Simplex(f,A,b)
仿真結果如下：
