算法及编程语言 - 声振论坛 - 振动,动力学,声学,信号处理,故障诊断

frogfish 发表于 2005-11-16 15:31

单纯型法解线性规划问题（两阶段法）

/************************************************************************* 单纯型法解线性规划问题（两阶段法） 编程环境：VC++6.0 方程组输入说明： 变量非负，按提示输入相关参数。 *************************************************************************/ #include <stdio.h> #include <stdlib.h> #define MAX 100 #define STP 100
 int stop=1; //迭代记数变量 int status; //iterative迭代返回值：1唯一最优，0无界解，-1无穷多最优解 -2迭代超过限制次数 int step=1; //目前阶段
double a,b,c,temp_c,max=0; //方程组相关系数 int num_x; //变量个数 int num_st; //约束方程数 int num_ar=0; //人工变量个数 int arti; //人工变量下标 int base; //基变量下标 int ma_mi; //1为求最大值，2为求最小值
void create(); //建立方程组 void iterative(); //单纯型法迭代 void output(); //输出结果 void banner(); //打印程序标题 void exchange(); //交换两阶段价值系数 void show(); //输出方程组
void main() { int i,j,k; banner(); create(); //保存原价值系数，转换为第一阶段价值系数 for(i=1;i<=num_x;i++) { k=0; for(j=1;j<=num_ar;j++) if(i==arti) k=1; if(k==1) temp_c=-1; else temp_c=0; } exchange(c,temp_c);
printf("\n\n第一阶段问题为：\n\n"); show(); step++; printf("\n\n按回车开始第一阶段迭代"); getchar(); getchar(); iterative(); if(status==-2) { puts("迭代超过限制次数强行终止！\n"); puts("\n按回车结束"); getchar(); exit(0); } output();
if(max!=0) { puts("\n\n原问题无可行解。\n"); puts("\n按回车结束"); getchar(); exit(0); }
//转换为第二阶段价值系数 exchange(c,temp_c); //把人工变量列全设为0 for(i=1;i<=num_ar;i++) { c]=0; for(j=1;j<=num_st;j++) a]=0; }
puts("\n\n第二阶段问题为：\n\n"); show(); puts("\n\n按回车开始第二阶段迭代"); getchar(); iterative(); switch(status) { case 1: output(); puts("\n\n原问题有唯一最优解。\n"); puts("\n按回车结束"); getchar(); exit(0); case 0: puts("\n\n原问题为无界解。\n"); puts("\n按回车结束"); getchar(); exit(0); case -1: output(); puts("\n\n原问题有无穷多最优解。\n"); puts("\n按回车结束"); getchar(); exit(0); case -2: puts("迭代超过限制次数强行终止！\n"); puts("\n按回车结束"); getchar(); exit(0); }//switch }
void banner() { printf("\t\t****************************************\n"); printf("\t\t 单纯型法解线性规划问题\n"); printf("\t\t 作者：Thunder\n"); printf("\t\t****************************************\n"); printf("\n"); }
void show() { //对方程组以自然的格式输出，系数为零的x不显示 //为1的不显示系数1，-1系数只显示负号 int i,j,k; switch(step) { case 1: printf("min z= "); printf("x[%d]",arti); for(i=2;i<=num_ar;i++) printf(" + x[%d]",arti); break; case 2: printf("max z= "); printf("%lg x[%d]",c,1); for(i=2;i<=num_x;i++) { if(c==1) printf(" + x[%d]",i); else if(c==-1) printf(" - x[%d]",i); else if(c>=0) printf(" +%lg x[%d]",c,i); else printf(" %lg x[%d]",c,i); } break; }
printf("\nst:\n"); for(i=1;i<=num_st;i++) { k=0; for(j=1;j<=num_x;j++) { if(a!=0) { if(a==1&&k!=0) printf(" + x[%d]",j); else if(a==1&&k==0) printf("x[%d]",j); else if(a==-1) printf(" - x[%d]",j); else if(a>=0&&k!=0) printf(" +%lg x[%d]",a,j); else if(a>=0&&k==0) printf(" %lg x[%d]",a,j); else printf(" %lg x[%d]",a,j); k=1; } } printf(" == %lg\n",b); } printf(" x~x[%d]>=0",num_x); }
void exchange() { int i; double temp; for(i=1;i<=num_x;i++) { temp=temp_c; temp_c=c; c=temp; } }
void create() { //输入方程组系数，每个方程输完后回显确认 int i,j,k,re_st,tnum_x,num_addv=0,num_ba=0; char confirm; while(1) { printf("请选择：1、求最大值，2、求最小值：(1/2)"); scanf("%d",&ma_mi); if(ma_mi!=1&&ma_mi!=2) printf("输入错误，重新选择。"); else break; } while(1) { printf("指定变量个数："); scanf("%d",&num_x); printf("输入价值系数c1-c%d:\n",num_x); for(i=1;i<=num_x;i++) { printf("c%d=",i); scanf("%lf",&c); } if(ma_mi==1) printf("max z= "); else printf("min z= "); printf("%lg x[%d]",c,1); for(i=2;i<=num_x;i++) { if(c>=0) printf(" +%lg x[%d]",c,i); else printf("%lg x[%d]",c,i); } printf("\n正确吗？：(y/n)"); getchar(); confirm=getchar(); if (confirm=='y') break; else if(confirm=='n') continue; }
printf("输入约束方程组个数："); scanf("%d",&num_st); for(i=1;i<=num_st;i++) { printf("st.%d:\n",i); while(1) { printf("请选择：1、==，2、>=，3、<= ：(1/2/3)"); scanf("%d",&re_st); if(re_st!=1&&re_st!=2&&re_st!=3) printf("输入错误，请重新选择。"); else break; } printf("输入技术系数：\n"); for(j=1;j<=num_x;j++) { printf("a%d=",j); scanf("%lf",&a); } printf("输入资源拥有量：\nb%d=",i); scanf("%lf",&b); printf("st.%i:\n",i); printf("%lg x[%d]",a,1); for(j=2;j<=num_x;j++) { if(a>=0) printf(" +%lg x[%d]",a,j); else printf(" %lg x[%d]",a,j); } switch(re_st) { case 1: printf(" == %lg",b); break; case 2: printf(" >= %lg",b); break; case 3: printf(" <= %lg",b); break; }
while(1) { printf("\n正确吗？(y/n)"); getchar(); confirm=getchar(); if (confirm=='y') break; else if(confirm=='n') {i-=1; break;} } }
//显示输入的方程组 printf("\n原问题为：\n\n"); if(ma_mi==1) printf("max z= "); else printf("min z= "); printf("%lg x[%d]",c,1); for(i=2;i<=num_x;i++) { if(c==1) printf(" + x[%d]",i); else if(c==-1) printf(" - x[%d]",i); else if(c>=0) printf(" +%lg x[%d]",c,i); else printf(" %lg x[%d]",c,i); }
printf("\nst:\n"); for(i=1;i<=num_st;i++) { k=0; for(j=1;j<=num_x;j++) { if(a!=0) { if(a==1&&k!=0) printf(" + x[%d]",j); else if(a==1&&k==0) printf("x[%d]",j); else if(a==-1) printf(" - x[%d]",j); else if(a>=0&&k!=0) printf(" +%lg x[%d]",a,j); else if(a>=0&&k==0) printf(" %lg x[%d]",a,j); else printf(" %lg x[%d]",a,j); k=1; } } switch(re_st) { case 1: printf(" == %lg\n",b); break; case 2: printf(" >= %lg\n",b); break; case 3: printf(" <= %lg\n",b); break; } } printf(" x~x[%d]>=0\n",num_x);
tnum_x=num_x; for(i=1;i<=num_st;i++) { switch(re_st) { case 1: case 3: num_x+=1; break; case 2: num_x+=2; break; } }
//化为标准形式 if(ma_mi==2) for(i=1;i<=tnum_x;i++) c*=-1; //求最小值时，系数变相反数 for(i=1;i<=num_st;i++) { switch(re_st) { case 1: num_addv++; num_ba++; num_ar++; c=0; base=arti=tnum_x+num_addv; for(j=tnum_x+1;j<=num_x;j++) if(j==tnum_x+num_addv) a=1; else a=0; break; case 2: num_addv++; c=0; num_addv++; num_ba++; num_ar++; c=0; base=arti=tnum_x+num_addv; for(j=tnum_x+1;j<=num_x;j++) if(j==tnum_x+num_addv-1) a=-1; else if(j==tnum_x+num_addv) a=1; else a=0; break; case 3: num_addv++; num_ba++; c=0; base=tnum_x+num_addv; for(j=tnum_x+1;j<=num_x;j++) if(j==tnum_x+num_addv) a=1; else a=0; break; }//switch }//增加松弛变量、剩余变量、人工变量、确定基变量
//显示标准化后的方程组 printf("\n化为标准形式后：\n\n"); if(ma_mi==1) printf("max z= "); else printf("max z'= "); printf("%lg x[%d]",c,1); for(i=2;i<=num_x;i++) { k=0; for(j=1;j<=num_ar;j++) if(i==arti) k=1; if(k==1) printf(" -M x[%d]",i); else if(c==1) printf(" + x[%d]",i); else if(c==-1) printf(" - x[%d]",i); else if(c>=0) printf(" +%lg x[%d]",c,i); else printf(" %lg x[%d]",c,i); }
printf("\nst:\n"); for(i=1;i<=num_st;i++) { k=0; for(j=1;j<=num_x;j++) { if(a!=0) { if(a==1&&k!=0) printf(" + x[%d]",j); else if(a==1&&k==0) printf("x[%d]",j); else if(a==-1) printf(" - x[%d]",j); else if(a>=0&&k!=0) printf(" +%lg x[%d]",a,j); else if(a>=0&&k==0) printf(" %lg x[%d]",a,j); else printf(" %lg x[%d]",a,j); k=1; } } printf(" == %lg\n",b); } printf(" x~x[%d]>=0",num_x); }
void iterative() { int i,j,k,k_a,k_f,l; //k_a,k_f值为0或1，记录当前下标在arti[]或base[]里的搜索结果 int base_elem; int base_out,base_in; double sigma,temp; double value_be; //高斯消元里保存主元素值
printf("\n\n第%d次迭代：\n\n",stop); for(i=1;i<=num_st;i++) { printf("c%d=%lg\t",base,c]); printf("b%d=%lg\t",i,b);
switch(step) { case 1: for(j=1;j<=num_x;j++) { printf("a[%d][%d]=%lg\t",i,j,a); } printf("\n"); break; case 2: for(j=1;j<=num_x;j++) { k_a=0; for(l=1;l<=num_ar;l++) if(j==arti)k_a=1; if(k_a!=1) printf("a[%d][%d]=%lg\t",i,j,a); } printf("\n"); break; } } //求检验数sigma for(i=1;i<=num_x;i++) { sigma=c; for(j=1;j<=num_st;j++) sigma-=c]*a; for(j=1;j<=num_st;j++) if(i==base) sigma=0; switch(step) { case 1: printf("sigma[%d]=%lg\t",i,sigma); break; case 2: k_a=0; for(l=1;l<=num_ar;l++) if(i==arti) k_a=1; if(k_a!=1) printf("sigma[%d]=%lg\t",i,sigma); break; } } putchar('\n'); //检验检验数sigma是否全小于等于0 k=0; for(i=1;i<=num_x;i++) { if(sigma>0) k=1; } if(k==0) { //sigma是全小于等于0时，检查是否为无穷多最优解 for(i=1;i<=num_x;i++) { k_f=k_a=0; for(j=1;j<=num_ar;j++) if(i==arti) k_a=1; if(sigma==0&&k_a!=1) { for(j=1;j<=num_st;j++) if(i==base) k_f=1; if(k_f==0) {status=-1; return;} } } status=1; return; } //检查是否为无界解 for(i=1;i<=num_x;i++) { k_f=0; if(sigma>0) { for(j=1;j<=num_st;j++) if(a>0) k_f=1; if(k_f!=1) {status=0; return;} } }
//确定换入变量 for(i=1;i<=num_x;i++) { k=0; for(j=1;j<=num_st;j++) if(i==base) k=1; if(k==0&&sigma>0) temp=sigma-1; }//temp赋初值 for(i=1;i<=num_x;i++) { k=0; for(j=1;j<=num_st;j++) if(i==base) k=1; if(k==0) if(sigma>temp&&sigma>0) { base_in=i; temp=sigma; } }
//确定换出变量 for(i=1;i<=num_st;i++) if(a>0) { temp=b/a+1; break; }//temp赋初值 for(i=1;i<=num_st;i++) { if(b/a<=temp&&a>0) { for(j=1;j<=num_ar;j++) if(base==arti) { base_out=base; base_elem=i; temp=b/a; break; } }//人工变量优先换出 if(b/a<temp&&a>0) { base_out=base; base_elem=i; temp=b/a; } }
printf(" 基变量："); for(i=1;i<=num_st;i++) printf("x[%d] ",base); printf("换入变量：x[%d] 换出变量：x[%d]",base_in,base_out); //基变量变换，进行新方程初始化后迭代 for(i=1;i<=num_st;i++) { if(base==base_out) base=base_in; } //初始化主元素行系数 value_be=a; b/=value_be; for(i=1;i<=num_x;i++) a/=value_be;
for(i=1;i<=num_st;i++) { if(i!=base_elem) { b-=b*a; value_be=a; for(j=1;j<=num_x;j++) a-=a*value_be; } } stop++; if(stop>STP) {status=-2; return;} iterative(); }
void output() { int i,j; double X; printf("\n结果如下：\n"); printf("\nX=("); for(i=1;i<=num_x;i++) { for(j=1;j<=num_st;j++) if(i==base) {X=b;break;} else X=0; printf("%lg ",X); } printf(")"); for(i=1;i<=num_x;i++) max+=c*X; if(ma_mi==1) printf("\nMax z= %lf\n",max); else printf("\nMin z= %lf\n",-max); }

幻觉发表于 2006-7-20 18:42

调试过吗

页: [1]

声振论坛's Archiver

单纯型法解线性规划问题（两阶段法）