/*本程序是单纯形法,参考书目:河海大学<吴凤平>运筹学方法与应用*/
#include<stdio.h>
#include<math.h>
#define X 5
#define Y 7
void xi_max(int *m2,int *mn1,float *c,int *is,int *ir,int *j0,float (*a)[X][Y])
{
  int j;
  *c=0;
  for(j=1;j<=*is;j++) 
    {
      if((*a)[*ir][j]-*c>0)
    {
      *c=(*a)[*ir][j];
      *j0=j;
    }
     }
}      /***************** 参数说明 **********************/
      /*   m_约束方程个数(基变量个数),n_非基变量个数   */
      /*   m2-m+2整个变量,(*a)[X][Y]存放初始数据       */
      /*   (*k)[]存放基变量脚标,(*x)存放基变量最优值   */
      /********* http://happyyangxu.home.sunbo.net*******/ 
int xi_sm(int m,int n,int m2,int mn1,int l1,float (*a)[X][Y],
       int(*k)[],float(*x)[])
{
  int i,m1,mn,j0,i0,j;
  float c,g;
  m1=m+1;
  mn=m+n;
  for(i=1;i<=m;i++)
    (*k)[i]=n+i;
leap1:
  if(l1-1==0)
     xi_max(&m2,&mn1,&c,&mn,&m1,&j0,a);
  else
     {
leap2:   if(l1-50==0)
      xi_max(&m2,&mn,&c,&n,&m2,&j0,a);
       else
      xi_max(&m2,&mn1,&c,&mn,&m2,&j0,a);
     }
   c-=1e-8;
   if((c<=0)&&(l1==1)&&((*a)[m1][mn1]-1e-8>0))
      {
    printf("\n\t*********Not min&&No solution**********\n");
    return(-1);
      }
   if((c<0)&&(l1==1))
       {
     l1=50;
     goto leap2;
       }
    if(c<=0)
       {
      for(i=1;i<=m;i++)
         (*x)[i]=(*a)[i][mn1];
         printf("\n\t********Optimal solution********\n");
      for(i=1;i<=m;i++)
         printf("\ni=%d,x[i]=%f\n",(*k)[i],(*x)[i]);
      printf("\nF=%f\n",(*a)[m2][mn1]);
      return 0;
    }
     c=1e8;
     for(i=1;i<=m;i++)
       {
     if((*a)[i][j0]>1e-8)
        {
          g=(*a)[i][mn1]/(*a)[i][j0];
          if(g-c<0);
         {
           c=g;
           i0=i;
          }
         }
    }
      if(c==1e+8)
    {
      printf("\n\t*******Lp no solution********\n");
      return -3;
     }
       (*k)[i0]=j0;
       for(j=1;j<=mn1;j++)
      {
        if((j==j0)||(l1==50)&&(n<j)&&(j<mn1))
           continue;
         g=(*a)[i0][j]/(*a)[i0][j0];
           (*a)[i0][j]=g;
         for(i=1;i<=m2;i++)
           {
         if((i==i0)||(!(l1==1)&(i==m1)))
            continue;
         (*a)[i][j]=(*a)[i][j]-(*a)[i][j0]*g;
           }
       }
       for(i=1;i<=m2;i++)
      (*a)[i][j0]=0;
       (*a)[i0][j0]=1;
       goto leap1;
}  
main()
{
  float a[5][7]={{0,0,0, 0,  0,0,0},
                 {0,3,-4,3,  1,0,12},
                 {0,3,0, 6,  0,1,12},
                 {0,0,0, 0,  0,0,0},
                 {0,69,0,144,0,0,300}};
  int k[3];
  float x[3];
  clrscr();
  xi_sm(2,3,4,6,0,&a,&k,&x);
}