數(shù)學(xué)建模常用算法程序

1、假設(shè)圖G權(quán)的鄰接矩陣為A,0aa，a1112aaa2122.,aaan1n21n2n0nn來(lái)存放各邊長(zhǎng)度，其中：iii=1,2,n;ijijij是i,j之間邊的長(zhǎng)度，i,j=1,2,對(duì)于無(wú)向圖，A是對(duì)稱(chēng)矩陣，ijjiiji,j之間沒(méi)有邊，在程序中以各邊都不可能達(dá)到的充分大的數(shù)代替；Floyd算法的基本思想是:遞推產(chǎn)生一個(gè)矩陣序列A,A，A，A，其中A(i,j)表01knk示從頂點(diǎn)v到頂點(diǎn)v的路徑上所經(jīng)過(guò)的頂點(diǎn)序號(hào)不大于k的最短路徑長(zhǎng)度。ij計(jì)算時(shí)用迭代公式:A(i,j)=min(A(i,j),A(i,k)+A(k,j)kk-1k-1k-1k是迭代次數(shù)，i,j,k=1,2,n。最后，當(dāng)k=n時(shí)，

2、A即是各頂點(diǎn)之間的最短通路值。n例10用Floyd算法求解例1。矩陣path用來(lái)存放每對(duì)頂點(diǎn)之間最短路徑上所經(jīng)過(guò)的頂點(diǎn)的序號(hào)。Floyd算法的Matlab程序如下:clear;clc;M=10000;a(1,:)=0,50,M,40,25,10;a(2,:)=zeros(1,2),15,20,M,25;a(3,:)=zeros(1,3),10,20,M;a(4,:)=zeros(1,4),10,25;a(5,:)=zeros(1,5),55;a(6,:)=zeros(1,6);b=a+a;path=zeros(length(b);fork=1:6fori=1:6forj=1:6ifb(i,j)

3、b(i,k)+b(k,j)b(i,j)=b(i,k)+b(k,j);path(i,j)=k;endendendendb,pathFloyd最短路算法的MATLAB程序%floyd.m%采用floyd算法計(jì)算圖a中每對(duì)頂點(diǎn)最短路%d是矩離矩陣%r是路由矩陣functiond,r=floyd(a)n=size(a,1);d=a;fori=1:nforj=1:nr(i,j)=j;endendrfork=1:nfori=1:nforj=1:nifd(i,k)+d(k,j)d(i,j)d(i,j)=d(i,k)+d(k,j);r(i,j)=r(i,k)endendendkdrend兩個(gè)指定頂點(diǎn)之間的最短

4、路徑問(wèn)題如下：給出了一個(gè)連接若干個(gè)城鎮(zhèn)的鐵路網(wǎng)絡(luò)，在這個(gè)網(wǎng)絡(luò)的兩個(gè)指定城鎮(zhèn)間，找一條最短鐵路線。以各城鎮(zhèn)為圖G的頂點(diǎn)，兩城鎮(zhèn)間的直通鐵路為圖G相應(yīng)兩頂點(diǎn)間的邊，得圖G。對(duì)G的每一邊e，賦以一個(gè)實(shí)數(shù)w(e)直通鐵路的長(zhǎng)度，稱(chēng)為e的權(quán)，得到賦權(quán)圖G。G的子圖的權(quán)是指子圖的各邊的權(quán)和。問(wèn)題就是求賦權(quán)圖G中指定的兩個(gè)頂點(diǎn)u,v間的具最小00權(quán)的軌。這條軌叫做u,v間的最短路，它的權(quán)叫做u,v間的距離，亦記作d(u,v)。000000求最短路已有成熟的算法：迪克斯特拉(Dijkstra)算法，其基本思想是按距u從近到0遠(yuǎn)為順序，依次求得u到G的各頂點(diǎn)的最短路和距離，直至v(或直至G的所有頂點(diǎn))，00算法

5、結(jié)束。為避免重復(fù)并保留每一步的計(jì)算信息，采用了標(biāo)號(hào)算法。下面是該算法。令l(u)=0，對(duì)v豐u，令l(v)=g，S=u，i=0。0000對(duì)每個(gè)veS(S=VS)，用iiminl(v),l(u)+w(uv)ueSi代替l(v)。計(jì)算minl(v)，把達(dá)到這個(gè)最小值的一個(gè)頂點(diǎn)記為u，令S=SUu。i+1i+1ii+1veSi.若i=IVI-1，停止；若iIVI-1，用i+1代替i，轉(zhuǎn)(ii)。算法結(jié)束時(shí)，從u到各頂點(diǎn)v的距離由v的最后一次的標(biāo)號(hào)l(v)給出。在v進(jìn)入S之0i前的標(biāo)號(hào)l(v)叫T標(biāo)號(hào)，v進(jìn)入S時(shí)的標(biāo)號(hào)l(v)叫P標(biāo)號(hào)。算法就是不斷修改各項(xiàng)點(diǎn)的Ti標(biāo)號(hào)，直至獲得P標(biāo)號(hào)。若在算法運(yùn)行過(guò)程

6、中，將每一頂點(diǎn)獲得P標(biāo)號(hào)所由來(lái)的邊在圖上標(biāo)明，則算法結(jié)束時(shí)，u至各項(xiàng)點(diǎn)的最短路也在圖上標(biāo)示出來(lái)了。0例9某公司在六個(gè)城市c,c,c中有分公司，從c到c的直接航程票價(jià)記在下述126ij矩陣的(i,j)位置上。(g表示無(wú)直接航路)，請(qǐng)幫助該公司設(shè)計(jì)一張城市c到其它城市間的1票價(jià)最便宜的路線圖。050g4025105001520g25g1501020g4020100102525g20100551025g25550用矩陣a(n為頂點(diǎn)個(gè)數(shù))存放各邊權(quán)的鄰接矩陣，行向量pb、index、index、nxn12d分別用來(lái)存放P標(biāo)號(hào)信息、標(biāo)號(hào)頂點(diǎn)順序、標(biāo)號(hào)頂點(diǎn)索引、最短通路的值。其中分量1當(dāng)?shù)趇頂點(diǎn)已標(biāo)號(hào)pb

7、(i)=;0當(dāng)?shù)趇頂點(diǎn)未標(biāo)號(hào)index(i)存放始點(diǎn)到第i點(diǎn)最短通路中第i頂點(diǎn)前一頂點(diǎn)的序號(hào);2d(i)存放由始點(diǎn)到第i點(diǎn)最短通路的值。求第一個(gè)城市到其它城市的最短路徑的Matlab程序如下：clear;clc;M=10000;a(1,:)=0,50,M,40,25,10;a(2,:)=zeros(1,2),15,20,M,25;a(3,:)=zeros(1,3),10,20,M;a(4,:)=zeros(1,4),10,25;a(5,:)=zeros(1,5),55;a(6,:)=zeros(1,6);a=a+a;pb(1:length(a)=0;pb(1)=1;index1=1;index

8、2=ones(1,length(a);d(1:length(a)=M;d(1)=0;temp=1;whilesum(pb)=2index=index(1);endindex2(temp)=index;endd,index1,index24.2.1prim算法構(gòu)造最小生成樹(shù)設(shè)置兩個(gè)集合P和Q，其中P用于存放G的最小生成樹(shù)中的頂點(diǎn)，集合Q存放G的最小生成樹(shù)中的邊。令集合P的初值為P=v（假設(shè)構(gòu)造最小生成樹(shù)時(shí)，從頂點(diǎn)v出發(fā)），11集合Q的初值為Q=。prim算法的思想是，從所有peP，veV-P的邊中，選取具有最小權(quán)值的邊pv，將頂點(diǎn)v加入集合P中，將邊pv加入集合Q中，如此不斷重復(fù)，直到P=V時(shí)，

9、最小生成樹(shù)構(gòu)造完畢，這時(shí)集合Q中包含了最小生成樹(shù)的所有邊。prim算法如下：(i)P=v，Q=；1(ii)whileP=Vpv=minW,peP,veV一PpvP=P+vQ=Q+pvend例11用prim算法求右圖的最小生成樹(shù)。Matlab我們用result的第一、二、三行分別表示生成樹(shù)邊的起點(diǎn)、終點(diǎn)、權(quán)集合3xn程序如下：clc;clear;M=1000;a(1,2)=50;a(1,3)=60;a(2,4)=65;a(2,5)=40;a(3,4)=52;a(3,7)=45;a(4,5)=50;a(4,6)=30;a(4,7)=42;a(5,6)=70;a=a;zeros(2,7);a=a+a

10、;a(find(a=0)=M;result=;p=1;tb=2:length(a);whilelength(result)=length(a)-1temp=a(p,tb);temp=temp(:);d=min(temp);jb,kb=find(a(p,tb)=d);j=p(jb(1);k=tb(kb(1);result=result,j;k;d;p=p,k;tb(find(tb=k)=;endresult4.2.1Kruskal算法構(gòu)造最小生成樹(shù)科茹斯克爾(Kruskal)算法是一個(gè)好算法。Kruskal算法如下：選eeE(G)，使得w(e)=min。11(ii)若e,e，e已選好，則從E(G

11、)-e,e，,e中選取e，使得12i12ii+1Ge,e,e,e中無(wú)圈，且12ii+1w(e)=min。i+1(iii)直到選得e為止。V-1例12用Kruskal算法構(gòu)造例3的最小生成樹(shù)。我們用index存放各邊端點(diǎn)的信息，當(dāng)選中某一邊之后，就將此邊對(duì)應(yīng)的頂點(diǎn)序號(hào)中xn較大序號(hào)u改記為此邊的另一序號(hào)v，同時(shí)把后面邊中所有序號(hào)為u的改記為v。此方法的幾何意義是：將序號(hào)u的這個(gè)頂點(diǎn)收縮到v頂點(diǎn)，u頂點(diǎn)不復(fù)存在。后面繼續(xù)尋查時(shí)，發(fā)現(xiàn)某邊的兩個(gè)頂點(diǎn)序號(hào)相同時(shí)，認(rèn)為已被收縮掉，失去了被選取的資格。Matlab程序如下：clc;clear;M=1000;a(1,2)=50;a(1,3)=60;a(2,4

12、)=65;a(2,5)=40;a(3,4)=52;a(3,7)=45;a(4,5)=50;a(4,6)=30;a(4,7)=42;a(5,6)=70;i,j=find(a=0)&(a=M);b=a(find(a=0)&(a=M);data=i;j;b;index=data(1:2,:);loop=max(size(a)-1;result=;whilelength(result)v2index(find(index=v1)=v2;elseindex(find(index=v2)=v1;enddata(:,flag)=;index(:,flag)=;endresult一名推銷(xiāo)員準(zhǔn)備前往若干城市推銷(xiāo)

13、產(chǎn)品，然后回到他的出發(fā)地。如何為他設(shè)計(jì)一條最短的旅行路線(從駐地出發(fā)，經(jīng)過(guò)每個(gè)城市恰好一次，最后返回駐地)？這個(gè)問(wèn)題稱(chēng)為旅行商問(wèn)題。用圖論的術(shù)語(yǔ)說(shuō)，就是在一個(gè)賦權(quán)完全圖中，找出一個(gè)有最小權(quán)的Hamilton圈。稱(chēng)這種圈為最優(yōu)圈。與最短路問(wèn)題及連線問(wèn)題相反，目前還沒(méi)有求解旅行商問(wèn)題的有效算法。所以希望有一個(gè)方法以獲得相當(dāng)好(但不一定最優(yōu))的解。一個(gè)可行的辦法是首先求一個(gè)Hamilton圈C，然后適當(dāng)修改C以得到具有較小權(quán)的另一個(gè)Hamilton圈。修改的方法叫做改良圈算法。設(shè)初始圈C=vvvv。12n1對(duì)于1i+1jn，構(gòu)造新的Hamilton圈：C=vvvvvvvvvvv,j12ijj-1j-

14、2i+1j+1j+2n1它是由C中刪去邊vv和vv,添加邊vv和vv而得到的。若ii+1jj+1iji+1j+1w(vv)+w(vv)0flag=0;form=1:L-3forn=m+2:L-1ifa(c1(m),c1(n)+a(c1(m+1),c1(n+1)0flag=0;form=1:L-3forn=m+2:L-1ifa(c1(m),c1(n)+a(c1(m+1),c1(n+1).a(c1(m),c1(m+1)+a(c1(n),c1(n+1)flag=1;c1(m+1:n)=c1(n:-1:m+1);endendendendsum1=0;fori=1:L-1sum1=sum1+a(c1(i

15、),c1(i+1);endifsum1sumsum=sum1;circle=c1;endcircle,sum已知非線性整數(shù)規(guī)劃為：Maxz=x2+x2+3x2+4x2+2x2123458x2x3xx2x123450 x99(i=1,，5)ix+x+x+x+x400TOC o 1-5 h z12345s.t.x+2x+2x+x+6x800123452x+x+6x200123x+x+5x20045對(duì)該題，目前尚無(wú)有效方法求出準(zhǔn)確解。如果用顯枚舉法試探，共需計(jì)算(100)5=1010個(gè)點(diǎn)，其計(jì)算量非常之大。然而應(yīng)用蒙特卡洛去隨機(jī)計(jì)算106個(gè)點(diǎn)，便可找到滿意解，那么這種方法的可信度究竟怎樣呢？下面就分

16、析隨機(jī)取樣采集106個(gè)點(diǎn)計(jì)算時(shí)，應(yīng)用概率理論來(lái)估計(jì)一下可信度。不失一般性，假定一個(gè)整數(shù)規(guī)劃的最優(yōu)點(diǎn)不是孤立的奇點(diǎn)。假設(shè)目標(biāo)函數(shù)落在高值區(qū)的概率分別為0.01，0.00001，則當(dāng)計(jì)算106個(gè)點(diǎn)后，有任一個(gè)點(diǎn)能落在高值區(qū)的概率分別為1-0.9910000000.9999(100多位)，1-0.9999910000000.999954602。解(i)首先編寫(xiě)M文件mente.m定義目標(biāo)函數(shù)f和約束向量函數(shù)g,程序如下：functionf,g=mengte(x);f=x(l廠2+x(2廠2+3*x(3廠2+4*x(4廠2+2*x(5)-8*x(l)-2*x(2)-3*x(3).-x(4)-2*x(5);g(1)=sum(x)-400;g(2)=x(1)+2*x(2)+2*x(3)+x(4)+6*x(5)-800;g(3)