運籌學 matlab 表上作業法 運輸問題

文章目錄

• • 運籌學 matlab 表上作業法 運輸問題
  • • 一，實驗內容
    • 二，建立運輸問題的數學模型
    • 三，實驗方法與步驟
    • • 1，定義輸入，利用表上作業法求解
      • 2，判斷是何種運輸問題，即是否產銷平衡
      • 3，求出初始基可行解
      • 4，用位勢法求檢驗數
      • 5，判斷是否是最優解（無負檢驗數）
      • 6，閉回路調整法改進
      • 7，完整matlab實現
    • 四，實驗結果
    • 五，實驗結果分析
    • 六，附錄：部分其余樣例

一，實驗內容

二，建立運輸問題的數學模型

根據題意是求預期盈利最大的采購方案，而我們所學的運輸問題模型是已知運價求最小運價的運輸方案，則我們需要將表中的利潤轉換為運價，需滿足利潤越高的則運價越低，可以定義單位運價表中的運價為表中最大的利潤減去原來的利潤，滿足利潤越高則運價越低。
單位運價表為：

ABCD產量

城市1	0	5	4	3	2500
城市2	2	8	3	4	2500
城市3	1	7	6	2	5000
銷量	1500	2000	3000	3500

$$minz=0x_{11}+5x_{12}+4x_{13}+3x_{14}+2x_{21}+8x_{22}+3x_{23}+4x_{24}+1x_{31}+7x_{32}+6x_{33}+2x_{34}$$

$$?\begin{array}{cc}{x}_{11}+x_{12}+x_{13}+x_{14}& =2500\\ x_{21}+x_{22}+x_{23}+x_{24}& =2500\\ x_{31}+x_{32}+x_{33}+x_{34}& =5000\\ x_{11}+x_{21}+x_{31}& =1500\\ x_{12}+x_{22}+x_{32}& =2000\\ x_{13}+x_{23}+x_{33}& =3000\\ x_{14}+x_{24}+x_{34}& =3500\\ x_{11},x_{12},x_{13},x_{14},x_{21},x_{22},x_{23},x_{24},x_{31},x_{32},x_{33},x_{34}& \ge 0\end{array}$$

三，實驗方法與步驟

1，定義輸入，利用表上作業法求解

function [sigma]=Transport(N,out,in) %N:運價表 out:每個產地產量(按行輸入) in:每個銷地銷量(按行輸入) %x:最優運輸方案

2，判斷是何種運輸問題，即是否產銷平衡

sum_out=sum(out); sum_in=sum(in);if sum_out>sum_in %產大于銷的情況，轉換為產銷平衡問題[old_row,old_col]=size(N);new=zeros(old_row,1);N=[N,new];in=[in,sum_out-sum_in];disp("該問題產大于銷,方案最后一列為虛擬銷地") elseif sum_in>sum_out %銷大于產的情況，轉換為產銷平衡問題[old_row1,old_col1]=size(N);new1=zeros(1,old_col1);N=[N;new1];out=[out,sum_in-sum_out];disp("該問題銷大于產,方案最后一行為虛擬產地") elsedisp("該問題為產銷平衡問題") end

3，求出初始基可行解

%求初始方案(西北角法) for i=1:rowfor j=1:colif out1(i)==0out1(i)=-1;%若該產地產量為0，則置為-1，即劃掉該行endif in1(j)==0in1(j)=-1; %若該銷地銷量為0，則置為-1，即劃掉該列endif out1(i)>0&&in1(j)>0if out1(i)>in1(j) x(i,j)=in1(j); %若產量大于銷量，則該銷地銷量填入方案out1(i)=out1(i)-x(i,j); %該產地剩余的產量in1(j)=in1(j)-x(i,j); %該銷地剩余的銷量else x(i,j)=out1(i); %若銷量大于產量，則直接將產地產量填入方案in1(j)=in1(j)-x(i,j); %該產地剩余的產量out1(i)=out1(i)-x(i,j); %該銷地剩余的銷量break;endendend end

4，用位勢法求檢驗數

%求檢驗數(位勢法) for i=2:row %初始化ui，除了第一個元素為0，其他所有ui元素變為infui(i)=inf; end for i=1:col %初始化vj,所有vj元素變為infvj(i)=inf; end while 1 %求ui,vjchecku=find(ui==inf);checkv=find(vj==inf);if isempty(checku)&&isempty(checkv) %當ui,vj都不為inf時,ui和vj計算完成,跳出循環break;endfor i=1:rowfor j=1:colif x(i,j)~=-1if ui(i)==inf&&vj(j)==inf%若ui,vj全為inf,則先跳過計算continue;elseif vj(j)==infvj(j)=N(i,j)-ui(i);%若vj為inf,ui不為inf,計算vjelseui(i)=N(i,j)-vj(j);%若ui為inf,vj不為inf,計算uiendendend end end for i=1:row %初始化檢驗數表，所有元素變為inffor j=1:colsigma(i,j)=inf; end end for i=1:row %計算檢驗數表for j=1:colif x(i,j)==-1sigma(i,j)=N(i,j)-ui(i)-vj(j); endend end

5，判斷是否是最優解（無負檢驗數）

if sigma>0 disp("有唯一最優方案,最優方案為(表中-1表示空格):")

6，閉回路調整法改進

%找閉回路 while 1 for i=1:row %找第c列中有無未被訪問的點if visit(i,c)==0 %若該點未被訪問,則訪問該點,進行標記并存入路徑表visit(i,c)=1;r=i; %記錄該點行標circle(p,1)=i; circle(p,2)=c;circle(p,3)=x(i,c);p=p+1;break;endendfor j=1:col %找第r行中有無未被訪問的點if visit(r,j)==0 %若該點未被訪問,則訪問該點,進行標記并存入路徑表visit(r,j)=1;c=j; %記錄該點列標circle(p,1)=r;circle(p,2)=j;circle(p,3)=x(r,j);p=p+1;break;endenda=find(visit(r,:)==0);b=find(visit(:,c)==0);a1=find(visit(r,:)==2);b1=find(visit(:,c)==2);if isempty(a)&&isempty(b) %判斷該點所在行和列中有無未被訪問的點if ~isempty(a1)||~isempty(b1) %判斷最后訪問點是否與起始點在同一行或同一列，若是則跳出循環，找到閉回路break;else %若不是,則將該點置為-1(此路不通,下次循環不走此路),重置訪問表和路徑表,開始下一次循環visit(r,c)=-1;for i=1:rowfor j=1:colif visit(i,j)==1visit(i,j)=0;endendendr=r1;c=c1;circle=-ones(row+col,3);circle(1,1)=r1;circle(1,2)=c1;p=2;endend end[rows,cols]=size(circle); %定義circle表時我們給了足夠大的維度,而由于閉回路可能無法包括所有點,現在要將多余行刪去 for i=1:rows if circle(i,1)==-1break;end end circle(i:rows,:)=[];%開始找閉回路中的頂點表 add=circle(1,:); circle=[circle;add]; %將起始點填入circle表,形成完整的閉回路 [row1,col1]=size(circle); i=1;%若閉回路中同一行或同一列中有兩個以上的點,則刪去中間點,只留下頂點 while 1 if i==row1-1break; endi=i+1;if (circle(i-1,1)==circle(i+1,1))circle(i,:)=[];row1=row1-1;i=i-1;endif (circle(i-1,2)==circle(i+1,2))circle(i,:)=[];row1=row1-1; i=i-1;end end%由于起始點在閉回路中出現了兩次,因此刪去第二次出現的起始點 k=find(circle(:,3)==-1); circle(k(2),:)=[];[row2,col2]=size(circle);%將頂點中奇數編號和偶數編號分開 %頂點表的結構:[行標 列標 運量] if mod(row2,2)==0 %若頂點總數是偶數single=zeros(row2/2,3); %定義奇數編號頂點表double=zeros(row2/2,3); %定義偶數編號頂點表 end if mod(row2,2)==1 %若頂點總數是奇數single=zeros(row2/2+0.5,3); %定義奇數編號頂點表double=zeros(row2/2-0.5,3); %定義偶數編號頂點表 endj1=1; j2=1; %將閉回路中的點按編號分別存入奇數頂點表和偶數頂點表 for i=1:row2if mod(i,2)==1single(j1,:)=circle(i,:);j1=j1+1;endif mod(i,2)==0double(j2,:)=circle(i,:);j2=j2+1;end end%更新運量表 value=double(:,3); [min_x,index]=min(value(value>=0));%找到偶數頂點的運量最小值及其位置 x(r1,c1)=min_x; %把該運量最小值填入閉回路起始點的運量 x(double(index,1),double(index,2))=-1; %把該偶數頂點的運量標記為-1 double(index,:)=[]; %將該偶數頂點從偶數頂點表中刪去,以免影響后續計算[row3,col3]=size(single); [row4,col4]=size(double); %將奇數編號頂點的運量加上min_x for i=2:row3x(single(i,1),single(i,2))=single(i,3)+min_x; end %將偶數編號頂點的運量減去min_x for i=1:row4x(double(i,1),double(i,2))=double(i,3)-min_x; end end

7，完整matlab實現

function [sigma]=Transport(N,out,in) %N:運價表 out:每個產地產量(按行輸入) in:每個銷地銷量(按行輸入) %x:最優運輸方案sum_out=sum(out); sum_in=sum(in);if sum_out>sum_in %產大于銷的情況，轉換為產銷平衡問題[old_row,old_col]=size(N);new=zeros(old_row,1);N=[N,new];in=[in,sum_out-sum_in];disp("該問題產大于銷,方案最后一列為虛擬銷地") elseif sum_in>sum_out %銷大于產的情況，轉換為產銷平衡問題[old_row1,old_col1]=size(N);new1=zeros(1,old_col1);N=[N;new1];out=[out,sum_in-sum_out];disp("該問題銷大于產,方案最后一行為虛擬產地") elsedisp("該問題為產銷平衡問題") end[row,col]=size(N); sigma=zeros(row,col);%定義檢驗數表 x=-ones(row,col);%定義初始運輸表 ui=zeros(row,1);%定義ui vj=zeros(1,col);%定義vj out1=out; %求初始方案的運算過程的產量out1 in1=in; %求初始方案的運算過程的銷量in1%求初始方案(西北角法) for i=1:rowfor j=1:colif out1(i)==0out1(i)=-1;%若該產地產量為0，則置為-1，即劃掉該行endif in1(j)==0in1(j)=-1; %若該銷地銷量為0，則置為-1，即劃掉該列endif out1(i)>0&&in1(j)>0if out1(i)>in1(j) x(i,j)=in1(j); %若產量大于銷量，則該銷地銷量填入方案out1(i)=out1(i)-x(i,j); %該產地剩余的產量in1(j)=in1(j)-x(i,j); %該銷地剩余的銷量else x(i,j)=out1(i); %若銷量大于產量，則直接將產地產量填入方案in1(j)=in1(j)-x(i,j); %該產地剩余的產量out1(i)=out1(i)-x(i,j); %該銷地剩余的銷量break;endendend end%迭代過程 while 1 %求檢驗數(位勢法) for i=2:row %初始化ui，除了第一個元素為0，其他所有ui元素變為infui(i)=inf; end for i=1:col %初始化vj,所有vj元素變為infvj(i)=inf; end while 1 %求ui,vjchecku=find(ui==inf);checkv=find(vj==inf);if isempty(checku)&&isempty(checkv) %當ui,vj都不為inf時,ui和vj計算完成,跳出循環break;endfor i=1:rowfor j=1:colif x(i,j)~=-1if ui(i)==inf&&vj(j)==inf%若ui,vj全為inf,則先跳過計算continue;elseif vj(j)==infvj(j)=N(i,j)-ui(i);%若vj為inf,ui不為inf,計算vjelseui(i)=N(i,j)-vj(j);%若ui為inf,vj不為inf,計算uiendendend end end for i=1:row %初始化檢驗數表，所有元素變為inffor j=1:colsigma(i,j)=inf; end end for i=1:row %計算檢驗數表for j=1:colif x(i,j)==-1sigma(i,j)=N(i,j)-ui(i)-vj(j); endend end%判斷是否得到最優方案 if sigma>0 disp("有唯一最優方案,最優方案為(表中-1表示空格):")xsum_min=0;for i=1:row %計算最小運價for j=1:colif x(i,j)~=-1sum_min=sum_min+x(i,j)*N(i,j);endendenddisp("最小運價為：")sum_minbreak; elseif sigma>=0disp("最優方案不唯一,其中一個為(表中-1表示空格):")xsum_min=0;for i=1:row %計算最小運價for j=1:colif x(i,j)~=-1sum_min=sum_min+x(i,j)*N(i,j);endendend disp("最小運價為：")sum_minbreak; end%閉回路調整法 visit=x; for i=1:row %初始化訪問表,可以被訪問的點標為0for j=1:colif visit(i,j)~=-1visit(i,j)=0; endend endm=min(sigma(sigma<0));%找到小于零的最小檢驗數m [r2,c2]=find(sigma==m);%找到m的位置%記錄m的行標和列標,由于可能出現檢驗數相同的情況,我們取其中第一個 r=r2(1); c=c2(1); r1=r2(1); c1=c2(1); visit(r,c)=2; %標記m已被訪問,記為2 circle=-ones(row+col+1,3); %定義閉回路路徑表 %circle表的結構:[行標 列標 運量]%將起點(m點)填入路徑表 circle(1,1)=r1; circle(1,2)=c1; p=2;%找閉回路 while 1 for i=1:row %找第c列中有無未被訪問的點if visit(i,c)==0 %若該點未被訪問,則訪問該點,進行標記并存入路徑表visit(i,c)=1;r=i; %記錄該點行標circle(p,1)=i; circle(p,2)=c;circle(p,3)=x(i,c);p=p+1;break;endendfor j=1:col %找第r行中有無未被訪問的點if visit(r,j)==0 %若該點未被訪問,則訪問該點,進行標記并存入路徑表visit(r,j)=1;c=j; %記錄該點列標circle(p,1)=r;circle(p,2)=j;circle(p,3)=x(r,j);p=p+1;break;endenda=find(visit(r,:)==0);b=find(visit(:,c)==0);a1=find(visit(r,:)==2);b1=find(visit(:,c)==2);if isempty(a)&&isempty(b) %判斷該點所在行和列中有無未被訪問的點if ~isempty(a1)||~isempty(b1) %判斷最后訪問點是否與起始點在同一行或同一列，若是則跳出循環，找到閉回路break;else %若不是,則將該點置為-1(此路不通,下次循環不走此路),重置訪問表和路徑表,開始下一次循環visit(r,c)=-1;for i=1:rowfor j=1:colif visit(i,j)==1visit(i,j)=0;endendendr=r1;c=c1;circle=-ones(row+col,3);circle(1,1)=r1;circle(1,2)=c1;p=2;endend end[rows,cols]=size(circle); %定義circle表時我們給了足夠大的維度,而由于閉回路可能無法包括所有點,現在要將多余行刪去 for i=1:rows if circle(i,1)==-1break;end end circle(i:rows,:)=[]; %開始找閉回路中的頂點表 add=circle(1,:); circle=[circle;add]; %將起始點填入circle表,形成完整的閉回路 [row1,col1]=size(circle); i=1;%若閉回路中同一行或同一列中有兩個以上的點,則刪去中間點,只留下頂點 while 1 if i==row1-1break; endi=i+1;if (circle(i-1,1)==circle(i+1,1))circle(i,:)=[];row1=row1-1;i=i-1;endif (circle(i-1,2)==circle(i+1,2))circle(i,:)=[];row1=row1-1; i=i-1;end end %由于起始點在閉回路中出現了兩次,因此刪去第二次出現的起始點 k=find(circle(:,3)==-1); circle(k(2),:)=[]; [row2,col2]=size(circle); %將頂點中奇數編號和偶數編號分開 %頂點表的結構:[行標 列標 運量] if mod(row2,2)==0 %若頂點總數是偶數single=zeros(row2/2,3); %定義奇數編號頂點表double=zeros(row2/2,3); %定義偶數編號頂點表 end if mod(row2,2)==1 %若頂點總數是奇數single=zeros(row2/2+0.5,3); %定義奇數編號頂點表double=zeros(row2/2-0.5,3); %定義偶數編號頂點表 end j1=1; j2=1; %將閉回路中的點按編號分別存入奇數頂點表和偶數頂點表 for i=1:row2if mod(i,2)==1single(j1,:)=circle(i,:);j1=j1+1;endif mod(i,2)==0double(j2,:)=circle(i,:);j2=j2+1;end end %更新運量表 value=double(:,3); [min_x,index]=min(value(value>=0));%找到偶數頂點的運量最小值及其位置 x(r1,c1)=min_x; %把該運量最小值填入閉回路起始點的運量 x(double(index,1),double(index,2))=-1; %把該偶數頂點的運量標記為-1 double(index,:)=[]; %將該偶數頂點從偶數頂點表中刪去,以免影響后續計算[row3,col3]=size(single); [row4,col4]=size(double); %將奇數編號頂點的運量加上min_x for i=2:row3x(single(i,1),single(i,2))=single(i,3)+min_x; end %將偶數編號頂點的運量減去min_x for i=1:row4x(double(i,1),double(i,2))=double(i,3)-min_x; end end

8，matlab函數輸入

N=[ 0 5 4 32 8 3 41 7 6 2]; in=[1500 2000 3000 3500]; out=[2500 2500 5000]; [sigma]=Transport(N,out,in)

四，實驗結果

結果為：

該問題為產銷平衡問題 有唯一最優方案,最優方案為(表中-1表示空格): x =0 2000 500 -1-1 -1 2500 -11500 -1 -1 3500 最小運價為： sum_min = 28000sigma =Inf Inf Inf 23 4 Inf 4Inf 1 1 Inf

以下是matlab截圖：

轉換為原問題：

最優解為：(預期贏利最大的采購方案)

ABCD

城市1		2000	500
城市2			2500
城市3	1500			3500

最大贏利為：10*(2500+2500+5000)-28000=72000(元)

五，實驗結果分析

經過計算，實驗結果正確，并使用多個樣例測試程序皆有正確結果。

六，附錄：部分其余樣例

%P94運籌學第四版（清華大學出版社） N=[ 3 11 3 101 9 2 87 4 10 5]; in=[3 6 5 6]; out=[7 4 9]; [sigma]=Transport(N,out,in)

總結

以上是生活随笔為你收集整理的运筹学 matlab实现运输问题（表上作业法）的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。