//---------------------------------------------------------------------------
#include "gomsolv.h"
#include "fruc.h"
#include <Classes.hpp>

//---------------------------------------------------------------------------
#pragma package(smart_init)

void Jstep( PSimplexTable st, PIntVector basic,
            int n, int nb, int ncut,
            bool EachStep, TStrings *S )
{
 int i,j;
 int col,row;
 int haveneg;
 TFruc f;

 do{
  // Найти отрицательный элемент в 1ом столбце
  haveneg = 0;
  for( i=1; (i <= n)&&(!haveneg); i++ )
   haveneg = (st[(nb+1)*i].sign() < 0);
  row = i-1;

  if(haveneg)
  {
   // Ищем положительное отношение
   for( i=1; i <= nb; i++ )
   {
    f = st[(nb+1)*row+i];
    if( f.x == 0 )continue;
    if( ((-st[i])/f).sign() >= 0 )break;
   }
   if( i > nb )throw Exception( "Целевая функция неограничена" );

   // Ищем направляющий элемент в строке row
   TFruc min = (-st[i])/st[(nb+1)*row+i];

   col = i++;
   for( ; i <= nb; i++ )
   {
    if ( st[(nb+1)*row+i].sign() < 0 )
    {
     f = st[(nb+1)*row+i];
     if( f.x == 0 )continue;
     f = (-st[i])/f;
     if( f < min )
     {
      min = f;
      col = i;
     }
    }
   }

   if( EachStep )
    if( ncut > 0 )PrintTable( st, basic, n-1, nb, ncut, row, col, S );else
     PrintTable( st, basic, n, nb, 0, row, col, S );

   f = st[(nb+1)*row+col];

   // Делаем шаг Жорданова исключения
   for( i = 0; i <= n; i++ )
    for( j = 0; j <= nb; j++ )
    {
     if( i == row || j == col )continue;

     st[(nb+1)*i+j] = st[(nb+1)*i+j] - (st[(nb+1)*i+col]*st[(nb+1)*row+j])/f;
    }

   for( i = 0; i <= n; i++ )
    st[(nb+1)*i+col] = st[(nb+1)*i+col]/(-f);

   for( i = 0; i <= nb; i++ )
    st[(nb+1)*row+i] = 0;
   st[(nb+1)*row+col] = -1;

   if( ncut > 0 )
   {
    basic[col-1] = ncut;
    ncut = 0;
    n--;
   }else basic[col-1] = row;
  }
 }while(haveneg);
}

void gomory( PSimplexTable st, PBoolVector Z,
             int n, int nb, int alg,
             bool EachStep, TStrings *S )
{
 int i,j;
 int col,row;

 int haveneg;
 int haveint;
 int NVar=n;
 TFruc f;

 PIntVector basic = new int[nb];
try
{
 for( i = 0; i < nb; i++ )basic[i] = i+1;

 do{
  // Найти отрицательный элемент в 1й строке
  haveneg = 0;
  for( i=1; (i <= nb)&&(!haveneg); i++ )
   haveneg = (st[i].sign() < 0);
  col = i-1;

  if(haveneg)
  {
   // Ищем в столбце col положительный элемент
   for( i=1; (i <= n)&&(st[(nb+1)*i+col].sign() <= 0); i++ );
   if( i > n )throw Exception( "Система ограничений несовместна" );
   row = i;
   for( i=row+1; (i <= n)&&(st[(nb+1)*i+col].sign() <= 0); i++ );
   if( i <= n )row=i;

   // Ищем направляющий элемент в строке row
   TFruc min = (-st[col])/st[(nb+1)*row+col];

   for( i=1; i <= nb; i++ )
   {
    f = st[(nb+1)*row+i];
    if( f.x == 0 )continue;
    f = (-st[i])/f;
    if( (f.sign() > 0)&&(f < min) )
    {
     min = f;
     col = i;
    }
   }

   if( EachStep )
    PrintTable( st, basic, n, nb, 0, row, col, S );

   // Делаем шаг Жорданового исключения
   f = st[(nb+1)*row+col];
   for( i = 0; i <= n; i++ )
    for( j = 0; j <= nb; j++ )
    {
     if( i == row || j == col )continue;
     st[(nb+1)*i+j] = st[(nb+1)*i+j] - (st[(nb+1)*i+col]*st[(nb+1)*row+j])/f;
    }

   for( i = 0; i <= n; i++ )
    st[(nb+1)*i+col] = (-st[(nb+1)*i+col])/f;

   for( i = 0; i <= nb; i++ )
    st[(nb+1)*row+i] = 0;
   st[(nb+1)*row+col] = -1;

   basic[col-1] = row;
  }
 }while(haveneg);

 Jstep( st, basic, n, nb, 0, EachStep, S );
 do{
// Проверка решения на целочисленность
  haveint = 1;
  for( i = 0; i < nb; i++ )haveint &= Z[i];
  if( haveint )haveint = (st[0].x%st[0].y == 0);else haveint = 1;

  for( i = 1; (i <= nb)&&(haveint); i++ )
  {
   f = st[(nb+1)*i];
   haveint = (Z[i-1]&&(f.x%f.y == 0))||(!Z[i-1]);
  }
  i--;

  if( !haveint )
  {
   // Отсечение
   switch(alg)
   {
    case algGomoryI:
     for( j = 0; j <= nb; j++ )
      st[(n+1)*(nb+1)+j]=-(!(st[(nb+1)*i+j]));
     break;
    case algGomoryII:
     st[(n+1)*(nb+1)]=-(!(st[(nb+1)*i]));
     for( j = 1; j <= nb; j++ )
     {
      if( (basic[j-1]<=nb&&Z[basic[j-1]]) ) /*basic[j-1] > n*/
       //Для целочисленных переменных
       if((!(st[(nb+1)*i+j]))>(!(st[(nb+1)*i])))
        st[(n+1)*(nb+1)+j]=-(!(st[(nb+1)*i]))/(-(!(st[(nb+1)*i]))+1)*(-(!st[(nb+1)*i+j])+1);
       else
        st[(n+1)*(nb+1)+j]=-(!(st[(nb+1)*i+j]));
      else
       //Для нецелочисленных переменных
       if( st[(nb+1)*i+j].sign() < 0 )
        st[(n+1)*(nb+1)+j]=-(!(st[(nb+1)*i]))/(-(!(st[(nb+1)*i]))+1)*(-st[(nb+1)*i+j]);
       else
        st[(n+1)*(nb+1)+j]=-(st[(nb+1)*i+j]);
     }
   }
   NVar++;

   Jstep( st, basic, n+1, nb, NVar, EachStep, S );
  }
 }while(!haveint);

 if( EachStep )
    PrintTable( st, basic, n, nb, 0, -1, -1, S );

 PrintSolve( st, nb, S );
}
__finally
{
 delete[] basic;
}
}

void PrintTable( PSimplexTable st, PIntVector basic, int n, int nb, int ncut,
                 int r, int s, TStrings *S )
{
 int i,j;
 AnsiString str1;

 str1 = AnsiString("        1    ");
 for( i=0; i<nb; i++ )
  str1 += Format( "   -x%d    ", ARRAYOFCONST(((int)basic[i])) );
 S->Add( str1 );

 for( i=0; i<=n; i++ )
 {
  str1 = Format( "x%d ", ARRAYOFCONST(((int)i)) );
  for( j=0; j<=nb; j++ )
   str1 += Format( " %4d/%-3d ", ARRAYOFCONST(((int)st[i*(nb+1)+j].x,(int)st[i*(nb+1)+j].y)) );
  S->Add( str1 );
 }

 if( ncut > 0 )
 {
  S->Add( "Отсечение" );
  str1 = Format( "x%d ", ARRAYOFCONST(((int)ncut)) );
  for( j=0; j<=nb; j++ )
   str1 += Format( " %4d/%-3d ", ARRAYOFCONST(((int)st[(n+1)*(nb+1)+j].x,(int)st[(n+1)*(nb+1)+j].y)) );
  S->Add( str1 );
 }

 if( r >= 0 && s >= 0 )
  S->Add( Format( "Разрешающий элемент (%d, %d)", ARRAYOFCONST(((int)r, (int)s)) ) );
  
 S->Add( "" );
}

void PrintSolve( PSimplexTable st, int nb, TStrings *S )
{
 AnsiString str;

 S->Add( "Ответ" );
 S->Add( Format( " целевая функция = %d/%d", ARRAYOFCONST(((int)st[0].x,(int)st[0].y)) ) );

 str = " (";
 for( int i=1; i<=nb; i++ )
 {
  str += Format( "%d/%d", ARRAYOFCONST(((int)st[i*(nb+1)].x,(int)st[i*(nb+1)].y)) );
  if( i < nb )str += ",";
 }
 str += ")";
 S->Add( str );
}
