#include <iostream.h>
#include <math.h>
#include <stdlib.h>
void main()
{
int n;
float *X,*Y,*V,*R;
cout<<"请输入结点个数：";
cin>>n;
X=new float[n];
Y=new float[n];
V=new float[n];
R=new float[n];
cout<<"请输入各结点的横坐标：";
for (int i=0;i<n;i++)
cin>>X[i];
cout<<"请输入各结点的纵坐标：";
for (i=0;i<n;i++)
cin>>Y[i];
cout<<"请输入各结点的供货量或需求量：";
for (i=0;i<n;i++)
cin>>V[i];
cout<<"请输入各结点的内/外向运输费率：";
for (i=0;i<n;i++)
cin>>R[i];
float sumx=0.0,sumy=0.0,sumqr=0.0;
float avex,avey;
for (i=0;i<n;i++)
{
sumx+=V[i]*R[i]*X[i];
sumy+=V[i]*R[i]*Y[i];
sumqr+=V[i]*R[i];
}
avex=sumx/sumqr;
avey=sumy/sumqr;
cout<<"初始仓库地址为："<<"X="<<avex<<' '<<"Y="<<avey<<endl;
float *d;
d=new float[i];
float xzsumx=0.0,xzsumy=0.0,xzsumqr=0.0; //xz表示修正
float xzavex,xzavey;
int j=0;
while (1)
{
for (i=0;i<n;i++)
d[i]=sqrt((avex-X[i])*(avex-X[i])+(avey-Y[i])*(avey-Y[i]));
for (i=0;i<n;i++)
{
xzsumx+=V[i]*R[i]*X[i]/d[i];
xzsumy+=V[i]*R[i]*Y[i]/d[i];
xzsumqr+=V[i]*R[i]/d[i];
}
xzavex=xzsumx/xzsumqr;
xzavey=xzsumy/xzsumqr;
if (j==50)
break;
else
{
cout<<"第"<<j+1<<"次迭代："<<"tempX="<<xzavex<<' '<<"tempY="<<xzavey<<endl;
avex=xzavex;
avex=xzavey;
xzsumx=0.0;
xzsumy=0.0;
xzsumqr=0.0;
j++;
}
}
cout<<"仓库最佳地址为："<<"X="<<xzavex<<' '<<"Y="<<xzavey<<endl;
}
对不起，发的时候忘记写了：
这个程序运行结果有误，是在后面的迭代部分，请大家看一看错误在哪里
公式：x0 = ( ∑ xiwi ) / ( ∑ wi )
y0 = ( ∑ yiwi ) / ( ∑ wi )
单设施选址模型有不同的名称，如：精确重心法、网格法和重心法等，因为选址因素只包括运输费率和该点的货物运输量，所以此方法较为简单。
我们以预计仓库点位到达各个目的点位的成本最低为考察目标：
即：
Min TC=∑ViRidi i 遍历所有目的地 …（0）
其中： TC——总运输成本；
Vi——i 点运输量；
Ri——到 i 点的运输费率；
di——从待定仓库位置到i 点的距离。
在坐标平面中设待定仓库坐标为（X0,Y0）
∑ViRiXi/di
X0= ...(1）
∑ViRi/di
∑ViRiYi/di
Y0= …（2）
∑ViRi/di
其中di可以坐标间公式
di =SQR((X0-Xi)2+(Y0-Yi)2) …（3）
该方法求解过程：
1、确定各目的地点的坐标，同时确定各点货物运输量和直线距离运费；
2、不考虑距离因素，用重心公式估算初始选址点：
∑ViRiXi
X0= …（4）
∑ViRi
∑ViRiYi
Y0 …（5）
∑ViRi
3、根据公式（3），用步骤2得到的（X0,Y0）计算di；
4、将di代入公式（1），（2），得到修正的（X0,Y0）坐标；
5、根据修正的（X0,Y0）坐标，再重新计算di；
6、 重复步骤4和步骤5直至（X0,Y0）坐标在连续迭代过程中都不再变化或变化在误差范围内；
7、最后，如果需要，利用公式（0）计算最优选址的总成本。
#include <stdio.h>
#include <math.h>
struct xuanzhi
{
double x,y,v,r;
}xz[5]={{3,8,2000,0.050},{8,2,3000,0.050},{2,5,2500,0.075},{6,4,1000,0.075},{8,8,1500,0.075}};
main()
{
struct xuanzhi;
int i;
double d[5];
double x0=0.0,y0=0.0,min=0.0,TC=0.0,t=0.0,m=0.0,w=0.0,x1=0.0,y1=0.0,m1=0.0,m2=0.0,w2=0.0,t1=0.0,m3=0.0,w3=0.0;
for(i=0;i<5;i++)
{
m+=xz[i].v*xz[i].r*xz[i].x;
w+=xz[i].v*xz[i].r;
m1+=xz[i].v*xz[i].r*xz[i].y;
}
x0=m/w;
y0=m1/w;
for(i=0;i<5;i++)
{
d[i]=sqrt((xz[i].x-x0)*(xz[i].x-x0)+(xz[i].y-y0)*(xz[i].y-y0));
m3+=xz[i].v*xz[i].r*xz[i].x/d[i];
t1+=xz[i].v*xz[i].r*xz[i].y/d[i];
w3+=xz[i].v*xz[i].r/d[i];
}
x1=m3/w3;
y1=t1/w3;
do
{
x0=x1;
y0=y1;
for(i=0;i<5;i++)
{
d[i]=sqrt((xz[i].x-x0)*(xz[i].x-x0)+(xz[i].y-y0)*(xz[i].y-y0));
m2+=xz[i].v*xz[i].r*xz[i].x/d[i];
t+=xz[i].v*xz[i].r*xz[i].y/d[i];
w2+=xz[i].v*xz[i].r/d[i];
}
x1=m2/w2;
y1=t/w2;
m2=0.0;
t=0.0;
w2=0.0;
}while(fabs(x0-x1)-0.00000001>0&&fabs(y0-y1)-0.00000001>0);
printf("x0=%12lf,y0=%12lf\n",x0,y0);
for(i=0;i<5;i++)
{
d[i]=sqrt((xz[i].x-x0)*(xz[i].x-x0)+(xz[i].y-y0)*(xz[i].y-y0));
min=xz[i].v*xz[i].r*d[i];
TC+=min;
}
printf("min TC=%12lf\n",TC);
getch();
}
