19565数学建模附录.docx

19565数学建模附录.docx

《19565数学建模附录.docx》由会员分享，可在线阅读，更多相关《19565数学建模附录.docx（16页珍藏版）》请在冰点文库上搜索。

19565数学建模附录

附录：

（1）：

数据的准备：

原始数据文件a:

\input.dat如下：

16118181191711819124201212112222123

23124241252412626127261282812929130

30131291322913333134331353513737136

37138381391104223013375044606551946

62057720188680994801010300111122012

12210131342014145001523163217460019

51020652171022731238122494226107028

1110301210341362351411036143037152038

15203916450181811501917801819110024

20306212119522222023232022424120025

2472026266902726520282817029298830

3046231296903229160333370343332035

351603737703637290383830393013033

301903434260353510036

数据准备程序：

a:

\mind\mind.cpp如下：

/******************************************************/

/*filename：

a:

\mind\mind.cpp*/

/*功能：

为解决钢管订购问题准备数据，数据存放在*/

/*c:

\cij0.dat（解问题

（1）所需数据）中和*/

/*c:

\cij1.dat（解问题（3）所需数据）中*/

/*说明：

本程序必须在执行主程序前执行*/

/*作者：

丁林阁*/

/******************************************************/

#include

#include

#include

#include

#defineSIZE39

#defineINFILENAME"a:

\\input.dat"

#defineOUTFILENAME1"c:

\\cij0.dat"

#defineOUTFILENAME2"c:

\\cij1.dat"

#defineINF100000

staticdoublecc[SIZE][SIZE];

staticdoubleisrail[SIZE][SIZE];

staticintisend=0;

intmind（intl,intn,double*d,double*c）;//最小路求解函数

intcost（intr）;//铁路费用计算函数

voidmain（）

{

doublec[SIZE*SIZE],c1[SIZE][SIZE],c2[SIZE][SIZE];

doubledd[SIZE*SIZE];

intn,i,j,rtn,m;

ints[7]={23,25,27,32,31,36,39};

inta[2][21],counter,flag,l;

intvalue;

charfil1[]=OUTFILENAME1;

FILE*fp,*fp1;

flag=54+23;

l=15;

if（（fp=fopen（INFILENAME,"r"））==NULL）//打开原始数据文input.dat

{

printf（"can'topeninputfile!

!

!

\n"）;

exit

（1）;

}

for（i=0;i<21;i++）

{

a[1][i]=i+1;

a[0][i]=i+1;

}

a[1][15]=26;a[1][16]=30;a[1][17]=33;a[1][18]=34;a[1][19]=35;a[1][20]=36;

for（m=0;m<2;m++）//m=0则生成cij0.dat,m=1则生成cij1.dat

{

if（（fp1=fopen（fil1,"w"））==NULL）

{

printf（"can'topeninputfile!

!

!

\n"）;

exit

（1）;

}

counter=0;//读取数据计数器以区分图

（1）数据和图

（2）数据

n=SIZE;

for（i=0;i

{

for（j=0;j

{

cc[i][j]=（double）（（i==j）?

0:

INF）;

c1[i][j]=（double）（（i==j）?

0:

INF）;

c2[i][j]=（double）（（i==j）?

0:

INF）;

israil[i][j]=0;

dd[i*SIZE+j]=0;

}

}

fseek（fp,0,SEEK_SET）;//读input.dat

while（!

feof（fp）&&counter

{

fscanf（fp,"%d",&i）;

fscanf（fp,"%d",&value）;

fscanf（fp,"%d",&j）;

if（value==1）//value=1表示i到j由铁路连接

{

israil[i-1][j-1]=1;

counter++;

continue;

}

if（israil[i-1][j-1]==1）

{

c1[i-1][j-1]=value;

c1[j-1][i-1]=value;

}

else

{

c2[i-1][j-1]=value;

c2[j-1][i-1]=value;

}

counter++;

}

for（i=0;i

{

for（j=0;j

{

c[i*n+j]=c1[i][j];

}

}

rtn=mind（1,SIZE,c,dd）;//生成铁路完全图

for（i=0;i

{

for（j=0;j

{

c1[i][j]=（double）cost（（int）c[i*n+j]）;

}

}

for（i=0;i

{

for（j=0;j

{

c[i*n+j]=c2[i][j];

dd[i*SIZE+j]=0;

}

}

rtn=mind（1,SIZE,c,dd）;//生成公路完全图

for（i=0;i

{

for（j=0;j

{

c2[i][j]=c[i*n+j]*0.1;

}

}

for（i=0;i

{

for（j=0;j

{

c[i*n+j]=（c1[i][j]

c1[i][j]:

c2[i][j]）;//合并成运费完全图

dd[i*SIZE+j]=0;

}

}

rtn=mind（1,SIZE,c,dd）;//生成最小运费完全图

for（i=0;i

{

for（j=0;j

{

cc[i][j]=c[i*n+j];

}

}

for（i=0;i<7;i++）

{

for（j=0;j

{

fprintf（fp1,"%f",cc[s[i]-1][a[m][j]-1]）;

}

fprintf（fp1,"\n"）;

}

fclose（fp1）;

fil1[6]++;

flag=58+23;

l=21;

}

fclose（fp）;

}

intmind（intl,intn,double*d,double*c）

{

inti,ii,j,k,n1;

intj1;

doublex;

if（l!

=0）

{

n1=n-1;

for（i=0;i

{

d[i*n+i]=0;

ii=i+1;

for（j=ii;j

d[i*n+j]=d[j*n+i];

}

}/*endofif*/

for（i=0;i

{

c[i*n+i]=0;

}

j1=0;

if（l!

=1）n1=n;

for（k=0;k

{

for（i=0;i

{

if（l==1）j1=i+1;

for（j=j1;j

{

x=d[i*n+k]+d[k*n+j];

c[i*n+j]=x

x:

d[i*n+j];

if（l==1）c[j*n+i]=c[i*n+j];

}

}

for（i=0;i

{

for（j=0;j

{

d[i*n+j]=c[i*n+j];

}

}

}/*endoffor*/

return（0）;

}

intcost（intr）

{

if（r>=0&&r<=300）return（20）;

else

if（r>300&&r<=350）return（23）;

else

if（r>350&&r<=400）return（26）;

else

if（r>400&&r<=450）return（29）;

else

if（r>450&&r<=500）return（32）;

else

if（r>500&&r<=600）return（37）;

else

if（r>600&&r<=700）return（44）;

else

if（r>700&&r<=800）return（50）;

else

if（r>800&&r<=900）return（55）;

else

if（r>900&&r<=1000）return（60）;

else

if（（r-1000）%100>0）

return（60+（（r-1000）/100+1）*5）;

else

return（60+（r-1000）/100*5）;

}

（2）：

解问题

（1）的程序：

主程序：

a:

\solut2\getmin.m如下：

clear

globalcij;

globalb;

globals;

globalcounter;

globalminc;

globalminx;

getcij;

loadcij;

counter=0;

b=[80080010002000200020003000];

s=[104301750606194205201680480300220210420500];

%x0=rand（7,15）*234;

x0=ones（7,15）*234;

options=foptions;

options（14）=3301;

x=constr（'cost',x0,options）;

minx

minc

数据处理函数a:

\solut2\getcij.m如下：

functionf=getcij（）

fid=0;

whilefid<1

filename='c:

\cij0.dat';

[fid,message]=fopen（filename,'r'）;

iffid==-1

disp（message）;

end

end

cij=fscanf（fid,'%g'）;

fclose（fid）;

p=[160155155160155150160];

cij=cij';

%cij=[cij（1:

15）+p

（1）;cij（16:

30）+p

（2）;cij（31:

45）+p（3）;cij（46:

%60）+p（4）;cij（61:

75）+p（5）;cij（76:

90）+p（6）;cij（91:

105）+p（7）];

cij=[cij（1:

15）+p

（1）;cij（16:

30）+p

（2）;cij（31:

45）+p（3）;…

cij（46:

60）+p（4）;cij（61:

75）+p（5）;cij（76:

90）+p（6）;…

cij（91:

105）+p（7）];

savecijcij;

目标函数a:

\solut2\cost.m如下：

function[f,g]=cost（x）

globalcij;

globalb

globals

globalcounter;

globalminc;

globalminx;

counter=counter+1;

a=sum（x）;

%k=sum（x'）;

%m=find（k~=0）;

l

（1）=0;r

（1）=0;

fori=2:

15

l（i）=s（i-1）-r（i-1）;

r（i）=a（i）-l（i）;

end

f=sum（sum（cij.*x））+0.05*sum（l.*（l-1）+r.*（r+1））

g=-[x;l;r;zeros（1,15）];

%ll=length（m）;

%g（10,1:

2*ll）=[k（m）-b（m）500-k（m）];

g（10,1:

4）=-[800-sum（x（1,:

））800-sum（x（2,:

））1000-sum（x（3,:

））sum（x（7,:

））-500];

ifcounter>3300

minc=f;

minx=round（x）;

end

（3）：

解问题（3）的程序：

主程序：

a:

\solut3\getmin.m如下：

clear

globalcij;

globalb;

globals;

globalcounter;

globalminc;

globalminx;

get2cij;

loadcij;

counter=0;

b=[80080010002000200020003000];

s=[1043017506061942052016804803002202104205004210130190260100];

%x0=rand（7,21）*234;

%x0=ones（7,21）*234;

x0=zeros（7,21）;

options=foptions;

options（14）=3601;

x=constr（'cost',x0,options）;

minx（6,20）=minx（6,20）+1;

minx（5,6）=minx（5,6）+1;

minx（7,6）=minx（7,6）+1;

minx（3,8）=minx（3,8）+1;

minx（1,7）=minx（1,7）+1;

counter=0;

minc=cost（minx）;

minx

minc

数据处理函数a:

\solut3\get2cij.m如下：

functionf=get2cij（）

fid=0;

whilefid<1

filename='c:

\cij1.dat';

[fid,message]=fopen（filename,'r'）;

iffid==-1

disp（message）;

end

end

cij=fscanf（fid,'%g'）;

fclose（fid）;

p=[160155155160155150160];

cij=cij';

%cij=[cij（1:

15）+p

（1）;cij（16:

30）+p

（2）;cij（31:

45）+p（3）;cij（46:

60）+p（4）;cij（61:

75）+p（5）;cij（76:

90）+p（6）;cij（91:

105）+p（7）];

cij=[cij（1:

21）+p

（1）;cij（22:

42）+p

（2）;cij（43:

63）+p（3）;cij（64:

84）+p（4）;cij（85:

105）+p（5）;cij（106:

126）+p（6）;cij（127:

147）+p（7）];

%fori=1:

7

%cij（（i-1）*15+1:

（i-1）*15+15）=cij（（i-1）*15+1:

（i-1）*15+15）+p（i）;

%end

savecijcij;

目标函数a:

\solut3\cost.m如下：

function[f,g]=cost（x）

globalcij;

globalb

globals

globalcounter;

globalminc;

globalminx;

counter=counter+1;

a=sum（x）;

%k=sum（x'）;

%m=find（k~=0）;

l

（1）=0;r

（1）=a

（1）;

fori=2:

8

l（i）=s（i-1）-r（i-1）;

r（i）=a（i）-l（i）;

end

l（9）=s（8）-r（8）;r（16）=0;l（16）=a（16）;u

（1）=s（15）-l（16）;r（9）=a（9）-l（9）-u

（1）;

l（10）=s（9）-r（9）;r（10）=a（10）-l（10）;l（11）=s（10）-r（10）;r（15）=0;l（15）=a（15）;

fori=14:

-1:

12

r（i）=s（i）-l（i+1）;

l（i）=a（i）-r（i）;

end

r（11）=s（11）-l（12）;u

（2）=a（11）-l（11）-r（11）;l（17）=s（16）-u

（2）;r（18）=0;

l（18）=a（18）;u（3）=s（17）-l（18）;r（17）=a（17）-l（17）-u（3）;l（19）=s（18）-r（17）;

r（19）=a（19）-l（19）;l（20）=s（19）-r（19）;r（20）=a（20）-l（20）;l（21）=s（20）-r（20）;

r（21）=a（21）-l（21）;

f=sum（sum（cij.*x））+0.05*sum（l.*（l-1）+r.*（r+1））+0.05*sum（u.*（u-1））

g=-[x;l;r;zeros（1,21）];

g（10,1:

7）=-[u800-sum（x（1,:

））800-sum（x（2,:

））1000-sum（x（3,:

））sum（x（7,:

））-500];

ifcounter>3600

minc=f;

minx=round（x）;

end