单纯形解决线形规划问题的代码（c）
总算完成了，下面是代码。
又：因为运筹学并没有学过，这次凑巧帮管院一个同学写单纯形的代码，
所以就有了这个。是从昨天开始看书的，如果有错误，还请见谅------觉得那个运筹学还真不容易呢，虽比不上咱们的离散：）~~~~~~~~
测试结果：
Successful:
X[3]=77.142860 X[1]=7.142854 X[2]=28.571430
max Z=585.714294
Press any key to continue
/*
*data.txt
*/
5 2 1 0 0
2 3 0 1 0
1 5 0 0 1
170 100 150
10 18 0 0 0
1 0 0
0 1 0
0 0 1
1 0 0
0 1 0
0 0 1
0 0 0
2 3 4
/*
本题的例子是：
max Z=10*x1+18*x2;
s.t. 5*x1+2*x2<=170
2*x1+3*x2<=100
x1+5*x2<=150
x1,x2>=0
本文件的数据是上式经过加入松弛变量后的标准形式的数据
*/
/*
*fmain.cpp
*/
#include <stdio.h>
#include <stdlib.h>
#define M 3
#define N 5
/* 解决形如
max z=cx
s.t. Ax=b,
x>=0
的线形规划问题
_B[M][M]是B[M][M]的逆阵
*/
float A[M][N],b[M],
c[N],B[M][M],
cn[M],X[M],Y[M],
Kao[N],P[M],
E[M][M],_B[M][M];
/*
解序号
*/
int Xcount[M];
float sita;
/*
进基变量，出基变量
*/
int xin,xout;
int init(FILE *pf);
void computeE();
void compute_B();
void computeXY();
int computeKao();
int computePsita();
void computeBcn();
int succ();
int fall();
void main()
{
FILE *pf;
float result=0.0;
if(!init(pf))
{
printf("Error:can't open data file,please check...\n");
exit(-1);
}
computeXY();
xin=computeKao();
while(!succ())
{
xout=computePsita();
if(fall())
{
printf("Fall...already exit...\n");
exit(0);
}
computeBcn();
computeE();
compute_B();
computeXY();
xin=computeKao();
}
printf("Successful:\n");
for(int i=0;i<M;i++)
printf("X[%d]=%f\t",Xcount[i]+1,X[i]);
for(int i=0;i<M;i++)
result+=c[Xcount[i]]*X[i];
printf("\nmax Z=%f\n",result);
}
/*
用data.txt文件初始化，包括：
A,b,c,B,_B,cn,Xcount
*/
int init(FILE *pf)
{
pf=fopen("data.txt","r");
if(pf==NULL){
fclose(pf);
return 0;
}
/*input A[M][N]*/
for(int i=0;i<M;i++)
for(int j=0;j<N;j++){
fscanf(pf,"%f",&A[i][j]);
}
/*input b[M]*/
for(int i=0;i<M;i++){
fscanf(pf,"%f",&b[i]);
}
/*input c[N]*/
for(int i=0;i<N;i++){
fscanf(pf,"%f",&c[i]);
}
/*input B[M][M]*/
for(int i=0;i<M;i++)
for(int j=0;j<M;j++)
fscanf(pf,"%f",&B[i][j]);
/*input _B[M][M]*/
for(int i=0;i<M;i++)
for(int j=0;j<M;j++){
fscanf(pf,"%f",&_B[i][j]);
}
/*input cn[M]*/
for(int i=0;i<M;i++)
fscanf(pf,"%f",&cn[i]);
/* input Xcount[M] */
for(int i=0;i<M;i++){
fscanf(pf,"%d",&Xcount[i]);
}
/*其它*/
for(int i=0;i<M;i++)
P[i]=X[i]=Y[i]=0;
for(int i=0;i<N;i++)
Kao[i]=0;
fclose(pf);
return 1;
}
/*
计算初等变换矩阵E
*/
void computeE()
{
for(int i=0;i<M;i++)
for(int j=0;j<M;j++)
{
E[i][j]=0;
if(i==j)
E[i][j]=1;
}
for(int i=0;i<M;i++)
E[i][xout]=-P[i]/P[xout];
E[xout][xout]=1/P[xout];
}
/*
计算B逆矩阵_B
*/
void compute_B()
{
float __B[M][M];
for(int i=0;i<M;i++)
for(int j=0;j<M;j++)
__B[i][j]=0;
for(int i=0;i<M;i++)
for(int j=0;j<M;j++)
for(int k=0;k<M;k++)
__B[i][j]+=E[i][k]*_B[k][j];
for(int i=0;i<M;i++)
for(int j=0;j<M;j++){
_B[i][j]=__B[i][j];
}
}
/*
计算 X=_B*b ,Y=_B*cn
*/
void computeXY()
{
for(int i=0;i<M;i++)
P[i]=X[i]=Y[i]=0;
for(int i=0;i<M;i++){
for(int j=0;j<M;j++)
X[i]+=_B[i][j]*b[j];
}
for(int i=0;i<M;i++){
for(int j=0;j<M;j++)
Y[i]+=cn[j]*_B[j][i];
}
}
/* 计算检验数 kao=c-cn*_B*A
xin<-return
*/
int computeKao()
{
float midval=0.0;
int k;
for(int i=0;i<N;i++){
for(int j=0;j<M;j++)
midval+=Y[j]*A[j][i];
Kao[i]=c[i]-midval;
midval=0;
}
midval=Kao[0];
k=0;
for(int j=0;j<N;j++)
if(Kao[j]>midval){
midval=Kao[j];
k=j;
}
return k;
}
/*计算进基列向量P,判断出基变量
xout<-return
*/
int computePsita()
{
int k=0;
for(int i=0;i<M;i++){
for(int j=0;j<M;j++)
P[i]+=_B[i][j]*A[j][xin];
}
sita=X[0]/P[0];
for(int j=0;j<M;j++)
if(sita>X[j]/P[j])
{
sita=X[j]/P[j];
k=j;
}
return k;
}
/*
计算矩阵B,cn
*/
void computeBcn()
{
for(int i=0;i<M;i++)
B[i][xout]=A[i][xin];
cn[xout]=c[xin];
Xcount[xout]=xin;
}
/*
判断是否成功
*/
int succ()
{
int k=1;
for(int i=0;i<N;i++)
k=k&&(Kao[i]<=0.0);
return k;
}
/*
算法失败，退出
*/
int fall()
{
int k=1;
for(int i=0;i<N;i++)
k=k&&(P[i]<=0.0);
return k;
}
